Glassvision/JuliaFEM.jl
1
0
Fork 0
mirror of https://github.com/JuliaFEM/JuliaFEM.jl.git synced 2026-09-20 01:59:59 +00:00
Code Issues Packages Projects Releases Wiki Activity

...

Browse Source
This commit is contained in:
Jukka Aho
2016-01-01 17:57:24 +02:00
parent ac45c31ce3
commit 407b0c4849
1 changed files with 52 additions and 15 deletions
Download Patch File Download Diff File Expand all files Collapse all files
test/test_mortar_2d_contact.jl
+52 -15
View File
@@ -5,7 +5,29 @@ using JuliaFEM.Test
using JuliaFEM.Core: Node, Seg2, Tri3, update!, calculate_normal_tangential_coordinates!,
PlaneStressLinearElasticityProblem, DirichletProblem, MortarProblem,
get_elements, DirectSolver
get_elements, DirectSolver, calculate_nodal_vector, FieldAssembly,
FieldProblem, set_linear_system_solver!, set_nonlinear_max_iterations!
import JuliaFEM.Core: postprocess_assembly!, linear_system_solver_preprocess!
function postprocess_assembly!(assembly::FieldAssembly, problem::FieldProblem{PlaneStressLinearElasticityProblem}, time::Real)
f = full(assembly.force_vector)
info("force vector = $(f')")
end
function linear_system_solver_preprocess!(solver, iter, time, K, f, C1, C2, D, g, sol, la, ::Type{Val{:foo}})
info("stiffness matrix")
dump(round(full(K), 3))
info("constraint matrix C1")
dump(round(full(C1), 3))
info("constraint matrix C2")
dump(round(full(C2), 3))
info("force vector")
dump(round(full(f)', 3))
info("constraint vector")
dump(round(full(g)', 3))
end
macro debug(msg)
haskey(ENV, "DEBUG") || return
@@ -14,24 +36,28 @@ end
@testset "2d frictionless contact" begin
nodes = Node[
[3.0, 3.0],
[4.0, 6.0],
[6.0, 6.0],
[7.0, 8.0],
[0.0, 0.0],
[3.0, 0.0],
[4.0, 0.0],
[10.0, 0.0]]
fel1 = Tri3([3, 4, 1])
fel2 = Tri3([5, 6, 2])
[6.0, 0.0],
[7.0, 0.0],
[19.0, 0.0]]
fel1 = Tri3([1, 3, 4])
fel2 = Tri3([2, 5, 6])
force = Seg2([3, 1])
bnd1 = Seg2([3, 4])
bnd2 = Seg2([5, 6])
sel = Seg2([1, 4])
mel = Seg2([5, 2])
update!([fel1, fel2, sel, mel, bnd1, bnd2], "geometry", nodes)
update!([fel1, fel2, force, sel, mel, bnd1, bnd2], "geometry", nodes)
prob = PlaneStressLinearElasticityProblem()
push!(prob, fel1, fel2)
update!(get_elements(prob), "youngs modulus", 90.0)
update!(get_elements(prob), "poissons ratio", 0.25)
push!(prob, force)
update!([fel1, fel2], "youngs modulus", 90.0)
update!([fel1, fel2], "poissons ratio", 0.25)
update!([force], "displacement traction force 1", 6/sqrt(2))
fel1["displacement nodal load"] = Vector{Float64}[[0.0, 0.0], [0.0, 0.0], [18.0, 0.0]]
bc = DirichletProblem("displacement", 2)
push!(bc, bnd1, bnd2)
@@ -50,9 +76,20 @@ end
solver = DirectSolver()
push!(solver, prob)
push!(solver, bc)
push!(solver, cont)
solver.method = :UMFPACK
solver.max_iterations = 5
call(solver, 0.0)
# push!(solver, cont)
set_linear_system_solver!(solver, :UMFPACK)
set_nonlinear_max_iterations!(solver, 2)
push!(solver.linear_system_solver_preprocessors, :foo)
time = 0.0
call(solver, time)
@debug begin
u = calculate_nodal_vector("displacement", 2, get_elements(prob), time)
info("solution vector")
dump(reshape(round(u, 8), 2, 6))
la = calculate_nodal_vector("reaction force", 2, get_elements(prob), time)
info("reaction force")
dump(reshape(round(la, 8), 2, 6))
end
end