From 48fa2026b12e27fc128fc361da3c741f8e6f4e50 Mon Sep 17 00:00:00 2001 From: Jukka Aho Date: Sun, 2 Sep 2018 11:16:25 +0300 Subject: [PATCH] Make package compatible with Julia v0.7 (#6) --- .travis.yml | 19 +++++---- REQUIRE | 3 +- docs/deploy.jl | 2 +- src/MortarContact2D.jl | 3 +- src/contact2d.jl | 47 ++++++++++------------ src/mortar2d.jl | 44 ++++++++++---------- test/REQUIRE | 2 + test/runtests.jl | 9 +++-- test/test_contact2d.jl | 32 ++++++++------- test/test_contact_segmentation.jl | 4 +- test/test_mortar2d_calculate_projection.jl | 4 +- test/test_mortar2d_ex1.jl | 9 ++--- test/test_mortar2d_ex2.jl | 3 +- 13 files changed, 91 insertions(+), 90 deletions(-) create mode 100644 test/REQUIRE diff --git a/.travis.yml b/.travis.yml index 1879b29..24e3a34 100644 --- a/.travis.yml +++ b/.travis.yml @@ -1,12 +1,17 @@ language: julia + os: - linux + julia: - - 0.6 -before_script: - - julia --color=yes -e 'Pkg.clone("https://github.com/JuliaFEM/PkgTestSuite.jl.git")' - - julia --color=yes -e 'using PkgTestSuite; init()' -script: - - julia --color=yes -e 'using PkgTestSuite; test()' + - 0.7 + - 1.0 + - nightly + +matrix: + allow_failures: + - julia: 1.0 + - julia: nightly + after_success: - - julia --color=yes -e 'using PkgTestSuite; deploy()' + - julia -e 'cd(Pkg.dir("MortarContact2D")); Pkg.add("Coverage"); using Coverage; Coveralls.submit(Coveralls.process_folder())' diff --git a/REQUIRE b/REQUIRE index 67baa34..30b195f 100644 --- a/REQUIRE +++ b/REQUIRE @@ -1,3 +1,2 @@ -julia 0.6 +julia 0.7 FEMBase -Reexport diff --git a/docs/deploy.jl b/docs/deploy.jl index 7b9de9c..41a15df 100644 --- a/docs/deploy.jl +++ b/docs/deploy.jl @@ -5,7 +5,7 @@ using Documenter deploydocs( repo = "github.com/JuliaFEM/MortarContact2D.jl.git", - julia = "0.6", + julia = "0.7", target = "build", deps = nothing, make = nothing) diff --git a/src/MortarContact2D.jl b/src/MortarContact2D.jl index 413ed2b..a6fe34e 100644 --- a/src/MortarContact2D.jl +++ b/src/MortarContact2D.jl @@ -12,8 +12,7 @@ """ module MortarContact2D -using Reexport -@reexport using FEMBase +using FEMBase, SparseArrays, LinearAlgebra, Statistics include("mortar2d.jl") include("contact2d.jl") diff --git a/src/contact2d.jl b/src/contact2d.jl index 032b005..05b3757 100644 --- a/src/contact2d.jl +++ b/src/contact2d.jl @@ -161,7 +161,7 @@ function FEMBase.assemble_elements!(problem::Problem{Contact2D}, assembly::Assem continue end - Ae = eye(nsl) + Ae = Matrix(I, nsl, nsl) if props.dual_basis De = zeros(nsl, nsl) Me = zeros(nsl, nsl) @@ -172,7 +172,7 @@ function FEMBase.assemble_elements!(problem::Problem{Contact2D}, assembly::Assem xi = ip.coords[1] xi_s = dot([1/2*(1-xi); 1/2*(1+xi)], xi1) N1 = vec(get_basis(slave_element, xi_s, time)) - De += w*diagm(N1) + De += w*Matrix(Diagonal(N1)) Me += w*N1*N1' end Ae = De*inv(Me) @@ -293,8 +293,8 @@ function FEMBase.assemble_elements!(problem::Problem{Contact2D}, assembly::Assem C1 = sparse(problem.assembly.C1, ndofs, ndofs) C2 = sparse(problem.assembly.C2, ndofs, ndofs) D = sparse(problem.assembly.D, ndofs, ndofs) - g = full(problem.assembly.g, ndofs, 1) - c = full(problem.assembly.c, ndofs, 1) + g = Vector(sparse(problem.assembly.g, ndofs, 1)[:]) + c = Vector(sparse(problem.assembly.c, ndofs, 1)[:]) for j in S dofs = [2*(j-1)+1, 2*(j-1)+2] @@ -303,7 +303,7 @@ function FEMBase.assemble_elements!(problem::Problem{Contact2D}, assembly::Assem state = problem.properties.contact_state_in_first_iteration if problem.properties.iteration == 1 - info("First contact iteration, initial contact state = $state") + @info("First contact iteration, initial contact state = $state") if state == :AUTO avg_gap = mean([weighted_gap[j][1] for j in S]) @@ -313,7 +313,7 @@ function FEMBase.assemble_elements!(problem::Problem{Contact2D}, assembly::Assem else state = :UNKNOWN end - info("Average weighted gap = $avg_gap, std gap = $std_gap, automatically determined contact state = $state") + @info("Average weighted gap = $avg_gap, std gap = $std_gap, automatically determined contact state = $state") end end @@ -385,41 +385,36 @@ function FEMBase.assemble_elements!(problem::Problem{Contact2D}, assembly::Assem update!(slave_elements, "slip nodes", time => is_slip) end - info("# | A | I | St | Sl | gap | pres | comp") + @info("# | A | I | St | Sl | gap | pres | comp") for j in S str1 = "$j | $(is_active[j]) | $(is_inactive[j]) | $(is_stick[j]) | $(is_slip[j]) | " - str2 = "$(round(weighted_gap[j][1], 3)) | $(round(contact_pressure[j][1], 3)) | $(round(complementarity_condition[j][1], 3))" - info(str1 * str2) + str2 = "$(round(weighted_gap[j][1]; digits=3)) | $(round(contact_pressure[j][1]; digits=3)) | $(round(complementarity_condition[j][1]; digits=3))" + @info(str1 * str2) end # solve variational inequality - # constitutive modelling in tangent direction, frictionless contact - for j in S - dofs = [2*(j-1)+1, 2*(j-1)+2] - if (is_active[j] == 1) && (is_slip[j] == 1) - info("$j is in active/slip, removing tangential constraint $(dofs[2])") - C2[dofs[2],:] = 0.0 - g[dofs[2]] = 0.0 - D[dofs[2], dofs] = tangents[j] - end - end - # remove inactive nodes from assembly for j in S dofs = [2*(j-1)+1, 2*(j-1)+2] + n, t = dofs if is_inactive[j] == 1 - info("$j is inactive, removing dofs $dofs") - C1[dofs,:] = 0.0 - C2[dofs,:] = 0.0 - D[dofs,:] = 0.0 - g[dofs,:] = 0.0 + @debug("$j is inactive, removing dofs $dofs") + C1[dofs,:] .= 0.0 + C2[dofs,:] .= 0.0 + D[dofs,:] .= 0.0 + g[dofs,:] .= 0.0 + elseif (is_active[j] == 1) && (is_slip[j] == 1) + @debug("$j is in active/slip, removing tangential constraint $t") + C2[t,:] .= 0.0 + g[t] = 0.0 + D[t, dofs] .= tangents[j] end end problem.assembly.C1 = C1 problem.assembly.C2 = C2 problem.assembly.D = D - problem.assembly.g = g + problem.assembly.g = sparsevec(g) end diff --git a/src/mortar2d.jl b/src/mortar2d.jl index 1a0f8c1..3e644ad 100644 --- a/src/mortar2d.jl +++ b/src/mortar2d.jl @@ -11,8 +11,8 @@ function Mortar2D() return Mortar2D([]) end -function FEMBase.add_elements!(::Problem{Mortar2D}, ::Any) - error("use `add_slave_elements!` and `add_master_elements!` to add ", +function FEMBase.add_elements!(::Problem{Mortar2D}, ::Union{Tuple, Vector}) + error("use `add_slave_elements!` and `add_master_elements!` to add " * "elements to the Mortar2D problem.") end @@ -79,19 +79,19 @@ function project_from_master_to_slave(slave_element::Element{E}, x2, time; midpnt = mean(x1_) dist = norm(midpnt - x2) distval = dist/len - info("Projection from the master element to the slave failed.") - info("Arguments:") - info("Time: $time") - info("Point to project to the slave element: (x2): $x2") - info("Slave element geometry (x1): $x1_") - info("Slave element normal direction (n1): $n1_") - info("Midpoint of slave element: $midpnt") - info("Length of slave element: $len") - info("Distance between midpoint of slave element and x2: $dist") - info("Distance/Lenght-ratio: $distval") - info("Number of iterations: $max_iterations") - info("Wanted convergence tolerance: $tol") - info("Achieved convergence tolerance: $(norm(dxi1))") + @info("Projection from the master element to the slave failed.") + @info("Arguments:") + @info("Time: $time") + @info("Point to project to the slave element: (x2): $x2") + @info("Slave element geometry (x1): $x1_") + @info("Slave element normal direction (n1): $n1_") + @info("Midpoint of slave element: $midpnt") + @info("Length of slave element: $len") + @info("Distance between midpoint of slave element and x2: $dist") + @info("Distance/Lenght-ratio: $distval") + @info("Number of iterations: $max_iterations") + @info("Wanted convergence tolerance: $tol") + @info("Achieved convergence tolerance: $(norm(dxi1))") error("Failed to project point to slave surface in $max_iterations.") end @@ -234,7 +234,7 @@ function get_mortar_matrix_P(problem::Problem{Mortar2D}) C2 = sparse(problem.assembly.C2) D_ = C2[s,s] M_ = -C2[s,m] - P = ldltfact(1/2*(D_ + D_')) \ M_ + P = cholesky(1/2*(D_ + D_')) \ M_ return P end @@ -246,12 +246,12 @@ function FEMBase.eliminate_boundary_conditions!(problem::Problem{Mortar2D}, K, M P = get_mortar_matrix_P(problem) ndim = size(K, 1) Id = ones(ndim) - Id[s] = 0.0 - Q = spdiagm(Id) - Q[m,s] += P' - K[:,:] = Q*K*Q' - M[:,:] = Q*M*Q' - f[:] = Q*f + Id[s] .= 0.0 + Q = sparse(Diagonal(Id)) + Q[m,s] .+= P' + K[:,:] .= Q*K*Q' + M[:,:] .= Q*M*Q' + f[:] .= Q*f return true end diff --git a/test/REQUIRE b/test/REQUIRE new file mode 100644 index 0000000..b340ebe --- /dev/null +++ b/test/REQUIRE @@ -0,0 +1,2 @@ +Documenter +Literate diff --git a/test/runtests.jl b/test/runtests.jl index 40232ab..c5619e6 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -1,10 +1,11 @@ # This file is a part of JuliaFEM. # License is MIT: see https://github.com/JuliaFEM/MortarContact2D.jl/blob/master/LICENSE -using MortarContact2D -using Base.Test +using FEMBase, MortarContact2D, Test, SparseArrays, LinearAlgebra, Statistics -@testset "test MortarContact2D.jl" begin +include(joinpath("..", "docs", "make.jl")) + +@testset "MortarContact2D.jl" begin @testset "Projecting vertices between surfaces" begin include("test_mortar2d_calculate_projection.jl") @@ -27,3 +28,5 @@ using Base.Test end end + +include(joinpath("..", "docs", "deploy.jl")) diff --git a/test/test_contact2d.jl b/test/test_contact2d.jl index ca2b1ef..7ef72ff 100644 --- a/test/test_contact2d.jl +++ b/test/test_contact2d.jl @@ -1,9 +1,8 @@ # This file is a part of JuliaFEM. # License is MIT: see https://github.com/JuliaFEM/MortarContact2D.jl/blob/master/LICENSE -using MortarContact2D +using FEMBase, MortarContact2D, Test using MortarContact2D: get_slave_dofs, get_master_dofs -using Base.Test # Original problem: # @@ -83,6 +82,11 @@ update!(master, "geometry", ([1.0, 0.6], [0.0, 0.6])) update!(master, "displacement", 0.0 => (zeros(2), zeros(2))) update!(master, "displacement", 1.0 => ([-0.18125, -0.15], [-0.21875, -0.15])) problem = Problem(Contact2D, "LOWER_TO_UPPER", 2, "displacement") +push!(problem.properties.store_fields, "complementarity condition") +push!(problem.properties.store_fields, "active nodes") +push!(problem.properties.store_fields, "inactive nodes") +push!(problem.properties.store_fields, "stick nodes") +push!(problem.properties.store_fields, "slip nodes") add_slave_elements!(problem, [slave]) add_master_elements!(problem, [master]) um = master("displacement", (0.0,), 1.0) @@ -93,21 +97,21 @@ la = slave("lambda", (0.0,), 1.0) @test isapprox(la, [0.0, 30.375]) assemble!(problem, 0.0) -@test isempty(full(problem.assembly.K)) -@test isempty(full(problem.assembly.C1)) -@test isempty(full(problem.assembly.C2)) -@test isempty(full(problem.assembly.D)) -@test isempty(full(problem.assembly.f)) -@test isempty(full(problem.assembly.g)) +@test iszero(sparse(problem.assembly.K)) +@test iszero(sparse(problem.assembly.C1)) +@test iszero(sparse(problem.assembly.C2)) +@test iszero(sparse(problem.assembly.D)) +@test iszero(sparse(problem.assembly.f)) +@test iszero(sparse(problem.assembly.g)) empty!(problem.assembly) assemble!(problem, 1.0) -@test isempty(full(problem.assembly.K)) -@test isempty(full(problem.assembly.f)) -C1 = full(problem.assembly.C1, 4, 8) -C2 = full(problem.assembly.C2, 4, 8) -D = full(problem.assembly.D, 4, 8) -g = full(problem.assembly.g, 4, 1) +@test iszero(sparse(problem.assembly.K)) +@test iszero(sparse(problem.assembly.f)) +C1 = Matrix(sparse(problem.assembly.C1, 4, 8)) +C2 = Matrix(sparse(problem.assembly.C2, 4, 8)) +D = Matrix(sparse(problem.assembly.D, 4, 8)) +g = Vector(sparse(problem.assembly.g, 4, 1)[:]) C1_expected = [ 0.5 0.0 0.0 0.0 0.0 0.0 -0.5 0.0 0.0 0.5 0.0 0.0 0.0 0.0 0.0 -0.5 diff --git a/test/test_contact_segmentation.jl b/test/test_contact_segmentation.jl index d208950..c729ed2 100644 --- a/test/test_contact_segmentation.jl +++ b/test/test_contact_segmentation.jl @@ -1,11 +1,9 @@ # This file is a part of JuliaFEM. # License is MIT: see https://github.com/JuliaFEM/MortarContact2D.jl/blob/master/LICENSE -using MortarContact2D +using FEMBase, MortarContact2D, Test using MortarContact2D: create_contact_segmentation, calculate_normals -using FEMBase.Test - X = Dict(1 => [0.0, 0.0], 2 => [1.0, 0.0], 3 => [0.0, 1.0], diff --git a/test/test_mortar2d_calculate_projection.jl b/test/test_mortar2d_calculate_projection.jl index 3029604..5a16240 100644 --- a/test/test_mortar2d_calculate_projection.jl +++ b/test/test_mortar2d_calculate_projection.jl @@ -1,13 +1,11 @@ # This file is a part of JuliaFEM. # License is MIT: see https://github.com/JuliaFEM/MortarContact2D.jl/blob/master/LICENSE -using MortarContact2D +using FEMBase, MortarContact2D, Test using MortarContact2D: calculate_normals, project_from_master_to_slave, project_from_slave_to_master -using Base.Test - X = Dict( 1 => [0.0, 1.0], 2 => [5/4, 1.0], diff --git a/test/test_mortar2d_ex1.jl b/test/test_mortar2d_ex1.jl index f4bf7f1..1172336 100644 --- a/test/test_mortar2d_ex1.jl +++ b/test/test_mortar2d_ex1.jl @@ -3,11 +3,10 @@ # Simple coupling interface of four elements -using MortarContact2D +using FEMBase, MortarContact2D, Test using MortarContact2D: get_slave_dofs, get_master_dofs, get_mortar_matrix_D, get_mortar_matrix_M, get_mortar_matrix_P -using Base.Test # Slave side coordinates: Xs = Dict( @@ -54,6 +53,6 @@ P_expected = 1/320*[329 -24 15; 46 304 -30; -21 56 285] # Results: -@test isapprox(full(D), D_expected) -@test isapprox(full(M), M_expected) -@test isapprox(full(P), P_expected) +@test isapprox(D, D_expected) +@test isapprox(M, M_expected) +@test isapprox(P, P_expected) diff --git a/test/test_mortar2d_ex2.jl b/test/test_mortar2d_ex2.jl index 6be14d8..b2f3d16 100644 --- a/test/test_mortar2d_ex2.jl +++ b/test/test_mortar2d_ex2.jl @@ -4,8 +4,7 @@ # Construct discrete coupling operator between two elements and eliminate # boundary conditions from global matrix assembly. Thesis, p. 73. -using MortarContact2D -using Base.Test +using FEMBase, MortarContact2D, Test X = Dict( 1 => [1.0, 2.0],