diff --git a/geometry/3d_splitted_wedge/splitted_wedge.inp b/geometry/3d_splitted_wedge/splitted_wedge.inp new file mode 100644 index 0000000..368c58b --- /dev/null +++ b/geometry/3d_splitted_wedge/splitted_wedge.inp @@ -0,0 +1,117 @@ +**NSET COUNT = 48 +*NODE +1, 0.00000, 100.00000, 0.00000 +2, 100.00000, 0.00000, 0.00000 +3, 0.00000, 0.00000, 0.00000 +4, 50.00000, 50.00000, 0.00000 +5, 50.00000, 0.00000, 0.00000 +6, -0.00000, 50.00000, 0.00000 +7, 100.00000, 0.00000, 25.00000 +8, 0.00000, 100.00000, 25.00000 +9, 0.00000, 0.00000, 25.00000 +10, 50.00000, 50.00000, 25.00000 +11, -0.00000, 50.00000, 25.00000 +12, 50.00000, 0.00000, 25.00000 +13, 0.00000, 0.00000, 12.50000 +14, 50.00000, 0.00000, 12.50000 +15, 100.00000, 0.00000, 12.50000 +16, 0.00000, 100.00000, 12.50000 +17, -0.00000, 50.00000, 12.50000 +18, 50.00000, 50.00000, 12.50000 +19, 25.00000, 75.00000, 18.75000 +20, 25.00000, 75.00000, 6.25000 +21, 75.00000, 25.00000, 18.75000 +22, 75.00000, 25.00000, 6.25000 +23, 25.00000, 25.00000, 6.25000 +24, 25.00000, 25.00000, 18.75000 +25, 25.00000, 25.00000, 31.25000 +26, 25.00000, 25.00000, 43.75000 +27, 0.00000, 100.00000, 25.00000 +28, 100.00000, 0.00000, 25.00000 +29, 0.00000, 0.00000, 25.00000 +30, 50.00000, 50.00000, 25.00000 +31, 50.00000, 0.00000, 25.00000 +32, -0.00000, 50.00000, 25.00000 +33, 100.00000, 0.00000, 50.00000 +34, 0.00000, 100.00000, 50.00000 +35, 0.00000, 0.00000, 50.00000 +36, 50.00000, 50.00000, 50.00000 +37, -0.00000, 50.00000, 50.00000 +38, 50.00000, 0.00000, 50.00000 +39, 0.00000, 0.00000, 37.50000 +40, 50.00000, 0.00000, 37.50000 +41, 100.00000, 0.00000, 37.50000 +42, 0.00000, 100.00000, 37.50000 +43, -0.00000, 50.00000, 37.50000 +44, 50.00000, 50.00000, 37.50000 +45, 25.00000, 75.00000, 43.75000 +46, 25.00000, 75.00000, 31.25000 +47, 75.00000, 25.00000, 43.75000 +48, 75.00000, 25.00000, 31.25000 +** +**ELSET COUNT = 6 +**HWCOLOR COMP 98 0 +*ELEMENT, TYPE=C3D10, ELSET=LOWER + 26, 1, 3, 2, 18, 6, 5, 4, + 20, 23, 22 + 27, 1, 8, 9, 18, 16, 11, 17, + 20, 19, 24 + 28, 3, 9, 7, 18, 13, 12, 14, + 23, 24, 21 + 29, 7, 9, 8, 18, 12, 11, 10, + 21, 24, 19 + 30, 1, 9, 3, 18, 17, 13, 6, + 20, 24, 23 + 31, 2, 3, 7, 18, 5, 14, 15, + 22, 23, 21 +** +**ELSET COUNT = 6 +**HWCOLOR COMP 99 0 +*ELEMENT, TYPE=C3D10, ELSET=UPPER + 32, 27, 34, 35, 44, 42, 37, 43, + 46, 45, 26 + 33, 29, 35, 33, 44, 39, 38, 40, + 25, 26, 47 + 34, 33, 35, 34, 44, 38, 37, 36, + 47, 26, 45 + 35, 27, 35, 29, 44, 43, 39, 32, + 46, 26, 25 + 36, 28, 29, 33, 44, 31, 40, 41, + 48, 25, 47 + 37, 27, 29, 28, 44, 32, 31, 30, + 46, 25, 48 +** +**ELSET COUNT = 1 +*SURFACE,TYPE=ELEMENT,NAME=UPPER_FIXED +34,S1 +** +**ELSET COUNT = 1 +*SURFACE,TYPE=ELEMENT,NAME=LOWER_FIXED +26,S1 +** +**ELSET COUNT = 1 +*SURFACE,TYPE=ELEMENT,NAME=UPPER_TO_LOWER +37,S1 +** +**ELSET COUNT = 1 +*SURFACE,TYPE=ELEMENT,NAME=LOWER_TO_UPPER +29,S1 +** +**Property Definitions +** +*SOLID SECTION, ELSET=LOWER, MATERIAL=Def_Material +*SOLID SECTION, ELSET=UPPER, MATERIAL=Def_Material +** +**Material Definitions +** +**Material:Def_Material +*MATERIAL,NAME=Def_Material +*ELASTIC,TYPE=ISO +2.08000e+005,3.00000e-001 +*DENSITY +7.80000e-009, +*SPECIFIC HEAT +5.00000e-001 +*CONDUCTIVITY +4.98100e-002 +** \ No newline at end of file diff --git a/geometry/3d_splitted_wedge/tie_contact.jl b/geometry/3d_splitted_wedge/tie_contact.jl new file mode 100644 index 0000000..34686b1 --- /dev/null +++ b/geometry/3d_splitted_wedge/tie_contact.jl @@ -0,0 +1,107 @@ +using JuliaFEM +using JuliaFEM.Preprocess +using JuliaFEM.Postprocess +using JuliaFEM.Abaqus: create_surface_elements + +Logging.configure(level=Logging.DEBUG) + +function create_body(mesh, name, E, nu, rho) + body = Problem(mesh, Elasticity, name, 3) + update!(body, "youngs modulus", E) + update!(body, "poissons ratio", nu) + update!(body, "density", rho) + return body +end + +""" Create boundary condition from surface set. """ +function create_bc_from_surface_set(mesh, name, u1, u2, u3) + bc = Problem(Dirichlet, string(name), 3, "displacement") + bc.elements = create_surface_elements(mesh, name) + update!(bc, "displacement 1", u1) + update!(bc, "displacement 2", u2) + update!(bc, "displacement 3", u3) + return bc +end + +""" Create boundary condition from node set. """ +function create_bc_from_node_set(mesh, name, u1, u2, u3) + bc = Problem(Dirichlet, string(name), 3, "displacement") + bc.elements = [Element(Poi1, [nid]) for nid in mesh.node_sets[name]] + update!(bc, "geometry", mesh.nodes) + update!(bc, "displacement 1", u1) + update!(bc, "displacement 2", u2) + update!(bc, "displacement 3", u3) + return bc +end + +function create_bc(mesh, name, u1=0.0, u2=0.0, u3=0.0) + name = Symbol(name) + if haskey(mesh.surface_sets, name) + return create_bc_from_surface_set(mesh, name, u1, u2, u3) + elseif haskey(mesh.node_sets, name) + return create_bc_from_node_set(mesh, name, u1, u2, u3) + else + error("Mesh does not contain node or surface set $name") + end +end + +function create_interface(mesh, slave_surface::String, master_surface::String; +dual_basis=true) + interface = Problem(Mortar, "interface between $slave_surface and $master_surface", 3, "displacement") + interface.properties.dual_basis = dual_basis + slave_elements = create_surface_elements(mesh, Symbol(slave_surface)) + master_elements = create_surface_elements(mesh, Symbol(master_surface)) + nslaves = length(slave_elements) + nmasters = length(master_elements) + info("$nslaves slaves, $nmasters masters") + update!(slave_elements, "master elements", master_elements) + interface.elements = [slave_elements; master_elements] + return interface +end + +function create_interface(mesh, slave::Problem, master::Problem) + slave_surface = slave.name * "_TO_" * master.name + master_surface = master.name * "_TO_" * slave.name + return create_interface(mesh, slave_surface, master_surface) +end + +""" Convert Mesh object from quadratic to linear. """ +function to_linear!(mesh) + mapping = Dict(:Tet10 => :Tet4, :Tri6 => :Tri3) + nnodes = Dict(:Tet4 => 4, :Tri3 => 3) + for elid in keys(mesh.elements) + eltype = mesh.element_types[elid] + if haskey(mapping, eltype) + mesh.element_types[elid] = mapping[eltype] + nnodes_new = nnodes[mapping[eltype]] + mesh.elements[elid] = mesh.elements[elid][1:nnodes_new] + end + end +end + +# start of simulation +mesh = abaqus_read_mesh("splitted_wedge.inp") +#to_linear!(mesh) +info("element sets = ", collect(keys(mesh.element_sets))) +info("surface sets = ", collect(keys(mesh.surface_sets))) + +# parts +lower = create_body(mesh, "LOWER", 210.0e3, 0.3, 7.85e-9) +upper = create_body(mesh, "UPPER", 210.0e3, 0.3, 7.85e-9) + +# boundary conditions +bc1 = create_bc(mesh, "LOWER_FIXED") +bc2 = create_bc(mesh, "UPPER_FIXED") + +# load +load = Problem(Elasticity, "pressure load", 3) +load.elements = create_surface_elements(mesh, :UPPER_FIXED) +update!(load, "surface pressure", 20.0) + +contact = create_interface(mesh, "UPPER_TO_LOWER", "LOWER_TO_UPPER") + +# solution +solver = Solver(Linear, lower, upper, bc1, load, contact) +solver.xdmf = Xdmf("results"; overwrite=true) +solver() +