soli8e¶

Purpose:

Compute element matrices for an 8 node isoparametric solid element.

Syntax:

Ke = soli8e(ex, ey, ez, ep, D)
[Ke, fe] = soli8e(ex, ey, ez, ep, D, eq)

Description:

soli8e provides an element stiffness matrix Ke and an element load vector fe for an 8 node isoparametric solid element.

The element nodal coordinates \(x_1, y_1, z_1, x_2\) etc. are supplied to the function by \(\mathbf{ex}\), \(\mathbf{ey}\) and \(\mathbf{ez}\), and the number of Gauss points \(n\) are supplied by \(\mathbf{ep}\).

\((n \times n \times n)\) integration points, \(n=1,2,3\)

The material properties are supplied by the constitutive matrix D. Any arbitrary \(\mathbf{D}\)-matrix with dimensions \((6 \times 6)\) may be given. For an isotropic elastic material the constitutive matrix can be formed by the function hooke, see Section Material functions.

\[\begin{split}\begin{aligned} \mathbf{ex} &= [\,x_1,\, x_2,\, \dots,\, x_8\,] \\ \mathbf{ey} &= [\,y_1,\, y_2,\, \dots,\, y_8\,] \\ \mathbf{ez} &= [\,z_1,\, z_2,\, \dots,\, z_8\,] \end{aligned}\end{split}\]

\[\mathbf{ep} = [\, n \,]\]

\[\begin{split}\mathbf{D} = \begin{bmatrix} D_{11} & D_{12} & \cdots & D_{16} \\ D_{21} & D_{22} & \cdots & D_{26} \\ \vdots & \vdots & \ddots & \vdots \\ D_{61} & D_{62} & \cdots & D_{66} \end{bmatrix}\end{split}\]

If different \(\mathbf{D}_i\) matrices are used in the Gauss points these \(\mathbf{D}_i\) matrices are stored in a global vector D. For numbering of the Gauss points, see eci in soli8s.

\[\begin{split}\mathbf{D} = \begin{bmatrix} \mathbf{D}_1 \\ \mathbf{D}_2 \\ \vdots \\ \mathbf{D}_{n^3} \end{bmatrix}\end{split}\]

If uniformly distributed loads are applied to the element, the element load vector fe is computed. The input variable

\[\begin{split}\mathrm{eq} = \begin{bmatrix} b_x \\ b_y \\ b_z \end{bmatrix}\end{split}\]

containing loads per unit volume, \(b_x\), \(b_y\), and \(b_z\), is then given.

Theory:

The element stiffness matrix \(\mathbf{K}^e\) and the element load vector \(\mathbf{f}_l^e\), stored in Ke and fe, respectively, are computed according to

\[ \begin{align}\begin{aligned}\mathbf{K}^e = \int_V \mathbf{B}^{eT} \mathbf{D} \mathbf{B}^e \, dV\\\mathbf{f}_l^e = \int_V \mathbf{N}^{eT} \mathbf{b} \, dV\end{aligned}\end{align} \]

with the constitutive matrix \(\mathbf{D}\) defined by D, and the body force vector \(\mathbf{b}\) defined by eq.

The evaluation of the integrals for the isoparametric 8 node solid element is based on a displacement approximation \(\mathbf{u}(\xi, \eta, \zeta)\), expressed in a local coordinate system in terms of the nodal variables \(u_1, u_2, \dots, u_{24}\) as

\[ \begin{align}\begin{aligned}\mathbf{u}(\xi, \eta, \zeta) = \mathbf{N}^e \mathbf{a}^e\\\begin{split}\mathbf{u} = \begin{bmatrix} u_x \\ u_y \\ u_z \end{bmatrix} \quad \mathbf{N}^e = \begin{bmatrix} N^e_1 & 0 & 0 & N^e_2 & 0 & 0 & \dots & N^e_8 & 0 & 0 \\ 0 & N^e_1 & 0 & 0 & N^e_2 & 0 & \dots & 0 & N^e_8 & 0 \\ 0 & 0 & N^e_1 & 0 & 0 & N^e_2 & \dots & 0 & 0 & N^e_8 \\ \end{bmatrix} \quad \mathbf{a}^e = \begin{bmatrix} u_1 \\ u_2 \\ \vdots \\ u_{24} \end{bmatrix}\end{split}\end{aligned}\end{align} \]

The element shape functions are given by

\[\begin{split}\begin{aligned} N_1^e &= \frac{1}{8}(1-\xi)(1-\eta)(1-\zeta) & N_5^e &= \frac{1}{8}(1-\xi)(1-\eta)(1+\zeta) \\ N_2^e &= \frac{1}{8}(1+\xi)(1-\eta)(1-\zeta) & N_6^e &= \frac{1}{8}(1+\xi)(1-\eta)(1+\zeta) \\ N_3^e &= \frac{1}{8}(1+\xi)(1+\eta)(1-\zeta) & N_7^e &= \frac{1}{8}(1+\xi)(1+\eta)(1+\zeta) \\ N_4^e &= \frac{1}{8}(1-\xi)(1+\eta)(1-\zeta) & N_8^e &= \frac{1}{8}(1-\xi)(1+\eta)(1+\zeta) \end{aligned}\end{split}\]

The \(\mathbf{B}^e\)-matrix is obtained as

\[\mathbf{B}^e = \tilde{\nabla} \mathbf{N}^e\]

where

\[\begin{split}\tilde{\nabla} = \begin{bmatrix} \frac{\partial}{\partial x} & 0 & 0 \\ 0 & \frac{\partial}{\partial y} & 0 \\ 0 & 0 & \frac{\partial}{\partial z} \\ \frac{\partial}{\partial y} & \frac{\partial}{\partial x} & 0 \\ \frac{\partial}{\partial z} & 0 & \frac{\partial}{\partial x} \\ 0 & \frac{\partial}{\partial z} & \frac{\partial}{\partial y} \end{bmatrix} \qquad \begin{bmatrix} \frac{\partial}{\partial x} \\ \frac{\partial}{\partial y} \\ \frac{\partial}{\partial z} \end{bmatrix} = (\mathbf{J}^T)^{-1} \begin{bmatrix} \frac{\partial}{\partial \xi} \\ \frac{\partial}{\partial \eta} \\ \frac{\partial}{\partial \zeta} \end{bmatrix}\end{split}\]

\[\begin{split}\mathbf{J} = \begin{bmatrix} \frac{\partial x}{\partial \xi} & \frac{\partial x}{\partial \eta} & \frac{\partial x}{\partial \zeta} \\ \frac{\partial y}{\partial \xi} & \frac{\partial y}{\partial \eta} & \frac{\partial y}{\partial \zeta} \\ \frac{\partial z}{\partial \xi} & \frac{\partial z}{\partial \eta} & \frac{\partial z}{\partial \zeta} \end{bmatrix}\end{split}\]

Evaluation of the integrals is done by Gauss integration.