Choral
Data Types | Functions/Subroutines
elasticity Module Reference

MODULE FOR LINEAR ELASTICITY PROBLEMS More...

Data Types

interface  l2_product
 L2 product of a function with the basis functions. More...
 

Functions/Subroutines

subroutine, public elasticity_matrix_pattern (g, Y, qdm)
  Define the sparsity pattern for elasticity matrices More...
 
subroutine, public elasticity_massmat (mass, a, Y, qdm, dofToDof)
  Assemble the mass matrix of the bilinear product: More...
 
subroutine, public elasticity_stiffmat (stiff, lambda, mu, Y, qdm, dofToDof)
  Assemble the stiffness matrix of the bilinear product: More...
 
real(rp) function, dimension(2, 2) e2 (v, component)
 To compute the symmetrised gradient in dim 2. More...
 
real(rp) function, dimension(3, 3) e3 (v, component)
 To compute the symmetrised gradient in dim 3. More...
 
subroutine fespacexk_l2_product (FV, Y, qdm, f_1, f_2, f_3)
  L2 Product of a vector function \( f:~\R^3 \mapsto \R^d\) with the basis functions of a finite element space \( Y = [X_h]^d \)..
\( X_h \) is a scalar finite element space on a mesh \( \T \) with dimension \( d \). It us assumed that \( d \) =2, 3 and that \( \Omega \subset \R^d \). More...
 
subroutine, public elasticity_neumann_rhs (rhs, Y, quad_type, g_1, g_2, g_3, f)
  L2 scalar product of a vector function \( g:~\R^3 \mapsto \R^d \) with the basis functions fo the finite element space \( Y = [X_h]^d \) on \(\Gamma_f \subset \partial\Omega\) a part of the domain boundary.
\( X_h \) is a scalar finite element space on a mesh \( \T \) with dimension \( d \). It us assumed that \( d \) =2, 3 and that \( \Omega \subset \R^d \). More...
 
subroutine, public elasticity_dirichlet (K, rhs, Y, g1, g2, g3, f)
  DIRICHLET BOUNDARY CONDITION FOR AN ELASTICITY PROBLEM More...
 

Detailed Description

MODULE FOR LINEAR ELASTICITY PROBLEMS

Construction of finite element matrices for linear elasticity.

Linear elasticity problems: search for \(u:~\Omega \mapsto \R^2\)

\(~~~~~~~~~ -\dv(A(x) e(u)) = f ~~~\) on \(~~~ \Omega ~~~~~~~~~~~~~~~~~~~~~~~~~~~(1)\)


Hooke tensor: \( ~~A(x)\xi = \lambda(x) {\rm Tr}(\xi) Id + 2 \mu(x)\xi ~~\) with \(~~ \lambda,~\mu~: \Omega \mapsto \R \).


Symmetrised gradient: \(~~ e(u) = (\nabla u + ^T\nabla u)/2 \).

Tutorial examples:

Author
Charles Pierre

Function/Subroutine Documentation

◆ e2()

real(rp) function, dimension(2,2) elasticity::e2 ( real(rp), dimension(2), intent(in)  v,
integer, intent(in)  component 
)
private

To compute the symmetrised gradient in dim 2.

Definition at line 661 of file elasticity.F90.

Here is the call graph for this function:

◆ e3()

real(rp) function, dimension(3,3) elasticity::e3 ( real(rp), dimension(3), intent(in)  v,
integer, intent(in)  component 
)
private

To compute the symmetrised gradient in dim 3.

Definition at line 695 of file elasticity.F90.

Here is the call graph for this function:

◆ elasticity_dirichlet()

subroutine, public elasticity::elasticity_dirichlet ( type(csr), intent(inout)  K,
real(rp), dimension(:), intent(inout)  rhs,
type(fespacexk), intent(in)  Y,
procedure(r3tor)  g1,
procedure(r3tor)  g2,
procedure(r3tor), optional  g3,
procedure(r3tor), optional  f 
)

DIRICHLET BOUNDARY CONDITION FOR AN ELASTICITY PROBLEM

  • \( X_h \) is a scalar finite element space on a mesh \( \T \) with dimension \( d \).
  • It is assumed that \( d \) =2, 3 and that \( \Omega \subset \R^d \).

INPUT:

  • the boundary sub-domain \( \Gamma_f = \partial\Omega \cap \{f\ge 0\}, \)
  • the source term \( g \) on \( \Gamma_f \mapsto \R^d\). with \( d \) the mesh dimension and with:
    \( g = [g_1, g_2, g_3] \) if \( d=3\) or \( g = [g_1, g_2]\) if \( d=2\)
  • \( Y= [X_h]^d \) with \( [X_h] \) a scalar finite element space.
  • the linear ssytem

    \[ K U = \text{rhs}, \]

OUTPUT:
Modify the linear system on the lines \( l \) associated with a finite element node \( x_i\in \Gamma_f \):

  • for \( l=(d-1)i + j\), \( ~j=1 \dots d \):
  • \( \text{rhs}_l = g_j(x_i)\),
  • \( K_{l,m} = \delta_{lm} \)

Definition at line 992 of file elasticity.F90.

Here is the call graph for this function:

◆ elasticity_massmat()

subroutine, public elasticity::elasticity_massmat ( type(csr), intent(inout)  mass,
procedure(r3tor)  a,
type(fespacexk), intent(in)  Y,
type(quadmesh), intent(in)  qdm,
type(graph), intent(in), optional  dofToDof 
)

Assemble the mass matrix of the bilinear product:

\( (u,v) \mapsto \int_O a(x)~ u(x) \cdot v(x) \dx \)

  • \( O \) = integration domain, see quadmesh_mod for a definition
  • \( u, v \in Y = [X_h]^d\),
  • \( X_h =\) scalar finite element space on the mesh \( \T \),
  • \( d =\) dimension of the mesh \( \T \),

OUTPUT:

  • mass= mass matrix

INPUT :

  • \( a:~x \in\R^3 \mapsto a(x)\in\R \), density function
  • qdm = integration method on the mesh \( \T \)
  • dofToDof = connectivity graph DOF –> DOF

Definition at line 175 of file elasticity.F90.

Here is the call graph for this function:

◆ elasticity_matrix_pattern()

subroutine, public elasticity::elasticity_matrix_pattern ( type(graph), intent(inout)  g,
type(fespacexk), intent(in)  Y,
type(quadmesh), intent(in)  qdm 
)

Define the sparsity pattern for elasticity matrices


It is given by the connectivity graph Dof –> Dof:
between the degrees of freedom of the space \( Y = [X_h]^d \)
where \( X_h \) is a finite element space on the mesh \( \T \) which is of dimension \( d \).


DOF_i is in relation with DOF_j if:

  • they are asociated to a common cell of the mesh \( \T \),
  • and if that common cell is associated to a non-void quadrature rule in the integration method.

OUTPUT:

  • g = connevtivity graph

INPUT :

  • \( X_h \) finite element space on the mesh \( \T \)
  • qdm = integration method on the mesh \( \T \)

Definition at line 103 of file elasticity.F90.

Here is the call graph for this function:

◆ elasticity_neumann_rhs()

subroutine, public elasticity::elasticity_neumann_rhs ( real(rp), dimension(:), allocatable  rhs,
type(fespacexk), intent(in)  Y,
integer, intent(in)  quad_type,
procedure(r3tor)  g_1,
procedure(r3tor)  g_2,
procedure(r3tor), optional  g_3,
procedure(r3tor), optional  f 
)

L2 scalar product of a vector function \( g:~\R^3 \mapsto \R^d \) with the basis functions fo the finite element space \( Y = [X_h]^d \) on \(\Gamma_f \subset \partial\Omega\) a part of the domain boundary.
\( X_h \) is a scalar finite element space on a mesh \( \T \) with dimension \( d \). It us assumed that \( d \) =2, 3 and that \( \Omega \subset \R^d \).

INPUT:

  • \( g:~\R^3 \mapsto \R^d \) with \( d \) the mesh dimension and with:
    \( g = [g_1, g_2, g_3] \) if \( d=3\) or \( g = [g_1, g_2]\) if \( d=2\).
  • \( Y= [X_h]^d \) with \( [X_h] \) a scalar finite element space.
  • 'quad_type' quadrature rule for the boundary cells \( K \subset \Gamma_f\).
  • \( f:~\R^3 \mapsto \R \) a function to characterise the boundary sub-domain \(\Gamma_f = \partial\Omega \cap \{f\ge 0\} \).
    This is an optional argument, if not provided \( \Gamma_f = \partial\Omega \).

OUTPUT:

  • \( {\rm rhs} = (g_i)_{1\le i\le dN} \) with \( g_{(i-1)d+c} = \int_{\Gamma_f} g_c(x) v_i \dx \), where \( (v_i)_{1\le i\le N} \) are the basis functions of \( X_h\) and for \( 1\le c \le d \).

Definition at line 886 of file elasticity.F90.

Here is the call graph for this function:

◆ elasticity_stiffmat()

subroutine, public elasticity::elasticity_stiffmat ( type(csr), intent(inout)  stiff,
procedure(r3tor)  lambda,
procedure(r3tor)  mu,
type(fespacexk), intent(in)  Y,
type(quadmesh), intent(in)  qdm,
type(graph), intent(in), optional  dofToDof 
)

Assemble the stiffness matrix of the bilinear product:

\((u,v) \mapsto \int_O A(x)e(u)(x):e(v)(x) \dx~~\in \R\)

  • \( O \) = integration domain, see quadmesh_mod for a definition
  • \( u, v \in Y = [X_h]^d\),
  • \( X_h =\) scalar finite element space on the mesh \( \T \),
  • \( d =\) dimension of the mesh \( \T \),
  • \( e(u) = (\nabla(u) + ^T\nabla(u))/2 \in \R^d\times \R^d\) is the symmetrised gradient
  • \( A(x)\xi = \lambda(x) \Tr(\xi) Id + 2 \mu(x) \xi \in \R^d\times \R^d\)
    is the Hooke tensor
    with \( \lambda(x)\in \R \) and \( \mu(x)\in \R \) the Lame coefficients,
    and with \( \Tr(\xi) \) the trace of the matrix \( \xi \), \( Id \) the identity matrix.
  • \( M:N = \sum_{1 \le i,j \le d} M_{ij} N_{ij} \in \R \)
    for \( M,~ N\) two matrices in \( \R^d\times \R^d\).
  • Assumption: \( O \subset \R^d\)

OUTPUT:

INPUT :

  • \( \lambda:~x\in\R^3 \mapsto \lambda(x)\in\R \),
    \( \mu:~x\in\R^3 \mapsto \mu(x)\in\R \)
    the Lame coefficients (given as scalar functions),
  • dofToDof = connectivity graph DOF –> DOF (optional, of type graph_mod::graph)

Definition at line 402 of file elasticity.F90.

Here is the call graph for this function:

◆ fespacexk_l2_product()

subroutine elasticity::fespacexk_l2_product ( real(rp), dimension(:), allocatable  FV,
type(fespacexk), intent(in)  Y,
type(quadmesh), intent(in)  qdm,
procedure(r3tor)  f_1,
procedure(r3tor)  f_2,
procedure(r3tor), optional  f_3 
)
private

L2 Product of a vector function \( f:~\R^3 \mapsto \R^d\) with the basis functions of a finite element space \( Y = [X_h]^d \)..
\( X_h \) is a scalar finite element space on a mesh \( \T \) with dimension \( d \). It us assumed that \( d \) =2, 3 and that \( \Omega \subset \R^d \).

INPUT

  • \( f:~ \R^3 \mapsto \R^d\) with \( d \) the mesh dimension and with:
    \( f = [f_1, f_2, f_3] \) if \( d=3\) or \( f = [f_1, f_2]\) if \( d=2\)
  • \( Y = [X_h]^d \) with \( X_h \) a scalar finite element space on the mesh \( \T \).
  • qdm = integration method on the mesh \( \T \)

OUTPUT:

  • \( FV = (FV_i)_{1\le i\le dN} \) with \( FV_{(i-1)d+c} = \int_O f_c(x)\cdot v_i \dx = [ (f_c,v_i)_{{\rm L}^2(O)} ] \)
    where \( (v_i)_{1\le i\le N} \) are the basis functions of \( X_h\) and for \( 1\le c \le d \).
  • \( O \) = integration domain, see quadmesh_mod for a definition.

Definition at line 786 of file elasticity.F90.

Here is the call graph for this function: