跳转至

形函数的空间导数

形函数与有限元近似中引入的形函数 \(N_\alpha^e(\boldsymbol{r})\) 是单元局部坐标即自然坐标 \(\boldsymbol{r}\) 的函数。另一方面，通过应变-位移关系，弱形式的被积函数中会出现关于物理坐标（参考构形 \(\boldsymbol{X}\) 或当前构形 \(\boldsymbol{x}\)）的偏导数 \(\partial N_\alpha^e/\partial \boldsymbol{X}\) 或 \(\partial N_\alpha^e/\partial \boldsymbol{x}\)。本节整理从单元节点坐标和自然坐标导数计算这些形函数空间导数的步骤。

利用链式法则表示空间导数

在插值公式 \(\boldsymbol{X} = \sum_\alpha N_\alpha^e(\boldsymbol{r})\,\boldsymbol{X}^e_\alpha\)（上一节）下应用链式法则，并使用从自然坐标到物理坐标映射的Jacobian 矩阵 \(J_{ij} = \partial X_i/\partial r_j\)，可得

\[ \frac{\partial \boldsymbol{N}^e}{\partial \boldsymbol{X}} = \frac{\partial \boldsymbol{N}^e}{\partial \boldsymbol{r}}\,\left(\frac{\partial \boldsymbol{X}}{\partial \boldsymbol{r}}\right)^{-1} \]

其中，\(\partial \boldsymbol{N}^e/\partial \boldsymbol{X}\) 和 \(\partial \boldsymbol{N}^e/\partial \boldsymbol{r}\) 的第 \(\alpha\) 行表示节点 \(\alpha\) 的偏导数向量，二者均为 \(n_e \times 3\) 矩阵。

Jacobian 矩阵的构造

对插值公式关于 \(\boldsymbol{r}\) 求偏导，可得

\[ \frac{\partial \boldsymbol{X}}{\partial \boldsymbol{r}} = \left[\boldsymbol{X}^e_1, \ldots, \boldsymbol{X}^e_{n_e}\right]\frac{\partial \boldsymbol{N}^e}{\partial \boldsymbol{r}} \]

因此，可由单元节点坐标 \(\boldsymbol{X}^e_\alpha\) 和自然坐标导数 \(\partial \boldsymbol{N}^e/\partial \boldsymbol{r}\) 直接构造 Jacobian 矩阵。自然坐标导数由形函数的函数形式决定，并针对各单元类型预先实现。

Jacobian 矩阵的行列式 \(\det(\partial \boldsymbol{X}/\partial \boldsymbol{r})\) 用于数值积分中的体积微元变换 \(dV = \det(\partial \boldsymbol{X}/\partial \boldsymbol{r})\,d\boldsymbol{r}\)。当行列式为 \(0\) 时，单元发生退化，计算将失败。

计算流程

根据以上推导，空间导数 \(\partial \boldsymbol{N}^e/\partial \boldsymbol{X}\) 的计算可分为以下 4 个阶段。

  1. 计算自然坐标导数：以单元类型和评价点的自然坐标 \(\boldsymbol{r}\) 为输入，计算自然坐标导数矩阵 \(\partial \boldsymbol{N}^e/\partial \boldsymbol{r}\)。由于不同单元类型的形函数形式不同，因此按单元类型分支处理。
  2. 计算 Jacobian 矩阵：由单元节点坐标 \(\boldsymbol{X}^e\) 和 \(\partial \boldsymbol{N}^e/\partial \boldsymbol{r}\)，使用上一节的公式组装 Jacobian 矩阵 \(\partial \boldsymbol{X}/\partial \boldsymbol{r}\)。
  3. 计算逆矩阵和行列式：计算 Jacobian 矩阵的逆矩阵 \((\partial \boldsymbol{X}/\partial \boldsymbol{r})^{-1}\) 和行列式。行列式用于数值积分权重。
  4. 计算空间导数：由 \(\partial \boldsymbol{N}^e/\partial \boldsymbol{r}\) 与 \((\partial \boldsymbol{X}/\partial \boldsymbol{r})^{-1}\) 的乘积得到空间导数 \(\partial \boldsymbol{N}^e/\partial \boldsymbol{X}\)。

这 4 个阶段会根据单元类型和空间维数（二维、三维）分支，但整体流程相同。

参考构形与当前构形的通用处理

上述步骤只需将 \(\boldsymbol{X}^e\) 替换为当前构形的节点坐标 \(\boldsymbol{x}^e_\alpha = \boldsymbol{X}^e_\alpha + \boldsymbol{u}^e_\alpha\)，即可同样用于计算当前构形中的空间导数 \(\partial \boldsymbol{N}^e/\partial \boldsymbol{x}\)。因此，在 Total Lagrange 法和 Updated Lagrange 法中，只需切换输入节点坐标，即可共用同一计算步骤。

相关项目