diff --git a/.travis.yml b/.travis.yml index 49f5582..9af7731 100644 --- a/.travis.yml +++ b/.travis.yml @@ -4,8 +4,8 @@ os: - linux julia: - - 0.7 - - 1.0 + - 1.2 + - 1.3 - nightly addons: @@ -17,5 +17,14 @@ matrix: allow_failures: - julia: nightly +jobs: + include: + - stage: "Documentation" + os: linux + script: + - julia --project=docs/ -e 'using Pkg; Pkg.develop(PackageSpec(path=pwd())); Pkg.instantiate()' + - julia --project=docs/ docs/make.jl + after_success: skip + after_success: - julia -e 'cd(Pkg.dir("JuliaFEM")); Pkg.add("Coverage"); using Coverage; Coveralls.submit(Coveralls.process_folder())' diff --git a/Project.toml b/Project.toml index d010ff2..1b9a7f6 100644 --- a/Project.toml +++ b/Project.toml @@ -1,6 +1,7 @@ name = "JuliaFEM" uuid = "f80590ac-b429-510a-8a99-e7c46989f22d" -version = "0.5.0" +author = ["Jukka Aho ", "Tero Frondelius ", "Olli Väinölä"] +version = "0.5.1" [deps] AbaqusReader = "bc6b9049-e460-56d6-94b4-a597b2c0390d" @@ -19,6 +20,7 @@ LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e" MortarContact2D = "048d6160-1a0b-53cd-a5b3-316946cc8d80" MortarContact2DAD = "c1673bdb-6aff-560b-99da-c78ea6da9af3" Parameters = "d96e819e-fc66-5662-9728-84c9c7592b0a" +REPL = "3fa0cd96-eef1-5676-8a61-b3b8758bbffb" Reexport = "189a3867-3050-52da-a836-e630ba90ab69" SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf" Statistics = "10745b16-79ce-11e8-11f9-7d13ad32a3b2" @@ -29,5 +31,16 @@ Documenter = "e30172f5-a6a5-5a46-863b-614d45cd2de4" Pkg = "44cfe95a-1eb2-52ea-b672-e2afdf69b78f" Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40" +[compat] +HDF5 = "≥ 0.7.0" +LightXML = "≥ 0.4.0" +julia = "≥ 1.0" + +[extras] +LinearAlgebra = "37e2e46d-f89d-539d-b4ee-838fcccc9c8e" +SparseArrays = "2f01184e-e22b-5df5-ae63-d93ebab69eaf" +Statistics = "10745b16-79ce-11e8-11f9-7d13ad32a3b2" +Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40" + [targets] -test = ["Test", "Pkg", "Documenter"] +test = ["Test"] diff --git a/README.md b/README.md index 7df45f4..979fd8a 100644 --- a/README.md +++ b/README.md @@ -95,7 +95,7 @@ Errors should never pass silently. If you like using our package, please consider citing our [article](https://rakenteidenmekaniikka.journal.fi/article/view/64224/26397) ``` @article{frondelius2017juliafem, - title={JuliaFEM - open source solver for both industrial and academia usage}, + title={Julia{FEM} - open source solver for both industrial and academia usage}, volume={50}, url={https://rakenteidenmekaniikka.journal.fi/article/view/64224}, DOI={10.23998/rm.64224}, diff --git a/REQUIRE b/REQUIRE deleted file mode 100644 index 1c5b39d..0000000 --- a/REQUIRE +++ /dev/null @@ -1,17 +0,0 @@ -julia 0.7 -FEMBase -FEMBasis -FEMQuad -ForwardDiff -LightXML 0.4 -HDF5 0.7 -TimerOutputs -Reexport -Arpack -AbaqusReader -AsterReader -HeatTransfer -MortarContact2D -MortarContact2DAD -FEMBeam -Parameters diff --git a/docs/Project.toml b/docs/Project.toml new file mode 100644 index 0000000..cf645b5 --- /dev/null +++ b/docs/Project.toml @@ -0,0 +1,3 @@ +[deps] +Documenter = "e30172f5-a6a5-5a46-863b-614d45cd2de4" +Literate = "98b081ad-f1c9-55d3-8b20-4c87d4299306" diff --git a/docs/deploy.jl b/docs/deploy.jl index 2d5c1d6..df9a90d 100644 --- a/docs/deploy.jl +++ b/docs/deploy.jl @@ -5,7 +5,6 @@ using Documenter deploydocs( repo = "github.com/JuliaFEM/JuliaFEM.jl.git", - julia = "1.0", target = "build", deps = nothing, make = nothing) diff --git a/docs/make.jl b/docs/make.jl index 1ac2a68..79abe05 100644 --- a/docs/make.jl +++ b/docs/make.jl @@ -254,8 +254,9 @@ PAGES = [ @info("Pages in documentation", PAGES) makedocs(modules=[JuliaFEM], - format = :html, + format = Documenter.HTML(analytics="UA-83590644-1"), checkdocs = :all, - sitename = "JuliaFEM", - analytics = "UA-83590644-1", + sitename = "JuliaFEM.jl", pages = PAGES) + +include("deploy.jl") diff --git a/examples/2d_hertz_contact.jl b/examples/2d_hertz_contact.jl index 8ad1ad8..62daf37 100644 --- a/examples/2d_hertz_contact.jl +++ b/examples/2d_hertz_contact.jl @@ -23,7 +23,7 @@ # Substituting values, one gets accurate solution to be ``p_0 = 3585 \;\mathrm{MPa}`` and # ``a = 6.21 \;\mathrm{mm}``. -using JuliaFEM +using JuliaFEM, LinearAlgebra # Simulation starts by reading the mesh. Model is constructed and meshed using # SALOME, thus mesh format is .med. Mesh type is quite simple structure, @@ -33,7 +33,7 @@ using JuliaFEM # need to use `Mesh` in simulation anyway if we figure some other way to define # the geometry for elements. -datadir = Pkg.dir("JuliaFEM", "examples", "2d_hertz_contact") +datadir = abspath(joinpath(pathof(JuliaFEM), "..", "..", "examples", "2d_hertz_contact")) meshfile = joinpath(datadir, "hertz_2d_full.med") mesh = aster_read_mesh(meshfile) for (elset_name, element_ids) in mesh.element_sets @@ -153,6 +153,7 @@ Rt = 0.0 time = 0.0 for sel in contact_slave_elements for ip in get_integration_points(sel) + global Rn, Rt w = ip.weight*sel(ip, time, Val{:detJ}) n = sel("normal", ip, time) t = sel("tangent", ip, time) @@ -181,10 +182,11 @@ p0_acc = 3585.0 for (nid, n) in normal lan = dot(n, lambda[nid]) println("$nid => $lan") + global p0 p0 = max(p0, lan) end -p0 = round(p0, 2) -rtol = round(norm(p0-p0_acc)/max(p0,p0_acc)*100, 2) +p0 = round(p0, digits=2) +rtol = round(norm(p0-p0_acc)/max(p0,p0_acc)*100, digits=2) println("Maximum contact pressure p0 = $p0, p0_acc = $p0_acc, rtol = $rtol %") # To get rough approximation where does the contact open, we can find the element @@ -202,6 +204,7 @@ for element in contact_slave_elements println("Contact opening element lambda: la1 = $la1, la2 = $la2") x11, y11 = X1 x12, y12 = X2 + global a_rad a_rad = 1/2*abs(x11+x12) break end diff --git a/src/JuliaFEM.jl b/src/JuliaFEM.jl index 0254349..82c9676 100644 --- a/src/JuliaFEM.jl +++ b/src/JuliaFEM.jl @@ -119,7 +119,7 @@ using TimerOutputs export @timeit, print_timer import Base: getindex, setindex!, convert, length, size, isapprox, - similar, start, first, next, done, last, endof, vec, + similar, first, last, vec, ==, +, -, *, /, haskey, copy, push!, isempty, empty!, append!, read, copy diff --git a/src/io.jl b/src/io.jl index 1316a41..cc519e4 100644 --- a/src/io.jl +++ b/src/io.jl @@ -328,6 +328,29 @@ global const xdmf_element_mapping = Dict( "Wedge15" => "Wedge_15", "Hex20" => "Hex_20") +get_xdmf_element_code(::Element{Poi1}) = 1 +get_xdmf_element_code(::Element{Seg2}) = 2 +# get_xdmf_element_code(::Element{Polygon}) = 3 +get_xdmf_element_code(::Element{Tri3}) = 4 +get_xdmf_element_code(::Element{Quad4}) = 5 +get_xdmf_element_code(::Element{Tet4}) = 6 +get_xdmf_element_code(::Element{Pyr5}) = 7 +get_xdmf_element_code(::Element{Wedge6}) = 8 +get_xdmf_element_code(::Element{Hex8}) = 9 +# get_xdmf_element_code(::Element{Polyhedron}) = 16 + +get_xdmf_element_code(::Element{Seg3}) = 34 +get_xdmf_element_code(::Element{Quad9}) = 35 +get_xdmf_element_code(::Element{Tri6}) = 36 +get_xdmf_element_code(::Element{Quad8}) = 37 +get_xdmf_element_code(::Element{Tet10}) = 38 +# get_xdmf_element_code(::Element{Pyr13}) = 39 +get_xdmf_element_code(::Element{Wedge15}) = 40 +# get_xdmf_element_code(::Element{Wedge18}) = 41 +get_xdmf_element_code(::Element{Hex20}) = 48 +# get_xdmf_element_code(::Element{Hex24}) = 49 +get_xdmf_element_code(::Element{Hex27}) = 50 + """ get_spatial_collection() @@ -416,27 +439,48 @@ function update_xdmf!(xdmf::Xdmf, problem::Problem, time::Float64, fields::Vecto add_child(geometry, X_dataitem) # 5. save topology - all_elements = get_elements(problem) - nelements = length(all_elements) - element_types = unique(map(get_element_type, all_elements)) - nelement_types = length(element_types) - @debug("Xdmf: Saving topology of $nelements elements total, $nelement_types different element types.") - - for element_type in element_types - elements = collect(filter_by_element_type(element_type, all_elements)) - nelements = length(elements) - @debug("Xdmf: $nelements elements of type $element_type") - sort!(elements, by=get_element_id) - element_ids = map(get_element_id, elements) - element_conn = map(element -> [node_mapping[j]-1 for j in get_connectivity(element)], elements) - element_conn = hcat(element_conn...) - element_code = split(string(element_type), ".")[end] + mesh_type = "unstructured" + if mesh_type == "unstructured" + element_conn = Int64[] + for element in get_elements(problem) + xdmf_element_code = get_xdmf_element_code(element) + xdmf_element_code > 0 || continue + push!(element_conn, xdmf_element_code) + if xdmf_element_code == 2 + push!(element_conn, length(element)) + end + for j in get_connectivity(element) + push!(element_conn, node_mapping[j]-1) + end + end topology_dataitem = new_dataitem(xdmf, element_conn) - topology = new_child(frame, "Topology") - set_attribute(topology, "TopologyType", xdmf_element_mapping[element_code]) - set_attribute(topology, "NumberOfElements", length(elements)) + set_attribute(topology, "TopologyType", "Mixed") add_child(topology, topology_dataitem) + else + all_elements = get_elements(problem) + nelements = length(all_elements) + element_types = unique(map(get_element_type, all_elements)) + nelement_types = length(element_types) + @debug("Xdmf: Saving topology of $nelements elements total, $nelement_types different element types.") + if nelement_types != 1 + error("Xdmf: only single type of element supported by structured grid type!") + end + for element_type in element_types + elements = collect(filter_by_element_type(element_type, all_elements)) + nelements = length(elements) + @debug("Xdmf: $nelements elements of type $element_type") + sort!(elements, by=get_element_id) + element_ids = map(get_element_id, elements) + element_conn = map(element -> [node_mapping[j]-1 for j in get_connectivity(element)], elements) + element_conn = hcat(element_conn...) + element_code = split(string(element_type), ".")[end] + topology_dataitem = new_dataitem(xdmf, element_conn) + topology = new_child(frame, "Topology") + set_attribute(topology, "TopologyType", xdmf_element_mapping[element_code]) + set_attribute(topology, "NumberOfElements", length(elements)) + add_child(topology, topology_dataitem) + end end # 6. save requested fields @@ -444,7 +488,11 @@ function update_xdmf!(xdmf::Xdmf, problem::Problem, time::Float64, fields::Vecto field_dict = problem(field_name, time) field_center = "Node" field_node_ids = sort(collect(keys(field_dict))) - @assert node_ids == field_node_ids + if node_ids != field_node_ids + @error("geom node ids = $node_ids") + @error("field node ids = $field_node_ids") + error("!=, geometry does not match with field.") + end field_dim = length(field_dict[first(field_node_ids)]) if field_dim == 2 @debug("Xdmf: Field dimension = 2, extending to 3") diff --git a/src/problems_elasticity.jl b/src/problems_elasticity.jl index d0aab28..b63154b 100644 --- a/src/problems_elasticity.jl +++ b/src/problems_elasticity.jl @@ -154,9 +154,39 @@ function allocate_buffer(problem::Problem{Elasticity}, ::Vector{Element{El}}) wh nnodes = length(El) ndofs = dim*nnodes +<<<<<<< HEAD return Elasticity3DLocalBuffers(ndofs=ndofs, dim=dim, bi = BasisInfo(El)) end +======= + for element in elements + + u = element("displacement", time) + X = element("geometry", time) + + fill!(Km, 0.0) + fill!(Kg, 0.0) + fill!(f_int, 0.0) + fill!(f_ext, 0.0) + + for ip in get_integration_points(element) + eval_basis!(bi, X, ip) + w = ip.weight*bi.detJ + N = bi.N + dN = bi.grad # deriatives of basis functions w.r.t. X, i.e. ∂N/∂X + grad!(bi, gradu, u) # displacement gradient ∇u + + # calculate strain tensor and deformation gradient + fill!(strain, 0.0) + fill!(F, 0.0) + F[:,:] += I + if props.finite_strain + strain[:,:] = 1/2 * (gradu + gradu' + gradu'*gradu) + F[:,:] += gradu + else + strain[:,:] = 1/2 * (gradu + gradu') + end +>>>>>>> master function reset_element!(buf::Elasticity3DLocalBuffers) fill!(buf.Km, 0.0) @@ -175,6 +205,7 @@ function reset_integration_point!(buf::Elasticity3DLocalBuffers) return end +<<<<<<< HEAD function to_voigt!(strain_vec, strain) strain_vec[1] = strain[1,1] strain_vec[2] = strain[2,2] @@ -184,6 +215,16 @@ function to_voigt!(strain_vec, strain) strain_vec[6] = 2.0*strain[1,3] return end +======= + fill!(D, 0.0) + E = element("youngs modulus", ip, time)::Float64 + nu = element("poissons ratio", ip, time)::Float64 + la = E*nu/((1.0+nu)*(1.0-2.0*nu)) + mu = E/(2.0*(1.0+nu)) + D[1,1] = D[2,2] = D[3,3] = 2*mu + la + D[4,4] = D[5,5] = D[6,6] = mu + D[1,2] = D[2,1] = D[2,3] = D[3,2] = D[1,3] = D[3,1] = la +>>>>>>> master const u = ([0.0, 0.0, 0.0], [0.0, 0.0, 0.0], [0.0, 0.0, 0.0], [0.0, 0.0, 0.0], [0.0, 0.0, 0.0], [0.0, 0.0, 0.0], [0.0, 0.0, 0.0], [0.0, 0.0, 0.0], [0.0, 0.0, 0.0], [0.0, 0.0, 0.0]) const X = ([-93.7197, -93.7197, 150.883], [-91.657, -85.8251, 157.885], [-100.523, -88.8309, 157.883], [-91.6593, -88.8309, 157.883], [-92.6883, -89.7724, 154.384], [-96.0902, -87.328, 157.883], [-97.1216, -91.2753, 154.383], [-92.6895, -91.2753, 154.383], [-91.6581, -87.328, 157.883], [-96.0914, -88.8309, 157.883]) diff --git a/src/solvers_modal.jl b/src/solvers_modal.jl index 4ac9038..ca30cc4 100644 --- a/src/solvers_modal.jl +++ b/src/solvers_modal.jl @@ -289,7 +289,7 @@ function update_xdmf!(solver::Solver{Modal}) # save modes temporal_collection = get_temporal_collection(xdmf) - unknown_field_name = ucfirst(get_unknown_field_name(solver)) + unknown_field_name = uppercasefirst(get_unknown_field_name(solver)) frames = [] @timeit "save modes" for (j, eigval) in enumerate(real(solver.properties.eigvals)) diff --git a/test/REQUIRE b/test/REQUIRE deleted file mode 100644 index b340ebe..0000000 --- a/test/REQUIRE +++ /dev/null @@ -1,2 +0,0 @@ -Documenter -Literate diff --git a/test/runtests.jl b/test/runtests.jl index afba444..01a490d 100644 --- a/test/runtests.jl +++ b/test/runtests.jl @@ -3,8 +3,11 @@ using JuliaFEM, Test +<<<<<<< HEAD # include("../docs/make.jl") +======= +>>>>>>> master @testset "JuliaFEM.jl" begin @testset "test_dirichlet.jl" begin include("test_dirichlet.jl") @@ -163,5 +166,3 @@ using JuliaFEM, Test include("test_von_mises_material.jl") end end - -include("../docs/deploy.jl") diff --git a/test/test_elasticity_tet10_stiffness_matrix.jl b/test/test_elasticity_tet10_stiffness_matrix.jl index a60575f..b9b7ee1 100644 --- a/test/test_elasticity_tet10_stiffness_matrix.jl +++ b/test/test_elasticity_tet10_stiffness_matrix.jl @@ -34,4 +34,4 @@ eigs_expected = [8809.45, 4936.01, 2880.56, 2491.66, 2004.85, 195.832, 104.008, 72.7562, 64.4376, 53.8515, 23.8417, 16.6354, 9.54682, 6.93361, 2.22099, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0] -@test isapprox(eigs, eigs_expected; atol=1.0e-2) +@test isapprox(sort(eigs), sort(eigs_expected); atol=1.0e-2)