LFEM.jl

Solve the linear system K(x)U(x)=F

Solve_linear(mesh::Mesh, x::Vector{Float64}, kparam::Function; loadcase=1)

where

x is a ne x 1 vector of design varibles 
kparam(xe): R->R is the material parametrization (SIMP like)
loadcase is the loadcase

Returns:

U = displacement vector (dim*nn x 1)
F = Force vector (dim*nn x 1)
linsolve = LinearSolve object with factored linear problem

Solve the linear system KU=F

Solve_linear(mesh::Mesh; loadcase=1)

Returns:

U = displacement vector (dim*nn x 1)
F = Force vector (dim*nn x 1)
linsolve = LinearSolve object with factored linear problem

Solve the modal problem (M(x) - λK(x))ϕ(x) = 0

Solve_modal(mesh::Mesh, x::Vector{Float64}, kparam::Function; mparam::Function, nev=4, which=:SM, σ=1.0, loadcase::Int64=1)

where

x is a ne x 1 vector of design varibles 
kparam(xe): R->R is the material parametrization for K (SIMP like)
mparam(xe): R->R is the material parametrization for M (SIMP like)
nev is the number of eigenvalues and eigenvectors to compute
lumped is true for lumped mass matrix
loadcase is the loadcase

Returns:

λ = eigenvalues vector (nev x 1)
modes = matrix dim*nn x nev with the eigenvectors

Solve the modal problem (M - λK)ϕ = 0

Solve_modal(mesh::Mesh ;nev=4, loadcase=1)

where nev is the number of eigenvalues and eigenvectors to compute loadcase is the loadcase

Returns:

λ = eigenvalues vector (nev x 1)
modes = matrix dim*nn x nev with the eigenvectors

Solve the harmonic problem Kd(x,w)Ud(x,w) = F(w), where Kd(x,w)= K(x)-M(x)w^2 + imwC(x)

Solve_harmonic(mesh::Mesh, w::Float64, α_c::Float64, β_c::Float64, x::Vector{Float64}, 
               kparam::Function, mparam::Function ; 
               lumped=true,loadcase=1)

where

w is the angular frequency
α_c and  β_c are the parameters for proportional damping
x is a ne x 1 vector of design varibles 
kparam(xe): R->R is the material parametrization for K (SIMP like)
mparam(xe): R->R is the material parametrization for M (SIMP like)
lumped is true for lumped mass matrices
loadcase is the loadcase

Returns:

Ud = displacement vector (ComplexF64) of size dim*nn x 1
linsolve = LinearSolve object with factored linear problem

Solve the harmonic problem Kd(w)Ud(w) = F(w), where Kd(w)= K-Mw^2 + imwC

Solve_harmonic(mesh::Mesh, w::Float64, α_c::Float64, β_c::Float64 ; loadcase=1)

where

w is the angular frequency  
α_c and  β_c are the parameters for proportional damping
lumped=true is for lumped mass matrix
loadcase is the loadcase

Returns:

Ud = displacement vector (ComplexF64) of size dim*nn x 1  
linsolve = LinearSolve object with factored linear problem

Solve the transient problem M(x)A(x,t) + C(x)V(x,t) + K(x,t)U(x,t) = F(t), using Newmark-beta method.

Solve_newmark(mesh::Mesh, f!::Function, gls::Matrix{Int64}, 
              ts::Tuple{Float64, Float64}, Δt::Float64,
              x::Vector{Float64}, kparam::Function, mparam::Function,
              verbose=false;
              U0=Float64[], V0=Float64[],
              β=1/4, γ=1/2, 
              α_c=0.0, β_c=1E-6
              loadcase=1)

where

ts is Tupple with initial and end time (Ti,Tf)  
Δt is (fixed) time step  
x is a ne x 1 vector of design varibles   
kparam(xe): R->R is the material parametrization for K (SIMP like)  
mparam(xe): R->R is the material parametrization for M (SIMP like)  
verbose is false or true
U0 and V0 are the initial conditions  
β and γ are the parameters of the Newmark method
α_c and  β_c are the coefficients for proportional damping C=α_cM + β_c*K
lumped is true for lumped mass matrix
loadcase is the loadcase

f!(t,F,mesh,loadcase) must be a function of t, mesh and F where F is dim*nn x 1,
            Example: 

            function f!(t,F,mesh::Mesh,loadcase=1)
                        P = Point_load(mesh,loadcase)
                        F.= cos(2*t)*P 
            end   

gls is a matrix with [node gl ;
                      node gl ...] to monitor

Return three arrays of size ng x nt, where ng is size(gls,1) and nt is the number of time steps (length of t0:Δt:tf)

A_U displacements
A_V velocities
A_A accelerations

A_t is a vector of size nt x 1 with discrete times
dofs is a vector with the (global) monitored dofs

Solve the transient problem MA(t) + CV(t) + K(t)U(t) = F(t), using Newmark-beta method.

Solve_newmark(mesh::Mesh, f!::Function, gls::Matrix{Int64}, 
              ts::Tuple{Float64, Float64}, Δt::Float64,
              verbose=false;
              U0=Float64[], V0=Float64[], β=1/4, γ=1/2,loadcase=1)

where

ts is Tupple with initial and end time (Ti,Tf)
Δt is (fixed) time steps
verbose is false or true
U0 and V0 are the initial conditions  
β and γ are the parameters of the Newmar method
lumped is true for lumped mass matrix
loadcase is the loadcase

f!(t,F,mesh,loadcase) must be a function of t, mesh and F where F is dim*nn x 1,
            Example: 

            function f!(t,F,mesh::Mesh,loadcase=1)
                        P = Point_load(mesh,loadcase)
                        F.= cos(2*t)*P 
            end   

gls is a matrix with [node gl ;
                      node gl ...] to monitor

Return three arrays of size ng x nt, where ng is size(gls,1) and nt is the number of time steps (length of t0:Δt:tf)

A_U displacements
A_V velocities
A_A accelerations

A_t is a vector of size nt x 1 with discrete times
dofs is a vector with the (global) monitored dofs

Solve the transient problem M(x)A(x,t) + C(x)V(x,t) + K(x,t)U(x,t) = F(t), using Newmark-beta method. Matrizes are given and no mesh information is used

Solve_newmark(M::AbstractMatrix,C::AbstractMatrix,K::AbstractMatrix, f!::Function, gls::Matrix{Int64}, 
              ts::Tuple{Float64, Float64}, Δt::Float64,
              verbose=false;
              U0=Float64[], V0=Float64[],
              β=1/4, γ=1/2)

where

M, C and K are given matrices
ts is Tupple with initial and end time (Ti,Tf)  
Δt is (fixed) time step  
verbose is false or true
U0 and V0 are the initial conditions  
β and γ are the parameters of the Newmark method

f!(t,F) must be a function of t  and F where F is dim*nn x 1,
            Example: 

            function f!(t,F)
                        P = [0.0;10.0]
                        F.= cos(2*t)*P 
            end   

gls is a vector gls to monitor

Return three arrays of size ng x nt, where ng is size(gls,1) and nt is the number of time steps (length of t0:Δt:tf)

A_U displacements
A_V velocities
A_A accelerations

A_t is a vector of size nt x 1 with discrete times

Assembly the global stiffness matrix.

Global_K(mesh::Mesh, x::Vector{Float64}, kparam::Function)

where kparam(x) can be, for example

function kparam(x,p=1.0) x^p end

for a SIMP like material parametrization.

This function also considers entries :Stiffness in mesh.options with [node dof value;]

Assembly the global stiffness matrix.

Global_K(mesh::Mesh)

This function also considers entries :Stiffness in mesh.options with [node dof value;]

Assembly the global mass matrix.

 Global_M(mesh::Mesh, x=Float64[], mparam::Function; lumped=true)

where mparam(x) can be, for example

function param(x,p=2.0,cut=0.1) s = x if x<cut s = x^p end return s end

for a SIMP like material parametrization.

This function also considers entries :Mass in mesh.options with [node dof value;]

Assembly the global mass matrix.

 Global_M(mesh::Mesh)

This function also considers entries :Stiffness in mesh.options with [node dof value;]

Assembly the global damping matrix C = αc M + βc K

 Global_C(M,K,mesh::Mesh,α_c=0.0,β_c=0.0)

This function also considers entries :Damper in mesh.options with [node dof value;]

Return stresses for the entire mesh. This version evaluates only in the central point.

Stresses(mesh::Mesh,U::Vector{T};x=Float64[],sparam::Function; center=true)

The output is a matrix ne x ncol, where ncol is 1 for :truss2D and 3D, 3 for :solid2D and 6 for :solid3D if center = true and ncol = 34 for 2D and 68 for 3D if center = false since we return stresses for each superconvergent (Gauss) point in the incompatible elements

Function sparam(x) is used to parametrize stress. One possibility is

function sparam(x,p=1.0,q=0.0) x^(p=q) end

Return stresses for the entire mesh. This version evaluates only in the central point.

Stresses(mesh::Mesh,U::Vector{T}; center=true)

The output is a matrix ne x ncol, where ncol is 1 for :truss2D and 3D, 3 for :solid2D and 6 for :solid3D if center = true and ncol = 34 for 2D and 68 for 3D if center = false since we return stresses for each superconvergent (Gauss) point in the incompatible elements

Return (harmonic) stresses for the entire mesh. This version evaluates only in the central point.

Harmonicstresses(mesh::Mesh,U::Vector{T}, w::Float64, βc::Float64, x=Float64[],sparam::Function; center=true)

where w is the angular frequency and β_c is the damping parameter.

The output is a matrix ne x ncol, where ncol is 1 for :truss2D and 3D, 3 for :solid2D and 6 for :solid3D if center = true and ncol = 34 for 2D and 68 for 3D if center = false since we return stresses for each superconvergent (Gauss) point in the incompatible elements

Function sparam(x) is used to parametrize stress. One possibility is

function sparam(x,p=1.0,q=0.0) x^(p=q) end

Return stresses for the entire mesh. This version evaluates only in the central point.

Harmonic_stresses(mesh::Mesh,U::Vector{T}, w::Float64, β_c::Float64)

where w is the angular frequency and β_c is the damping parameter.

The output is a matrix ne x ncol, where ncol is 1 for :truss2D and 3D, 3 for :solid2D and 6 for :solid3D if center = true and ncol = 34 for 2D and 68 for 3D if center = false since we return stresses for each superconvergent (Gauss) point in the incompatible elements

Local stiffness matrix for truss2D K_truss2D(mesh::Mesh2D,ele::Int64)

Local mass matrix for truss2D M_truss2D(mesh::Mesh2D,ele::Int64; lumped=true)

Local B matrix for truss2D B_truss2D(mesh::Mesh2D,ele::Int64)

Local stress for truss2D Stress_truss2D(mesh::Mesh2D,ele::Int64,U::Vector{T})

It returns stress as [sxx] for compatibility with solid elements.

Local stiffness matrix for truss3D K_truss3D(mesh::Mesh3D,ele::Int64)

Local mass matrix for truss3D (lumped) M_truss3D(mesh::Mesh3D,ele::Int64; lumped=true)

Local B matrix for truss3D B_truss3D(mesh::Mesh3D,ele::Int64)

Local stress for truss3D Stress_truss3D(mesh::Mesh3D,ele::Int64,U::Vector{T})

It returns stress as [sxx] for compatibility with solid elements.

B Matrix (with additional bublle functions) for 2D elements B_solid2D(r::T,s::T,x::Vector{T},y::Vector{T}) where T

Stiffness Matrix for (incompatible) 2D element K_solid2D(m::Mesh2D,ele::Int64)

Consistent mass matrix for solid 2D M_solid2D(m::Mesh2D,ele::Int64,lumped=false)

Local stress for solid 2D ((Not expanding bubble DOFs) Stress_solid2D(r::Float64,s::Float64,mesh::Mesh2D,ele::Int64,U::Vector{T})

Return the volume of element ele

Volume_solid2D(mesh::Mesh,ele::Int64)

B Matrix (with additional bublle functions) for 3D elements B_solid3D(r::T,s::T,t::T,x::Vector{T},y::Vector{T},z::Vector{T}) where T

Stiffness Matrix for (incompatible) 3D element K_solid3D(m::Mesh3D,ele::Int64)

Consistent mass matrix for solid 3D M_solid3D(m::Mesh3D,ele::Int64,lumped=false)

Local stress for solid 3D ((Not expanding bubble DOFs) Stress_solid3D(r::Float64,s::Float64,t::Float64,mesh::Mesh2D,ele::Int64,U::Vector{T})

Return the volume of element ele

Volume_solid3D(mesh::Mesh,ele::Int64)

(Overloaded from BMesh) Initialize a gmsh mesh for post processing.

Gmsh_init(nome_arquivo::String,mesh::Mesh)

Export a nodal scalar view to gmsh

Gmsh_nodal_scalar(mesh::Mesh,escalares::Vector,nome_arquivo::String,
                  nome_vista::String,tempo=0.0)

Export an element (centroidal) scalar view to gmsh

Gmsh_element_scalar(mesh::Mesh,escalares::Vector,nome_arquivo::String,
                    nome_vista::String,tempo=0.0)

Export an nodal vectorial view to gmsh

Gmsh_nodal_vector(mesh::Mesh,vetor::Vector,nome_arquivo::String,
                  nome_vista::String,tempo=0.0)