solidDisplacementFoam解析

引言

大部分的固体应力分析都使用的是有限元。然而使用有限体积法计算应力应变以及位移也是可行的。同时有限体积法具有一些优良的优点：如守恒特性。在将力离散到网格单元面的时候，所有的力均被抵消仅仅在边界上存在受力。这与与流体的质量守恒类似。OpenFOAM中的solidDisplacementFoam则采用有限体积法来对位移矢量方程进行离散，进而获得应力应变张量。本文对其简述。

控制方程与算法

在不考虑体积力的情况下，solidDisplacementFoam中求解的为下述位移方程：

(1)\[ \frac{\p ^2 \bfD}{\p t ^2}-\nabla\cdot\left( \mu(\nabla \bfD+\nabla \bfD^T)+\lambda (\nabla\cdot\bfD)\bfI \right)=0 \]

(2)\[ \mu=\frac{E}{2(1 + \nu)} \]

(3)\[ \lambda=\frac{\nu E}{(1 + \nu)(1 - 2\nu)} \]

其中\(\bfD\)表示位移矢量，\(E\)表示杨氏模量，\(\lambda\)表示拉梅系数的第1个参数，\(\mu\)表示拉梅系数的第2个参数（在流体力学中理解为粘度，在固体力学中理解为剪切模量）。solidDisplacementFoam求解的为位移矢量。因此，在给定边界条件的时候，给定的是位移矢量的边界条件。例如：固定位移大小边界条件之类。在方程(1)求解获得\(\bfD\)之后，可以获得应变张量\(\boldsymbol{\varepsilon}\)：

(4)\[ \boldsymbol{\varepsilon}=\frac{1}{2}\left(\nabla\bfU+\nabla\bfU^T\right) \]

进一步的获取应变张量：

(5)\[ \boldsymbol{\sigma}=2\mu\boldsymbol{\varepsilon}+\lambda (\nabla\cdot\bfD)\bfI \]

仔细观察方程(1)，其实际上表示3个位移矢量的分量的方程。每个位移矢量分量的方程里面，涉及到其他的位移矢量分量。传统的FEM算法，是采用耦合求解的方法，离散成一个大的方程组同时求解。在CFD里面也有类似的求解算法，即为耦合算法。

本文讨论分离式求解，即将三个位移矢量分量的方程离散为三个分量方程。这与速度矢量的方程离散是一样的。采用这种离散方式可以得到3个稀疏矩阵，稀疏矩阵的好处之一是可以降低内存存储，另外是可以调用迭代矩阵求解器。然而分离式求解一些变量需要做显性离散处理，因此分离式求解不仅需要用稀疏线性系统求解器来迭代求解矩阵，还需要迭代求解位移矢量分量在三个方程之间的耦合关系。

具体的，在求解方程(1)的时候，其中显性与隐性离散的处理关系可以写为：

(6)\[\underset{\mathrm{implicit}}{\underbrace{\frac{\p ^2 \bfD}{\p t ^2}-\nabla\cdot ( \mu(\nabla \bfD)}} \underset{\mathrm{explicit}}{\underbrace{-\nabla\cdot ( \mu(\nabla \bfD^T)-\nabla\cdot\left(\lambda (\nabla\cdot\bfD)\bfI \right)}} =0 \]

然而在实际过程中，方程(6)的离散过程仅仅会部分收敛。这是因为显性项的贡献比隐性项的贡献还要大。因此，OpenFOAM将其改写为：

(7)\[\begin{split} \begin{multline} \underset{\mathrm{implicit}}{\underbrace{\frac{\p ^2 \bfD}{\p t ^2}-\nabla\cdot ( (2\mu+\lambda)\nabla \bfD)}} \\ \underset{\mathrm{explicit}}{\underbrace{-\nabla\cdot ( \mu(\nabla \bfD^T)-\nabla\cdot\left(\lambda (\nabla\cdot\bfD)\bfI \right) -\nabla\cdot((\mu+\lambda)\nabla\bfD)}} =0 \end{multline} \end{split}\]