new test geometry and test of quadratic tie contact using linear surface elements

This commit is contained in:
Jukka Aho
2016-12-15 09:19:00 +02:00
parent d24b1e178a
commit 57d9bd35b2
2 changed files with 224 additions and 0 deletions
@@ -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
**
+107
View File
@@ -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()