calculate stress, interpolate stress to nodes using lsq fitting

This commit is contained in:
Jukka Aho
2016-05-28 21:39:29 +03:00
parent ccbd2d7224
commit 11d8ebd970
9 changed files with 347 additions and 248 deletions
-51
View File
@@ -1,57 +1,6 @@
# This file is a part of JuliaFEM.
# License is MIT: see https://github.com/JuliaFEM/JuliaFEM.jl/blob/master/LICENSE.md
""" Calculate nodal vector from set of elements.
For example element 1 with dofs [1, 2, 3, 4] has [1, 1, 1, 1] and
element 2 with dofs [3, 4, 5, 6] has [2, 2, 2, 2] the result will
be sparse matrix with values [1, 1, 3, 3, 2, 2].
Parameters
----------
field_name
name of field, e.g. "geometry"
field_dim
degrees of freedom / node
elements
elements used to calculate vector
vec_dim
used to resize solution vector if given
time
"""
function calculate_nodal_vector(field_name, field_dim, elements::Vector{Element},
time, vec_dim=0)
A = SparseMatrixCOO()
b = SparseMatrixCOO()
for element in elements
haskey(element, field_name) || continue
gdofs = get_gdofs(element, 1)
for ip in get_integration_points(element, Val{3})
J = get_jacobian(element, ip, time)
w = ip.weight*norm(J)
f = element(field_name, ip, time)
N = element(ip, time)
add!(A, gdofs, gdofs, w*kron(N', N))
for dim=1:field_dim
add!(b, gdofs, w*f[dim]*N, dim)
end
end
end
A = sparse(A)
b = sparse(b)
nz = sort(unique(rowvals(A)))
x = zeros(size(b)...)
x[nz, :] = A[nz,nz] \ b[nz, :]
x = vec(transpose(x))
if vec_dim != 0
v = zeros(vec_dim)
v[1:length(x)] = x
return v
else
return x
end
end
function calculate_rotated_nodal_vector(field_name, field_dim, elements::Vector{Element},
time, vec_dim=0)
A = SparseMatrixCOO()