Damping matrix 3D

Let us proceed to computing the damping matrix for the three-dimensional case. The damping matrix is related to the integral of the square of the trial function. Consider the integral for a single tetrahedron with vertices

(xi,yi,zi),\displaystyle (x_i, y_i, z_i),(xi+1,yi+1,zi+1),\displaystyle (x_{i+1}, y_{i+1}, z_{i+1}),(xi+2,yi+2,zi+2),\displaystyle (x_{i+2}, y_{i+2}, z_{i+2}),(xi+3,yi+3,zi+3)\displaystyle (x_{i+3}, y_{i+3}, z_{i+3})
tetυ2dV.\int_{\text{tet}} \upsilon^2 \,dV.
(6.29)

The trial function on the tetrahedron has the form υ(i)(i+3)(x,y,z)\upsilon_{(i)(i+3)}(x, y, z)=qiϕi{} = q_i \cdot \phi_i+qi+1ϕi+1{} + q_{i+1} \cdot \phi_{i+1}+qi+2ϕi+2{} + q_{i+2} \cdot \phi_{i+2}+qi+3ϕi+3{} + q_{i+3} \cdot \phi_{i+3}. Analogously to the two-dimensional case, the hat functions for the tetrahedron are linear according to (6.13)

ϕi(x,y,z)\displaystyle \phi_i(x, y, z) =ai\displaystyle {} = a_i+bix\displaystyle {} + b_i \cdot x+ciy\displaystyle {} + c_i \cdot y+diz\displaystyle {} + d_i \cdot zϕi+1(x,y,z)\displaystyle \phi_{i+1}(x, y, z) =ai+1\displaystyle {} = a_{i+1}+bi+1x\displaystyle {} + b_{i+1} \cdot x+ci+1y\displaystyle {} + c_{i+1} \cdot y+di+1z\displaystyle {} + d_{i+1} \cdot zϕi+2(x,y,z)\displaystyle \phi_{i+2}(x, y, z) =ai+2\displaystyle {} = a_{i+2}+bi+2x\displaystyle {} + b_{i+2} \cdot x+ci+2y\displaystyle {} + c_{i+2} \cdot y+di+2z\displaystyle {} + d_{i+2} \cdot zϕi+3(x,y,z)\displaystyle \phi_{i+3}(x, y, z) =ai+3\displaystyle {} = a_{i+3}+bi+3x\displaystyle {} + b_{i+3} \cdot x+ci+3y\displaystyle {} + c_{i+3} \cdot y+di+3z\displaystyle {} + d_{i+3} \cdot z
(6.30)

Substitute the trial function into (6.29)

tetυ(i)(i+3)2dV\displaystyle \int_{\text{tet}} \upsilon_{(i)(i+3)}^2 \,dV =tet[qiϕi(x,y,z)+qi+1ϕi+1(x,y,z)+qi+2ϕi+2(x,y,z)+qi+3ϕi+3(x,y,z)]2dV\displaystyle {} = \int_{\text{tet}} \Big[ q_i \cdot \phi_i(x, y, z) + q_{i+1} \cdot \phi_{i+1}(x, y, z) + q_{i+2} \cdot \phi_{i+2}(x, y, z) + q_{i+3} \cdot \phi_{i+3}(x, y, z) \Big]^2 \,dV

For linear hat functions on the tetrahedron the following relations hold

tetϕmϕndV\displaystyle \int_{\text{tet}} \phi_m \cdot \phi_n \,dV ={Vtet10,m=nVtet20,mn\displaystyle {} = \begin{cases} \frac{\displaystyle V_{\text{tet}}}{\displaystyle 10}, & m = n\\ \frac{\displaystyle V_{\text{tet}}}{\displaystyle 20}, & m \neq n \end{cases}
(6.31)

where VtetV_{\text{tet}} is the volume of the tetrahedron, which is computed by formula (6.15).

Expand the square of the trial function

tetυ(i)(i+3)2dV\displaystyle \int_{\text{tet}} \upsilon_{(i)(i+3)}^2 \,dV =tet[qi2ϕi2+qi+12ϕi+12+qi+22ϕi+22+qi+32ϕi+32+2qiqi+1ϕiϕi+1+2qiqi+2ϕiϕi+2+2qiqi+3ϕiϕi+3+2qi+1qi+2ϕi+1ϕi+2+2qi+1qi+3ϕi+1ϕi+3+2qi+2qi+3ϕi+2ϕi+3]dV\displaystyle {} = \int_{\text{tet}} \Big[ q_i^2 \cdot \phi_i^2 + q_{i+1}^2 \cdot \phi_{i+1}^2 + q_{i+2}^2 \cdot \phi_{i+2}^2 + q_{i+3}^2 \cdot \phi_{i+3}^2 + 2 \cdot q_i \cdot q_{i+1} \cdot \phi_i \cdot \phi_{i+1} + 2 \cdot q_i \cdot q_{i+2} \cdot \phi_i \cdot \phi_{i+2} + 2 \cdot q_i \cdot q_{i+3} \cdot \phi_i \cdot \phi_{i+3} + 2 \cdot q_{i+1} \cdot q_{i+2} \cdot \phi_{i+1} \cdot \phi_{i+2} + 2 \cdot q_{i+1} \cdot q_{i+3} \cdot \phi_{i+1} \cdot \phi_{i+3} + 2 \cdot q_{i+2} \cdot q_{i+3} \cdot \phi_{i+2} \cdot \phi_{i+3} \Big] \,dV

Apply the formulas (6.31)

tetυ(i)(i+3)2dV\displaystyle \int_{\text{tet}} \upsilon_{(i)(i+3)}^2 \,dV =qi2Vtet10\displaystyle {} = q_i^2 \cdot \frac{\displaystyle V_{\text{tet}}}{\displaystyle 10}+qi+12Vtet10\displaystyle {} + q_{i+1}^2 \cdot \frac{\displaystyle V_{\text{tet}}}{\displaystyle 10}+qi+22Vtet10\displaystyle {} + q_{i+2}^2 \cdot \frac{\displaystyle V_{\text{tet}}}{\displaystyle 10}+qi+32Vtet10\displaystyle {} + q_{i+3}^2 \cdot \frac{\displaystyle V_{\text{tet}}}{\displaystyle 10}+2qiqi+1Vtet20\displaystyle {} + 2 \cdot q_i \cdot q_{i+1} \cdot \frac{\displaystyle V_{\text{tet}}}{\displaystyle 20}+2qiqi+2Vtet20\displaystyle {} + 2 \cdot q_i \cdot q_{i+2} \cdot \frac{\displaystyle V_{\text{tet}}}{\displaystyle 20}+2qiqi+3Vtet20\displaystyle {} + 2 \cdot q_i \cdot q_{i+3} \cdot \frac{\displaystyle V_{\text{tet}}}{\displaystyle 20}+2qi+1qi+2Vtet20\displaystyle {} + 2 \cdot q_{i+1} \cdot q_{i+2} \cdot \frac{\displaystyle V_{\text{tet}}}{\displaystyle 20}+2qi+1qi+3Vtet20\displaystyle {} + 2 \cdot q_{i+1} \cdot q_{i+3} \cdot \frac{\displaystyle V_{\text{tet}}}{\displaystyle 20}+2qi+2qi+3Vtet20\displaystyle {} + 2 \cdot q_{i+2} \cdot q_{i+3} \cdot \frac{\displaystyle V_{\text{tet}}}{\displaystyle 20}

Simplify the expression

tetυ(i)(i+3)2dV\displaystyle \int_{\text{tet}} \upsilon_{(i)(i+3)}^2 \,dV =Vtet20\displaystyle {} = \frac{\displaystyle V_{\text{tet}}}{\displaystyle 20}[2qi2+2qi+12+2qi+22+2qi+32+2qiqi+1+2qiqi+2+2qiqi+3+2qi+1qi+2+2qi+1qi+3+2qi+2qi+3]\displaystyle {} \cdot \Big[ 2 \cdot q_i^2 + 2 \cdot q_{i+1}^2 + 2 \cdot q_{i+2}^2 + 2 \cdot q_{i+3}^2 + 2 \cdot q_i \cdot q_{i+1} + 2 \cdot q_i \cdot q_{i+2} + 2 \cdot q_i \cdot q_{i+3} + 2 \cdot q_{i+1} \cdot q_{i+2} + 2 \cdot q_{i+1} \cdot q_{i+3} + 2 \cdot q_{i+2} \cdot q_{i+3} \Big]

We introduce the notation for the elements of the local damping matrix of the tetrahedron

c(i)(i)\displaystyle c_{(i)(i)} =Vtet10\displaystyle {} = \frac{\displaystyle V_{\text{tet}}}{\displaystyle 10}c(i+1)(i+1)\displaystyle c_{(i+1)(i+1)} =Vtet10\displaystyle {} = \frac{\displaystyle V_{\text{tet}}}{\displaystyle 10}c(i+2)(i+2)\displaystyle c_{(i+2)(i+2)} =Vtet10\displaystyle {} = \frac{\displaystyle V_{\text{tet}}}{\displaystyle 10}c(i+3)(i+3)\displaystyle c_{(i+3)(i+3)} =Vtet10\displaystyle {} = \frac{\displaystyle V_{\text{tet}}}{\displaystyle 10}c(m)(n)\displaystyle c_{(m)(n)} =c(n)(m)\displaystyle {} = c_{(n)(m)}=Vtet20,\displaystyle {} = \frac{\displaystyle V_{\text{tet}}}{\displaystyle 20},m\displaystyle mn,\displaystyle {} \neq n,m,n\displaystyle m, n{i,i\displaystyle {} \in \{i, i+1,i\displaystyle {} +1, i+2,i\displaystyle {} +2, i+3}\displaystyle {} +3\}
(6.32)

Thus, the local damping matrix for the tetrahedral element has the form

Ctet\displaystyle \mathbf{C}_{\text{tet}} =[c(i)(i)c(i)(i+1)c(i)(i+2)c(i)(i+3)c(i+1)(i)c(i+1)(i+1)c(i+1)(i+2)c(i+1)(i+3)c(i+2)(i)c(i+2)(i+1)c(i+2)(i+2)c(i+2)(i+3)c(i+3)(i)c(i+3)(i+1)c(i+3)(i+2)c(i+3)(i+3)]\displaystyle {} = \begin{bmatrix} c_{(i)(i)} & c_{(i)(i+1)} & c_{(i)(i+2)} & c_{(i)(i+3)}\\ c_{(i+1)(i)} & c_{(i+1)(i+1)} & c_{(i+1)(i+2)} & c_{(i+1)(i+3)}\\ c_{(i+2)(i)} & c_{(i+2)(i+1)} & c_{(i+2)(i+2)} & c_{(i+2)(i+3)}\\ c_{(i+3)(i)} & c_{(i+3)(i+1)} & c_{(i+3)(i+2)} & c_{(i+3)(i+3)} \end{bmatrix}=Vtet20\displaystyle {} = \frac{\displaystyle V_{\text{tet}}}{\displaystyle 20}[2111121111211112].\displaystyle \begin{bmatrix} 2 & 1 & 1 & 1\\ 1 & 2 & 1 & 1\\ 1 & 1 & 2 & 1\\ 1 & 1 & 1 & 2 \end{bmatrix}.
(6.33)

The local damping matrix is symmetric. The global damping matrix C\mathbf{C} is obtained by summing the contributions from all tetrahedral elements of the mesh using the assembly method: the elements of the local matrices are added to the corresponding elements of the global matrix according to the global node numbering. The dimension of the global damping matrix is N×NN \times N, where NN is the total number of mesh nodes.