Element Mass Matrix¶

We have

$\mathbb{M}^e(T) =\frac{|T|}{20}\begin{pmatrix} 2 & 1 & 1 & 1\\ 1 & 2 & 1 & 1\\ 1 & 1 & 2 & 1\\ 1 & 1 & 1 & 2 \end{pmatrix}$

Then with $\MAT{K}_g$ definition (see Section New Optimized assembling algorithm (version OptV2)) , we obtain

$\MAT{K}_g(4(i-1)+j,k)=|T_k|\frac{1+\delta_{i,j}}{20} \quad 1\le i,j \le 4,$

So the vectorized algorithm for $\mathbb{K}_g$ computation is simple and given in Algorithm 149.

Algorithm 149

Figure 149: Vectorized algorithm for $\mathbb{K}_g$ associated to 3d Mass matrix

Note

pyOptFEM.FEM3D.elemMatrixVec.ElemMassMat3DP1Vec(nme, volumes)[source]

Compute all the elementaries Mass matrices, $\mathbb{M}^e(T_k)$ for $k\in\{0,\hdots,\nme-1\}$

Parameters:	volumes ( $\nme$ numpy array of floats) – volumes of all the mesh elements.
Returns:	a one dimensional numpy array of size $16 \nme$

Element Stiffness Matrix¶

We have $\forall (\il,\jl)\in\ENS{1}{4}^2$

$\mathbb{S}_{\il,\jl}^e(T_k)= |T_k| \DOT{\GRAD\BasisFunc_{\jl}^k}{\GRAD\BasisFunc_{\il}^k}.$

Using vectorized algorithm function $\FNametxt{GradientVec3D}$ given in Algorithm 148, we obtain the vectorized algorithm 150 for $\mathbb{K}_g$ computation of the Stiffness matrix in 3d.

Algorithm 150

Figure 150: Vectorized algorithm for $\mathbb{K}_g$ associated to 3d Stiffness matrix

Note

pyOptFEM.FEM3D.elemMatrixVec.ElemStiffMat3DP1Vec(nme, q, me, volumes)[source]

Compute all the elementaries Stiff matrices, $\mathbb{S}^e(T_k)$ for $k\in\{0,\hdots,\nme-1\}$

Parameters:	nme (int) – number of mesh elements, q ( $3\times \nq$ numpy array of floats) – mesh vertices, me ( $4 \times\nme$ numpy array of integers) – mesh connectivity, areas ( $\nme$ numpy array of floats) – areas of all the mesh elements.
Returns:	a one dimensional numpy array of size $9 \nme$

Element Stiffness Elasticity Matrix¶

We define on tetrahedra $T_k$ the local alternate basis $\mathcal{B}_a^k$ by

$\begin{array}{c} \mathcal{B}_a^k=\{\BasisFuncTwoD_1^k,\hdots,\BasisFuncTwoD_{12}^k\}\\=\\ \left\{\tiny \begin{pmatrix} \BasisFunc_1^k \\ 0 \\ 0\end{pmatrix}, \begin{pmatrix} 0 \\ \BasisFunc_1^k \\0 \end{pmatrix}, \begin{pmatrix} 0 \\ 0 \\ \BasisFunc_1^k \end{pmatrix}, \begin{pmatrix} \BasisFunc_2^k \\ 0 \\ 0\end{pmatrix}, \begin{pmatrix} 0 \\ \BasisFunc_2^k \\0 \end{pmatrix}, \begin{pmatrix} 0 \\ 0 \\ \BasisFunc_2^k \end{pmatrix}, \begin{pmatrix} \BasisFunc_3^k \\ 0 \\ 0\end{pmatrix}, \begin{pmatrix} 0 \\ \BasisFunc_3^k \\0 \end{pmatrix}, \begin{pmatrix} 0 \\ 0 \\ \BasisFunc_3^k \end{pmatrix}, \begin{pmatrix} \BasisFunc_4^k \\ 0 \\ 0\end{pmatrix}, \begin{pmatrix} 0 \\ \BasisFunc_4^k \\0 \end{pmatrix}, \begin{pmatrix} 0 \\ 0 \\ \BasisFunc_4^k \end{pmatrix} \right\} \end{array}$

where $\BasisFunc_\il^k=\BasisFunc_{\me(\il,k)}.$ With notations of Presentation, we have $\forall (\il,\jl) \in \ENS{1}{12}^2$

$\StiffElasElem_{\il,\jl}(T_k)= \int_{T_k} \Odv^t(\BasisFuncTwoD^k_\jl) \mathbb{C}\Odv(\BasisFuncTwoD^k_\il)d\q= \int_{T_k} \mathcal{H}(\BasisFuncTwoD^k_\jl,\BasisFuncTwoD^k_\il)(\q)d\q$

with, $\forall \vecb{u}=(u_1,u_2,u_3)\in\HUnD{\DOM}^3,$ $\forall \vecb{v}=(v_1,v_2,v_3)\in\HUnD{\DOM}^3,$ by

$\begin{array}{c} \mathcal{H}(\vecb{u},\vecb{v})\\=\\ \tiny{ \DOT{\begin{pmatrix} \gamma & 0 &0\\ 0 & \mu &0\\ 0 & 0 &\mu\end{pmatrix}\GRAD u_1 }{\GRAD v_1} +\DOT{\begin{pmatrix} 0 & \lambda & 0\\ \mu & 0 & 0 \\ 0 & 0 & 0 \end{pmatrix}\GRAD u_2 }{\GRAD v_1} +\DOT{\begin{pmatrix} 0 & 0 & \lambda\\ 0 & 0 & 0 \\ \mu & 0 & 0 \end{pmatrix}\GRAD u_3 }{\GRAD v_1} }\\ %&+& \tiny{+ \DOT{\begin{pmatrix} 0 & \mu &0\\ \lambda & 0 &0\\ 0 & 0 &0\end{pmatrix}\GRAD u_1 }{\GRAD v_2} +\DOT{\begin{pmatrix} \mu & 0 & 0\\ 0 & \gamma & 0 \\ 0 & 0 & \mu \end{pmatrix}\GRAD u_2 }{\GRAD v_2} +\DOT{\begin{pmatrix} 0 & 0 & 0\\ 0 & 0 & \lambda \\ 0 & \mu & 0 \end{pmatrix}\GRAD u_3 }{\GRAD v_2} }\\ %&+& \tiny{+ \DOT{\begin{pmatrix} 0 & 0 &\mu\\ 0 & 0 &0\\ \lambda & 0 & 0\end{pmatrix}\GRAD u_1 }{\GRAD v_3} +\DOT{\begin{pmatrix} 0 & 0 & 0\\ 0 & 0 & \mu \\ 0 & \lambda & 0 \end{pmatrix}\GRAD u_2 }{\GRAD v_3} +\DOT{\begin{pmatrix} \mu & 0 & 0\\ 0 & \mu & 0 \\ 0 & 0 & \gamma \end{pmatrix}\GRAD u_3 }{\GRAD v_3} } \end{array}$

where $\lambda$ and $\mu$ are the Lame coefficients, and $\gamma=\lambda+2\mu.$

For example, we can compute explicitely the first two terms in the first column of $\StiffElasElem(T_k)$ which are given by

$\begin{array}{lcl} \StiffElasElem_{1,1}(T_k)&=&\int_{T_k} \mathcal{H}(\BasisFuncTwoD^k_{1},\BasisFuncTwoD^k_{1})(\q)d\q\\ &=&\int_{T_k} \mathcal{H}\left( \tiny\begin{pmatrix} \BasisFunc^k_{1}\\ 0\\ 0 \end{pmatrix}, \tiny\begin{pmatrix} \BasisFunc^k_{1}\\ 0\\ 0 \end{pmatrix} \right)(\q)d\q\\ &=&|T_k| \DOT{\tiny\begin{pmatrix} \gamma & 0 &0\\ 0 & \mu &0\\ 0 & 0 &\mu\end{pmatrix}\GRAD \BasisFunc^k_{1} }{\GRAD \BasisFunc^k_{1}} =|T_k|\left(\gamma\DP{\BasisFunc^k_{1}}{x}\DP{\BasisFunc^k_{1}}{x}+\mu(\DP{\BasisFunc^k_{1}}{y}\DP{\BasisFunc^k_{1}}{y}+\DP{\BasisFunc^k_{1}}{z}\DP{\BasisFunc^k_{1}}{z}) \right). \end{array}$

and

$\begin{array}{lcl} \StiffElasElem_{2,1}(T_k)&=&\int_{T_k} \mathcal{H}(\BasisFuncTwoD^k_{1},\BasisFuncTwoD^k_{2})(\q)d\q\\ &=&\int_{T_k} \mathcal{H}\left( \tiny\begin{pmatrix} \BasisFunc^k_{1}\\ 0\\0 \end{pmatrix}, \begin{pmatrix} 0\\ \BasisFunc^k_{1}\\ 0 \end{pmatrix} \right)(\q)d\q\\ &=&|T_k| \DOT{\tiny\begin{pmatrix} 0 & \mu &0\\ \lambda & 0 &0\\ 0 & 0 &0\end{pmatrix}\GRAD \BasisFunc^k_{1} }{\GRAD \BasisFunc^k_{1}} =|T_k|(\lambda+\mu)\DP{\BasisFunc^k_{1}}{x}\DP{\BasisFunc^k_{1}}{y}. \end{array}$

Using vectorized algorithm function $\FNametxt{GradientVec3D}$ given in Algorithm 148, we obtain the vectorized algorithm 151 for $\mathbb{K}_g$ computation of the Elasticity Stiffness matrix in 3d.

Algorithm 151

Figure 151: Vectorized algorithm for $\mathbb{K}_g$ associated to 3d Elasticity Stiffness matrix

Note

pyOptFEM.FEM3D.elemMatrixVec.ElemStiffElasMatBa3DP1Vec(nme, q, me, volumes, la, mu)[source]

Compute all the elementaries Stiffness elasticity matrices, $\mathbb{K}^e(T_k)$ for $k\in\{0,\hdots,\nme-1\}$ in local alternate basis.

Parameters:	nme (int) – number of mesh elements, q (`(3,nq)` numpy array of floats) – mesh vertices, me (`(4,nme)` numpy array of integers) – mesh connectivity, volumes (`(nme,)` numpy array of floats) – volumes of all the mesh elements. la (float) – the $\\lambda$ Lame parameter, mu (float) – the $\\mu$ Lame parameter.
Returns:	a `(144nme,)` numpy* array of floats.

We define on $T_k$ the local block basis $\mathcal{B}_b^k$ by

$\begin{array}{c} \mathcal{B}_b^k=\{\BasisFuncTwoDB_1^k,\hdots,\BasisFuncTwoDB_{12}^k\} \\ = \\ \left\{\tiny \begin{pmatrix} \BasisFunc_1^k \\ 0 \\ 0 \end{pmatrix}, \begin{pmatrix} \BasisFunc_2^k \\ 0 \\ 0 \end{pmatrix}, \begin{pmatrix} \BasisFunc_3^k \\ 0 \\ 0 \end{pmatrix}, \begin{pmatrix} \BasisFunc_4^k \\ 0 \\ 0 \end{pmatrix}, \begin{pmatrix} 0 \\ \BasisFunc_1^k \\ 0 \end{pmatrix}, \begin{pmatrix} 0 \\ \BasisFunc_2^k \\ 0 \end{pmatrix}, \begin{pmatrix} 0 \\ \BasisFunc_3^k \\ 0 \end{pmatrix}, \begin{pmatrix} 0 \\ \BasisFunc_4^k \\ 0 \end{pmatrix}, \begin{pmatrix} 0 \\ 0 \\ \BasisFunc_1^k \end{pmatrix}, \begin{pmatrix} 0 \\ 0 \\ \BasisFunc_2^k \end{pmatrix}, \begin{pmatrix} 0 \\ 0 \\ \BasisFunc_3^k \end{pmatrix}, \begin{pmatrix} 0 \\ 0 \\ \BasisFunc_4^k \end{pmatrix} \right\} \end{array}$

where $\BasisFunc_\il^k=\BasisFunc_{\me(\il,k)}.$

For example, using formula (?), we can explicitly compute the first two terms in the first column of $\StiffElasElem(T_k)$ which are given by

$\begin{array}{lcl} \StiffElasElem_{1,1}(T_k)&=&\int_{T_k} \mathcal{H}(\BasisFuncTwoDB^k_{1},\BasisFuncTwoDB^k_{1})(\q)d\q\\ &=&\int_{T_k} \mathcal{H}\left( \tiny\begin{pmatrix} \BasisFunc^k_{1}\\ 0\\0 \end{pmatrix}, \begin{pmatrix} \BasisFunc^k_{1}\\ 0\\0 \end{pmatrix} \right)(\q)d\q\\ &=&|T_k| \DOT{\tiny\begin{pmatrix} \gamma & 0 &0\\ 0 & \mu &0\\ 0 & 0 &\mu\end{pmatrix}\GRAD \BasisFunc^k_{1} }{\GRAD \BasisFunc^k_{1}} =|T_k|\left(\gamma\DP{\BasisFunc^k_{1}}{x}\DP{\BasisFunc^k_{1}}{x}+\mu(\DP{\BasisFunc^k_{1}}{y}\DP{\BasisFunc^k_{1}}{y} +\DP{\BasisFunc^k_{1}}{z}\DP{\BasisFunc^k_{1}}{z})\right). \end{array}$

and

$\begin{array}{lcl} \StiffElasElem_{2,1}(T_k)&=&\int_{T_k} \mathcal{H}(\BasisFuncTwoDB^k_{1},\BasisFuncTwoDB^k_{2})(\q)d\q\\ &=&\int_{T_k} \mathcal{H}\left( \tiny \begin{pmatrix} \BasisFunc^k_{1}\\ 0\\ 0 \end{pmatrix}, \begin{pmatrix} \BasisFunc^k_{2}\\ 0\\ 0 \end{pmatrix} \right)(\q)d\q\\ &=&|T_k| \DOT{\tiny\begin{pmatrix} \gamma & 0 &0\\ 0 & \mu &0\\ 0 & 0 &\mu\end{pmatrix}\GRAD \BasisFunc^k_{1} }{\GRAD \BasisFunc^k_{2}} =|T_k|\left(\gamma\DP{\BasisFunc^k_{1}}{x}\DP{\BasisFunc^k_{2}}{x}+\mu(\DP{\BasisFunc^k_{1}}{y}\DP{\BasisFunc^k_{2}}{y}+\DP{\BasisFunc^k_{1}}{z}\DP{\BasisFunc^k_{2}}{z}) \right). \end{array}$

Using vectorized algorithm function $\FNametxt{GradientVec3D}$ given in Algorithm 148, we obtain the vectorized algorithm 152 for $\mathbb{K}_g$ computation of the Elasticity Stiffness matrix in 3d.

Algorithm 152

Figure 152: Vectorized algorithm for $\mathbb{K}_g$ associated to 3d Elasticity Stiffness matrix

Note

pyOptFEM.FEM3D.elemMatrixVec.ElemStiffElasMatBb3DP1Vec(nme, q, me, volumes, L, M)[source]

Compute all the elementaries Stiffness elasticity matrices, $\mathbb{K}^e(T_k)$ for $k\in\{0,\hdots,\nme-1\}$ in local block basis.

Parameters:	nme (int) – number of mesh elements, q (`(3,nq)` numpy array of floats) – mesh vertices, me (`(4,nme)` numpy array of integers) – mesh connectivity, volumes (`(nme,)` numpy array of floats) – volumes of all the mesh elements. la (float) – the $\\lambda$ Lame parameter, mu (float) – the $\\mu$ Lame parameter.
Returns:	a `(144nme,)` numpy* array of floats.

Element Mass Matrix¶

Element Stiffness Matrix¶

Element Stiffness Elasticity Matrix¶

Table Of Contents

This Page

Navigation

Element Mass Matrix¶

Element Stiffness Matrix¶

Element Stiffness Elasticity Matrix¶

Table Of Contents

This Page

Quick search

Navigation