flw2teΒΆ
- Purpose:
Compute element stiffness matrix for a triangular heat flow element.
- Syntax:
Ke = flw2te(ex, ey, ep, D) [Ke, fe] = flw2te(ex, ey, ep, D, eq)- Description:
flw2teprovides the element stiffness (conductivity) matrixKeand the element load vectorfefor a triangular heat flow element.The element nodal coordinates \(x_1\), \(y_1\), \(x_2\) etc, are supplied to the function by
exandey, the element thickness \(t\) is supplied byepand the thermal conductivities (or corresponding quantities) \(k_{xx}\), \(k_{xy}\) etc are supplied byD.\[\begin{split}\begin{array}{l} \mathbf{ex} = [\, x_1 \;\; x_2 \;\; x_3\,] \\ \mathbf{ey} = [\, y_1 \;\; y_2 \;\; y_3\,] \end{array} \qquad \mathbf{ep} = [\, t \,] \qquad \mathbf{D} = \begin{bmatrix} k_{xx} & k_{xy} \\ k_{yx} & k_{yy} \end{bmatrix}\end{split}\]If the scalar variable
eqis given in the function, the element load vector \(\mathbf{fe}\) is computed, using\[\mathbf{eq} = [\, Q \,]\]where \(Q\) is the heat supply per unit volume.
- Theory:
The element stiffness matrix \(\mathbf{K}^e\) and the element load vector \(\mathbf{f}_l^e\), stored in
Keandfe, respectively, are computed according to\[\mathbf{K}^e = (\mathbf{C}^{-1})^T \int_A \bar{\mathbf{B}}^T \mathbf{D} \bar{\mathbf{B}}\, t\, dA\, \mathbf{C}^{-1}\]\[\mathbf{f}_l^e = (\mathbf{C}^{-1})^T \int_A \bar{\mathbf{N}}^T Q\, t\, dA\]with the constitutive matrix \(\mathbf{D}\) defined by
D.The evaluation of the integrals for the triangular element is based on the linear temperature approximation \(T(x, y)\) and is expressed in terms of the nodal variables \(T_1\), \(T_2\) and \(T_3\) as
\[T(x, y) = \mathbf{N}^e \mathbf{a}^e = \bar{\mathbf{N}}\, \mathbf{C}^{-1} \mathbf{a}^e\]where
\[\begin{split}\bar{\mathbf{N}} = [\, 1 \;\; x \;\; y\,] \qquad \mathbf{C} = \begin{bmatrix} 1 & x_1 & y_1 \\ 1 & x_2 & y_2 \\ 1 & x_3 & y_3 \end{bmatrix} \qquad \mathbf{a}^e = \begin{bmatrix} T_1 \\ T_2 \\ T_3 \end{bmatrix}\end{split}\]and hence it follows that
\[\begin{split}\bar{\mathbf{B}} = \nabla \bar{\mathbf{N}} = \begin{bmatrix} 0 & 1 & 0 \\ 0 & 0 & 1 \end{bmatrix} \qquad \nabla = \begin{bmatrix} \dfrac{\partial}{\partial x} \\ \dfrac{\partial}{\partial y} \end{bmatrix}\end{split}\]Evaluation of the integrals for the triangular element yields
\[\mathbf{K}^e = (\mathbf{C}^{-1})^T \bar{\mathbf{B}}^T \mathbf{D} \bar{\mathbf{B}} \mathbf{C}^{-1} t A\]\[\begin{split}\mathbf{f}_l^e = \frac{Q A t}{3} \begin{bmatrix} 1 \\ 1 \\ 1 \end{bmatrix}\end{split}\]where the element area \(A\) is determined as
\[A = \frac{1}{2} \det \mathbf{C}\]