From bf8e4cf03038a3d9052a85ad8a185d7d1d99a55f Mon Sep 17 00:00:00 2001 From: ovainola Date: Thu, 25 Jun 2015 19:13:28 +0300 Subject: [PATCH 01/10] Added type template to interpolate function and added test for interpolate --- src/elasticity_solver.jl | 4 ++-- test/test_elasticity_solver.jl | 1 + 2 files changed, 3 insertions(+), 2 deletions(-) diff --git a/src/elasticity_solver.jl b/src/elasticity_solver.jl index f9a429b..8d277a3 100644 --- a/src/elasticity_solver.jl +++ b/src/elasticity_solver.jl @@ -25,11 +25,11 @@ basis :: Function ip :: Array{Number, 1} Point to interpolate """ -> -function interpolate(field::Array{Float64,1}, basis::Function, ip) +function interpolate{T<:Real}(field::Array{T,1}, basis::Function, ip) result = dot(field, basis(ip)) return result end -function interpolate(field::Array{Float64,2}, basis::Function, ip) +function interpolate{T<:Real}(field::Array{T,2}, basis::Function, ip) m, n = size(field) bip = basis(ip) tmp = size(bip) diff --git a/test/test_elasticity_solver.jl b/test/test_elasticity_solver.jl index cc29d23..e51d306 100644 --- a/test/test_elasticity_solver.jl +++ b/test/test_elasticity_solver.jl @@ -154,6 +154,7 @@ facts("test interpolation of different field variables") do F3 = F2' F4 = [0.0 0.0; 10.0 0.0; 10.0 1.0; 0.0 1.0]' F5 = F4' + F6 = [36, 36, 36, 36] @fact interpolate(F1, N, [0.0, 0.0]) => 36.0 @fact interpolate(F2, N, [0.0, 0.0]) => 36.0 From 6249626227b0afaad3e774b3bd1402146aa4ca8e Mon Sep 17 00:00:00 2001 From: ovainola Date: Thu, 25 Jun 2015 19:20:15 +0300 Subject: [PATCH 02/10] Added test to interpolate function --- test/test_elasticity_solver.jl | 1 + 1 file changed, 1 insertion(+) diff --git a/test/test_elasticity_solver.jl b/test/test_elasticity_solver.jl index e51d306..202f67d 100644 --- a/test/test_elasticity_solver.jl +++ b/test/test_elasticity_solver.jl @@ -162,6 +162,7 @@ facts("test interpolation of different field variables") do @fact interpolate(F4, N, [0.0, 0.0]) => [5.0; 0.5] @fact interpolate(F5, N, [0.0, 0.0]) => [5.0; 0.5] @fact interpolate(F5, dNdξ, [0.0, 0.0]) => [5.0 0.0; 0.0 0.5] + @fact interpolate(F6, N, [0.0, 0.0]) => 36 end From 86745ada8303e6f43bb5783b3534d60da4f0c766 Mon Sep 17 00:00:00 2001 From: ovainola Date: Thu, 25 Jun 2015 19:50:34 +0300 Subject: [PATCH 03/10] Adding license comment to files missing it and removing non-ascii characters from elasticity-solver --- src/elasticity_solver.jl | 48 ++++++++++++++++++---------------- src/interfaces.jl | 6 +++-- src/shape_functions.jl | 5 +++- src/xdmf.jl | 6 +++-- test/test_elasticity_solver.jl | 3 +++ test/test_model.jl | 3 +++ test/test_xdmf.jl | 3 +++ 7 files changed, 46 insertions(+), 28 deletions(-) diff --git a/src/elasticity_solver.jl b/src/elasticity_solver.jl index 8d277a3..ce67724 100644 --- a/src/elasticity_solver.jl +++ b/src/elasticity_solver.jl @@ -1,4 +1,6 @@ -# This file is a part of JuliaFEM. License is MIT: https://github.com/ovainola/JuliaFEM/blob/master/README.md +# This file is a part of JuliaFEM. +# License is MIT: see https://github.com/JuliaFEM/JuliaFEM.jl/blob/master/LICENSE.md + module elasticity_solver using Logging @@ -63,7 +65,7 @@ end @doc """ Calculate local tangent stiffness matrix and residual force vector R = T - F """ -> -function calc_local_matrices!(X, u, R, Kt, N, dNdξ, λ_, μ_, ipoints, iweights) +function calc_local_matrices!(X, u, R, Kt, N, dNdchi, lambda_, mu_, ipoints, iweights) dim, nnodes = size(X) I = eye(dim) R[:,:] = 0.0 @@ -73,34 +75,34 @@ function calc_local_matrices!(X, u, R, Kt, N, dNdξ, λ_, μ_, ipoints, iweights for m = 1:length(iweights) w = iweights[m] - ξ = ipoints[m, :] + chi = ipoints[m, :] # interpolate material parameters from element node fields - #λ = (λ_*N(ξ))[1] - #μ = (μ_*N(ξ))[1] - # Jᵀ = X*dNdξ(ξ) - #@debug("Jt:\n",Jᵀ) - λ = interpolate(λ_, N, ξ) - μ = interpolate(μ_, N, ξ) - Jᵀ = interpolate(X, dNdξ, ξ) - detJ = det(Jᵀ) - ∇N = inv(Jᵀ)*dNdξ(ξ)' - ∇u = u*∇N' - F = I + ∇u # Deformation gradient - E = 1/2*(∇u' + ∇u + ∇u'*∇u) # Green-Lagrange strain tensor - S = λ*trace(E)*I + 2*μ*E # PK2 stress tensor + #lambda = (lambda_*N(chi))[1] + #mu = (mu_*N(chi))[1] + # Jt = X*dNdchi(chi) + #@debug("Jt:\n",Jt) + lambda = interpolate(lambda_, N, chi) + mu = interpolate(mu_, N, chi) + Jt = interpolate(X, dNdchi, chi) + detJ = det(Jt) + deltaN = inv(Jt)*dNdchi(chi)' + delta_u = u*deltaN' + F = I + delta_u # Deformation gradient + E = 1/2*(delta_u' + delta_u + delta_u'*delta_u) # Green-Lagrange strain tensor + S = lambda*trace(E)*I + 2*mu*E # PK2 stress tensor P = F*S # PK1 stress tensor - R[:,:] += w*P*∇N*detJ + R[:,:] += w*P*deltaN*detJ for p = 1:nnodes for i = 1:dim dF[:,:] = 0.0 - dF[i,:] = ∇N[:,p] + dF[i,:] = deltaN[:,p] dE = 1/2*(F'*dF + dF'*F) - dS = λ*trace(dE)*I + 2*μ*dE + dS = lambda*trace(dE)*I + 2*mu*dE dP = dF*S + F*dS for q = 1:nnodes for j = 1:dim - Kt[dim*(p-1)+i,dim*(q-1)+j] += w*(dP[j,:]*∇N[:,q])[1]*detJ + Kt[dim*(p-1)+i,dim*(q-1)+j] += w*(dP[j,:]*deltaN[:,q])[1]*detJ end end end @@ -274,7 +276,7 @@ end Solve one increment of elasticity problem """ -> function solve_elasticity_increment!(X, u, du, elmap, nodalloads, - dirichletbc, λ, μ, N, dNdξ, ipoints, + dirichletbc, lambda, mu, N, dNdchi, ipoints, iweights) if length(size(elmap)) == 1 # quick hack for just one element @@ -297,8 +299,8 @@ function solve_elasticity_increment!(X, u, du, elmap, nodalloads, # this can be parallelized for i in 1:nelements eldofs = elmap[:,i] - calc_local_matrices!(X[:, eldofs], u[:, eldofs], R, Kt, N, dNdξ, - λ[eldofs], μ[eldofs], ipoints, iweights) + calc_local_matrices!(X[:, eldofs], u[:, eldofs], R, Kt, N, dNdchi, + lambda[eldofs], mu[eldofs], ipoints, iweights) assemble!(Kt, eldofs, Imat, Jmat, Vmat) assemble!(R, eldofs, Ivec, Vvec) end diff --git a/src/interfaces.jl b/src/interfaces.jl index 860ae35..642fa39 100644 --- a/src/interfaces.jl +++ b/src/interfaces.jl @@ -1,4 +1,6 @@ -# This file is a part of JuliaFEM. License is MIT: https://github.com/ovainola/JuliaFEM/blob/master/README.md +# This file is a part of JuliaFEM. +# License is MIT: see https://github.com/JuliaFEM/JuliaFEM.jl/blob/master/LICENSE.md + module interfaces using Logging @@ -21,4 +23,4 @@ function solve_elasticity_interface!() return 0 end -end \ No newline at end of file +end diff --git a/src/shape_functions.jl b/src/shape_functions.jl index f28250e..ef9b8f0 100644 --- a/src/shape_functions.jl +++ b/src/shape_functions.jl @@ -1,10 +1,13 @@ +# This file is a part of JuliaFEM. +# License is MIT: see https://github.com/JuliaFEM/JuliaFEM.jl/blob/master/LICENSE.md + function test(a) return a*2 end function C2D4(xx) - θ = 1 + theta = 1 chi = xx[1] eta = xx[2] return = [(1 - chi) * (1 - eta), diff --git a/src/xdmf.jl b/src/xdmf.jl index 9126b2e..f285c4c 100644 --- a/src/xdmf.jl +++ b/src/xdmf.jl @@ -1,4 +1,6 @@ -# This file is a part of JuliaFEM. License is MIT: https://github.com/ovainola/JuliaFEM/blob/master/README.md +# This file is a part of JuliaFEM. +# License is MIT: see https://github.com/JuliaFEM/JuliaFEM.jl/blob/master/LICENSE.md + module xdmf using Logging @@ -103,4 +105,4 @@ function xdmf_new_field(grid, name, source, data) end -end \ No newline at end of file +end diff --git a/test/test_elasticity_solver.jl b/test/test_elasticity_solver.jl index 202f67d..74ff918 100644 --- a/test/test_elasticity_solver.jl +++ b/test/test_elasticity_solver.jl @@ -1,3 +1,6 @@ +# This file is a part of JuliaFEM. +# License is MIT: see https://github.com/JuliaFEM/JuliaFEM.jl/blob/master/LICENSE.md + using FactCheck using Logging @Logging.configure(level=INFO) diff --git a/test/test_model.jl b/test/test_model.jl index a8f3acc..4052fc9 100644 --- a/test/test_model.jl +++ b/test/test_model.jl @@ -1,3 +1,6 @@ +# This file is a part of JuliaFEM. +# License is MIT: see https://github.com/JuliaFEM/JuliaFEM.jl/blob/master/LICENSE.md + using FactCheck #using JuliaFEM using Logging diff --git a/test/test_xdmf.jl b/test/test_xdmf.jl index 148d3ec..fd0284f 100644 --- a/test/test_xdmf.jl +++ b/test/test_xdmf.jl @@ -1,3 +1,6 @@ +# This file is a part of JuliaFEM. +# License is MIT: see https://github.com/JuliaFEM/JuliaFEM.jl/blob/master/LICENSE.md + using FactCheck using Logging @Logging.configure(level=INFO) From 85be00bc0c985f6dfd129d9d58d9fc7efb82481c Mon Sep 17 00:00:00 2001 From: ovainola Date: Thu, 25 Jun 2015 19:58:09 +0300 Subject: [PATCH 04/10] Added '->' symbol at the end of docstring --- src/elasticity_solver.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/elasticity_solver.jl b/src/elasticity_solver.jl index ce67724..9881cdd 100644 --- a/src/elasticity_solver.jl +++ b/src/elasticity_solver.jl @@ -214,7 +214,7 @@ Raises Exception, if displacement boundary conditions given, i.e. DX=2 for some node, for example. -""" +""" -> function eliminate_boundary_conditions(dirichletbc, I, J, V) if any(dirichletbc .> 0) throw("displacement boundary condition not supported") From 03eb2f24d0f459d9a5613a4e5bfb284a83cc8e9f Mon Sep 17 00:00:00 2001 From: ovainola Date: Thu, 25 Jun 2015 20:27:47 +0300 Subject: [PATCH 05/10] new update --- abaqus_reader.jl | 133 +++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 133 insertions(+) create mode 100644 abaqus_reader.jl diff --git a/abaqus_reader.jl b/abaqus_reader.jl new file mode 100644 index 0000000..d7463a0 --- /dev/null +++ b/abaqus_reader.jl @@ -0,0 +1,133 @@ +# This file is a part of JuliaFEM. +# License is MIT: https://github.com/ovainola/JuliaFEM/blob/master/README.md +module abaqus_reader + +using Logging +@Logging.configure(level=DEBUG) + +VERSION < v"0.4-" && using Docile + +eldims = Dict("C3D10" => 10) +global handlers = Dict() + +""" +Register new handler for parser +""" +function add_handler(section, function_name) + handlers[section] = function_name +end + +function create_or_get(model, key) + if !(key in keys(model)) + model[key] = Dict() + end + return model[key] +end + +function parse_header(header_line) + args = map(s -> strip(s), split(header_line, ",")) + args[1] = strip(args[1], '*') + d = Dict({"section" => args[1]}) + options = Dict() + for k in args[2:end] + args2 = split(k, "=") + options[args2[1]] = args2[2] + end + d["options"] = options + return d +end + +function parse_node_section(model, header, data) + nodes = create_or_get(model, "nodes") + for line in split(data, "\n") + m = matchall(r"[-0-9.]+", line) + id = parse(Int, m[1]) + coords = float(m[2:end]) + nodes[id] = coords + end +end + +function parse_element_section(model, header, data) + eltype = header["options"]["TYPE"] + if !(eltype in keys(eldims)) + throw("Element $eltype dimension information missing") + end + eldim = eldims[eltype] + m = matchall(r"[0-9]+", data) + m = map(integer, m) + elements = create_or_get(model, "elements") + m = reshape(m, eldim+1, round(Int, length(m)/(eldim+1))) + nel = size(m)[2] + Logging.debug("$nel elements found") + for i=1:nel + elements[m[1,i]] = m[2:end,i] + end + if "ELSET" in keys(header["options"]) + elsets = create_or_get(model, "elsets") + elset_name = header["options"]["ELSET"] + Logging.info("Creating ELSET $elset_name") + elsets[elset_name] = Int64[] + for i=1:nel + push!(elsets[elset_name], m[1,i]) + end + end +end + + +function parse_nodeset_section(model, header, data) + nset_name = header["options"]["NSET"] + Logging.debug("Creating node set $nset_name") + m = matchall(r"[0-9]+", data) + node_ids = map(integer, m) + nsets = create_or_get(model, "nsets") + nsets[nset_name] = Int64[] + for j in node_ids + push!(nsets[nset_name], j) + end +end + + +function parse_abaqus(fid) + model = Dict() + section = None + header = None + data = "" + Logging.info("Registered handlers: $(keys(handlers))") + + function process_section(section) + if section == None + return + end + if !(section in keys(handlers)) + Logging.info("Don't know what to do with data in section $section") + Logging.info("Skipping $(length(data)) bytes of unknown data") + return + end + handlers[section](model, header, strip(data)) + data = "" + end + + for line in eachline(fid) + if beginswith(line, "**") + continue + end + if beginswith(line, "*") + process_section(section) + header = parse_header(line) + Logging.debug("Found ", header["section"], " section") + section = header["section"] + continue + end + data *= line + end + process_section(section) + return model +end + +# add handlers +add_handler("NODE", parse_node_section) +add_handler("ELEMENT", parse_element_section) +add_handler("NSET", parse_nodeset_section) + + +end From 0d81caf119d0a2d53c37df1708871ae58ff6fe63 Mon Sep 17 00:00:00 2001 From: ovainola Date: Thu, 25 Jun 2015 20:28:57 +0300 Subject: [PATCH 06/10] new commit --- abaqus_reader.jl => src/abaqus_reader.jl | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) rename abaqus_reader.jl => src/abaqus_reader.jl (99%) diff --git a/abaqus_reader.jl b/src/abaqus_reader.jl similarity index 99% rename from abaqus_reader.jl rename to src/abaqus_reader.jl index d7463a0..de5557d 100644 --- a/abaqus_reader.jl +++ b/src/abaqus_reader.jl @@ -10,9 +10,9 @@ VERSION < v"0.4-" && using Docile eldims = Dict("C3D10" => 10) global handlers = Dict() -""" +@doc """ Register new handler for parser -""" +""" -> function add_handler(section, function_name) handlers[section] = function_name end From 54e577f823ca92fec2b9a83565ada26569962b03 Mon Sep 17 00:00:00 2001 From: ovainola Date: Thu, 25 Jun 2015 21:04:39 +0300 Subject: [PATCH 07/10] own fault in aba reader --- src/abaqus_reader.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/abaqus_reader.jl b/src/abaqus_reader.jl index 749b6a9..480bade 100644 --- a/src/abaqus_reader.jl +++ b/src/abaqus_reader.jl @@ -14,7 +14,7 @@ VERSION < v"0.4-" && using Docile eldims = Dict("C3D10" => 10) global handlers = Dict() -<<<<<<< HEAD + @doc """ Register new handler for parser """ -> From bfde673714d223dbcf2cc7d12204dc63058b10bb Mon Sep 17 00:00:00 2001 From: ovainola Date: Thu, 25 Jun 2015 21:08:42 +0300 Subject: [PATCH 08/10] own fault in aba reader ver.2 --- src/abaqus_reader.jl | 2 -- 1 file changed, 2 deletions(-) diff --git a/src/abaqus_reader.jl b/src/abaqus_reader.jl index 480bade..e02523a 100644 --- a/src/abaqus_reader.jl +++ b/src/abaqus_reader.jl @@ -1,9 +1,7 @@ -<<<<<<< HEAD # This file is a part of JuliaFEM. # License is MIT: https://github.com/ovainola/JuliaFEM/blob/master/README.md ======= # This file is a part of JuliaFEM. License is MIT: https://github.com/ovainola/JuliaFEM/blob/master/README.md ->>>>>>> 3e48ae32abe7a17af89cd1e23d798532aed2960a module abaqus_reader using Logging From 3971cc462d66b5147690b61f39e9c7d5688f77ae Mon Sep 17 00:00:00 2001 From: ovainola Date: Thu, 25 Jun 2015 21:17:29 +0300 Subject: [PATCH 09/10] meneeko ikina lapi --- src/abaqus_reader.jl | 18 ++++-------------- src/shape_functions.jl | 5 ----- test/test_abaqus_reader.jl | 5 ++++- 3 files changed, 8 insertions(+), 20 deletions(-) diff --git a/src/abaqus_reader.jl b/src/abaqus_reader.jl index e02523a..6d89aa3 100644 --- a/src/abaqus_reader.jl +++ b/src/abaqus_reader.jl @@ -1,7 +1,6 @@ -# This file is a part of JuliaFEM. -# License is MIT: https://github.com/ovainola/JuliaFEM/blob/master/README.md -======= -# This file is a part of JuliaFEM. License is MIT: https://github.com/ovainola/JuliaFEM/blob/master/README.md +# This file is a part of JuliaFEM. +# License is MIT: see https://github.com/JuliaFEM/JuliaFEM.jl/blob/master/LICENSE.md + module abaqus_reader using Logging @@ -16,11 +15,6 @@ global handlers = Dict() @doc """ Register new handler for parser """ -> -======= -""" -Register new handler for parser -""" ->>>>>>> 3e48ae32abe7a17af89cd1e23d798532aed2960a function add_handler(section, function_name) handlers[section] = function_name end @@ -137,9 +131,5 @@ add_handler("NODE", parse_node_section) add_handler("ELEMENT", parse_element_section) add_handler("NSET", parse_nodeset_section) +end -<<<<<<< HEAD -end -======= -end ->>>>>>> 3e48ae32abe7a17af89cd1e23d798532aed2960a diff --git a/src/shape_functions.jl b/src/shape_functions.jl index ef9b8f0..968f7f3 100644 --- a/src/shape_functions.jl +++ b/src/shape_functions.jl @@ -1,11 +1,6 @@ # This file is a part of JuliaFEM. # License is MIT: see https://github.com/JuliaFEM/JuliaFEM.jl/blob/master/LICENSE.md -function test(a) - return a*2 -end - - function C2D4(xx) theta = 1 chi = xx[1] diff --git a/test/test_abaqus_reader.jl b/test/test_abaqus_reader.jl index cf0f91e..4274e13 100644 --- a/test/test_abaqus_reader.jl +++ b/test/test_abaqus_reader.jl @@ -1,3 +1,6 @@ +# This file is a part of JuliaFEM. +# License is MIT: see https://github.com/JuliaFEM/JuliaFEM.jl/blob/master/LICENSE.md + using FactCheck using Logging @Logging.configure(level=INFO) @@ -14,4 +17,4 @@ facts("test import abaqus model") do @fact length(model["nsets"]["SUPPORT"]) => 9 @fact length(model["nsets"]["LOAD"]) => 9 @fact length(model["nsets"]["TOP"]) => 83 -end \ No newline at end of file +end From b387070c8f59e7c800832eea2caac0d3c3362a14 Mon Sep 17 00:00:00 2001 From: ovainola Date: Thu, 25 Jun 2015 21:48:11 +0300 Subject: [PATCH 10/10] more fixes --- src/abaqus_reader.jl | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/src/abaqus_reader.jl b/src/abaqus_reader.jl index 6d89aa3..ce292d9 100644 --- a/src/abaqus_reader.jl +++ b/src/abaqus_reader.jl @@ -8,7 +8,7 @@ using Logging VERSION < v"0.4-" && using Docile -eldims = Dict("C3D10" => 10) +eldims = Dict({"C3D10" => 10}) global handlers = Dict()