make backwards compatible

This commit is contained in:
Kristoffer Carlsson
2018-11-29 17:40:55 -05:00
parent ee7532a749
commit 09c355b916
6 changed files with 414 additions and 32 deletions
+21 -19
View File
@@ -53,8 +53,6 @@ function get_formulation_type(problem::Problem{Elasticity})
return :incremental
end
using InteractiveUtils
"""
assemble!(assembly:Assembly, problem::Problem{Elasticity}, elements, time)
@@ -76,26 +74,23 @@ function assemble!(assembly::Assembly, problem::Problem{Elasticity},
elements::Vector{T}, time, formulation) where {T <: Element}
if problem.assemble_parallel
@assert problem.assemble_csc
# Threaded assembly
assemblers = [FEMSparse.start_assemble(assembly.K, assembly.f) for i in 1:Threads.nthreads()]
assemblers = [FEMSparse.start_assemble(assembly.K_csc, assembly.f_csc) for i in 1:Threads.nthreads()]
local_buffers = [allocate_buffer(problem, elements) for i in 1:Threads.nthreads()]
#TODO: We have to be a bit careful here, the index of the element is no longer
#
# should only loop over elements that exist in `elements` here
for (color, elements) in FEMBase.get_color_ranges(elements)
Threads.@threads for i in 1:length(elements)
for i in 1:length(elements)
element = elements[i]
tid = Threads.threadid()
assemble_element!(assembly, assemblers[tid], problem, element, local_buffers[tid], time, formulation)
assemble_element!(assembly, assemblers[tid], problem, element, local_buffers[tid], time, formulation, true)
end
end
else
# Normal assembly
local_buffer = allocate_buffer(problem, elements)
assembler = FEMSparse.start_assemble(assembly.K, assembly.f)
assembler = FEMSparse.start_assemble(assembly.K_csc, assembly.f_csc)
for i in 1:length(elements)
assemble_element!(assembly, assembler, problem, elements[i], local_buffer, time, formulation)
assemble_element!(assembly, assembler, problem, elements[i], local_buffer, time, formulation, problem.assemble_csc)
end
end
end
@@ -193,7 +188,8 @@ function assemble_element!(assembly::Assembly,
problem::Problem{Elasticity},
element::Element{El},
local_buffer::Elasticity3DLocalBuffers,
time, ::Type{Val{:continuum}}) where El<:Elasticity3DVolumeElements
time, ::Type{Val{:continuum}},
use_csc = false) where El<:Elasticity3DVolumeElements
props = problem.properties
dim = get_unknown_field_dimension(problem)
@@ -288,7 +284,6 @@ function assemble_element!(assembly::Assembly,
stress_vec[:] = Dtan * strain_vec
end
#=
if material_model == :ideal_plasticity
plastic_def = element("plasticity")[ip.id]
@@ -313,7 +308,6 @@ function assemble_element!(assembly::Assembly,
calculate_stress!(stress_vec, stress_last, dstrain_vec, plastic_strain, D, params, Dtan, yield_surface_, time, dt, Val{:type_3d})
end
=#
:strain in props.store_fields && update!(ip, "strain", time => strain_vec)
:stress in props.store_fields && update!(ip, "stress", time => stress_vec)
@@ -388,11 +382,19 @@ function assemble_element!(assembly::Assembly,
# Update f_ext in place to be f_ext - f_int
f_ext .-= f_int
# add contributions to K, Kg, f
FEMSparse.assemble_local!(assembler, gdofs, Km, f_ext)
if use_csc
# add contributions to K, Kg, f
FEMSparse.assemble_local!(assembler, gdofs, Km, f_ext)
if props.geometric_stiffness
FEMSparse.assemble_local_matrix!(assembler, gdofs, Kg)
if props.geometric_stiffness
FEMSparse.assemble_local_matrix!(assembler, gdofs, Kg)
end
else
add!(assembly.f, gdofs, f_ext)
add!(assembly.K, gdofs, gdofs, Km)
if props.geometric_stiffness
add!(assembly.Kg, gdofs, gdofs, Kg)
end
end
return nothing
@@ -445,7 +447,7 @@ function assemble!(assembly::Assembly,
end
gdofs = get_gdofs(problem, element)
FEMSparse.assemble_local_vector!(assembly.f, gdofs, f)
add!(assembly.f, gdofs, f)
end
end
+18 -11
View File
@@ -61,28 +61,35 @@ function get_field_assembly(solver::Solver)
M = problem.assembly.M
K = problem.assembly.K
f = problem.assembly.f
K_csc = problem.assembly.K_csc
f_csc = problem.assembly.f_csc
Kg = problem.assembly.Kg
fg = problem.assembly.fg
for problem in problems[2:end]
append!(M, problem.assembly.M)
K += problem.assembly.K
append!(K, problem.assembly.K)
append!(Kg, problem.assembly.Kg)
f += problem.assembly.f
append!(f, problem.assembly.f)
# Use in place addition with .+= ?
K_csc += problem.assembly.K_csc
f_csc += problem.assembly.f_csc
append!(fg, problem.assembly.fg)
end
N = size(K, 1)
M = sparse(M, N, N)
K = sparse(K, N, N)
if nnz(K) == 0
@warn("Field assembly seems to be empty. Check that elements are ",
"pushed to problem and formulation is correct.")
end
f = sparse(f, N, 1)
Kg = sparse(Kg, N, N)
fg = sparse(fg, N, 1)
return M, K, Kg, f, fg
return M, problem.assemble_csc ? K_csc : K, Kg, problem.assemble_csc ? f_csc : f, fg
end
""" Loop through boundary assemblies and check for possible overconstrain situations. """
@@ -141,11 +148,11 @@ function get_boundary_assembly(solver::Solver, N)
g = spzeros(N, 1)
for problem in get_boundary_problems(solver)
assembly = problem.assembly
# K_ = assembly.K
K_ = sparse(assembly.K, N, N)
C1_ = sparse(assembly.C1, N, N)
C2_ = sparse(assembly.C2, N, N)
D_ = sparse(assembly.D, N, N)
# f_ = assembly.f
f_ = sparse(assembly.f, N, 1)
g_ = sparse(assembly.g, N, 1)
for dof in assembly.removed_dofs
@info("$(problem.name): removing dof $dof from assembly")
@@ -166,12 +173,12 @@ function get_boundary_assembly(solver::Solver, N)
error("overconstrained dofs, not solving problem.")
end
#K += K_
C1 += C1_
C2 += C2_
D += D_
#f += f_
g += g_
K .+= K_
C1 .+= C1_
C2 .+= C2_
D .+= D_
f .+= f_
g .+= g_
end
return K, C1, C2, D, f, g
end