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 assembly algorithm (OptV2 version)) , 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 27.

Algorithm 27

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

Note

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

Computes all the element 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 26, we obtain the vectorized algorithm 28 for $\mathbb{K}_g$ computation for the Stiffness matrix in 3d.

Algorithm 28

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

Note

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

Computes all the element stiffness 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 Elastic Stiffness Matrix¶

We define on the tetrahedron $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 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}(\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 26, we obtain the vectorized algorithm 29 for $\mathbb{K}_g$ computation for the Elastic Stiffness matrix in 3d.

Algorithm 29

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

Note

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

Computes all the element elastic stiffness 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 26, we obtain the vectorized algorithm 30 for $\mathbb{K}_g$ computation for the Elastic Stiffness matrix in 3d.

Algorithm 30

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

Note

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

Compute all the element elastic stiffness 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 Elastic Stiffness Matrix¶

Table Of Contents

This Page

Navigation

Element Mass Matrix¶

Element Stiffness Matrix¶

Element Elastic Stiffness Matrix¶

Table Of Contents

This Page

Quick search

Navigation