弹性力学的核心任务是在给定的边界条件下，求解弹性体内部的应力、应变和位移场。这是一个典型的偏微分方程组边值问题。

弹性力学问题的求解由三大类方程支撑。

未知数 (15个):

• 应力分量 (6个): $\sigma_x, \sigma_y, \sigma_z, \tau_{xy}, \tau_{yz}, \tau_{zx}$
• 应变分量 (6个): $\varepsilon_x, \varepsilon_y, \varepsilon_z, \gamma_{xy}, \gamma_{yz}, \gamma_{zx}$
• 位移分量 (3个): $u, v, w$

基本方程 (15个): 原则上，联立以下15个方程可求解上述15个未知量。

偏微分方程的通解包含任意函数，必须引入边界条件才能确定具体问题的特解。

当物体边界上的面力（Surface Force）已知时： $$ \begin{cases} P_x = \sigma_x l + \tau_{xy} m + \tau_{xz} n \\ P_y = \tau_{yx} l + \sigma_y m + \tau_{yz} n \\ P_z = \tau_{zx} l + \tau_{zy} m + \sigma_z n \end{cases} $$ 简写为张量形式：$P_i = \sigma_{ij} n_j$。 其中 $l, m, n$ (或 $n_x, n_y, n_z$) 是边界外法线的方向余弦。

当物体边界上的位移已知时： $$ u = \bar{u}, \quad v = \bar{v}, \quad w = \bar{w} $$

• 第一类边值问题: 给定体力，且边界上全为面力已知（应力边界条件）。
• 第二类边值问题: 给定体力，且边界上全为位移已知（位移边界条件）。
• 第三类边值问题 (混合): 边界一部分给定面力，另一部分给定位移。

由于直接联立15个方程求解过于繁琐，实际中通常采用消元法，将方程组化简为以位移或应力为主的偏微分方程组。

核心思想： 选取位移 $u, v, w$ 为基本未知量。将本构方程代入平衡方程，利用几何方程消去应变和应力。

推导过程 (Lamé-Navier 方程):

1. 从平衡方程出发：

$$ \frac{\partial \sigma_x}{\partial x} + \frac{\partial \tau_{yx}}{\partial y} + \frac{\partial \tau_{zx}}{\partial z} + f_x = 0 $$

2. 利用本构方程（应力-应变关系）将应力用应变表示：

引入拉梅常数 $\lambda, \mu$ ($\mu = G$)：

$$ \sigma_x = \lambda \theta + 2\mu \varepsilon_x, \quad \tau_{xy} = \mu \gamma_{xy} $$

代入平衡方程。

3. 利用几何方程将应变用位移表示： $$ \varepsilon_x = \frac{\partial u}{\partial x}, \quad \theta = \frac{\partial u}{\partial x} + \frac{\partial v}{\partial y} + \frac{\partial w}{\partial z} = \nabla \cdot \mathbf{u} $$

$$ \gamma_{xy} = \frac{\partial v}{\partial x} + \frac{\partial u}{\partial y} $$

4. 代入整理 (以 $x$ 方向为例)：

$$ \frac{\partial}{\partial x}(\lambda \theta + 2\mu \frac{\partial u}{\partial x}) + \frac{\partial}{\partial y}[\mu(\frac{\partial v}{\partial x} + \frac{\partial u}{\partial y})] + \frac{\partial}{\partial z}[\mu(\frac{\partial w}{\partial x} + \frac{\partial u}{\partial z})] + f_x = 0 $$

整理各项：

$$ \lambda \frac{\partial \theta}{\partial x} + \mu \frac{\partial}{\partial x}(\frac{\partial u}{\partial x} + \frac{\partial v}{\partial y} + \frac{\partial w}{\partial z}) + \mu (\frac{\partial^2 u}{\partial x^2} + \frac{\partial^2 u}{\partial y^2} + \frac{\partial^2 u}{\partial z^2}) + f_x = 0 $$

注意到中间项即为 $\mu \frac{\partial \theta}{\partial x}$，且 $\nabla^2 = \frac{\partial^2}{\partial x^2} + \frac{\partial^2}{\partial y^2} + \frac{\partial^2}{\partial z^2}$。

最终结果 (Lamé-Navier 方程):

$$ \begin{cases} (\lambda + \mu) \frac{\partial \theta}{\partial x} + \mu \nabla^2 u + f_x = 0 \\ (\lambda + \mu) \frac{\partial \theta}{\partial y} + \mu \nabla^2 v + f_y = 0 \\ (\lambda + \mu) \frac{\partial \theta}{\partial z} + \mu \nabla^2 w + f_z = 0 \end{cases} $$

矢量形式：$(\lambda + \mu) \nabla (\nabla \cdot \mathbf{u}) + \mu \nabla^2 \mathbf{u} + \mathbf{f} = 0$

位移法将15个方程化简为3个关于位移的二阶偏微分方程。适用于第二类边值问题。若为第一类问题，需将应力边界条件也转化为位移形式。

核心思想： 选取应力为基本未知量。除了满足3个平衡方程外，还需满足变形协调条件（即6个相容方程用应力表示）。

推导过程 (Beltrami-Michell 方程):

1. 从圣维南应变协调方程出发 (以 $xy$ 面为例)：

$$ \frac{\partial^2 \varepsilon_x}{\partial y^2} + \frac{\partial^2 \varepsilon_y}{\partial x^2} = \frac{\partial^2 \gamma_{xy}}{\partial x \partial y} $$

2. 代入胡克定律 (应变用应力表示)：

$$ \varepsilon_x = \frac{1}{E}[\sigma_x - \nu(\sigma_y + \sigma_z)], \quad \gamma_{xy} = \frac{2(1+\nu)}{E}\tau_{xy} $$

代入上式，并利用 $\sigma_z = \Theta - \sigma_x - \sigma_y$ ($\Theta$ 为第一应力不变量)，经过繁琐的代数运算（详见手写推导图），可得：

$$ \frac{\partial^2}{\partial y^2}[(1+\nu)\sigma_x - \nu\Theta] + \frac{\partial^2}{\partial x^2}[(1+\nu)\sigma_y - \nu\Theta] = 2(1+\nu)\frac{\partial^2 \tau_{xy}}{\partial x \partial y} $$

3. 利用平衡方程消除切应力项：

由平衡方程 $\frac{\partial \tau_{xy}}{\partial x} = -f_y - \frac{\partial \sigma_y}{\partial y} - \frac{\partial \tau_{yz}}{\partial z}$ 等关系，对平衡方程求导并代入协调方程，消除混合偏导数项。

4. 引入不变量 $\Theta$ 的拉普拉斯关系：

通过将三个正应力的协调方程相加，可推导出：

$$ \nabla^2 \Theta = -\frac{1+\nu}{1-\nu} \nabla \cdot \mathbf{f} = -\frac{1+\nu}{1-\nu} (\frac{\partial f_x}{\partial x} + \frac{\partial f_y}{\partial y} + \frac{\partial f_z}{\partial z}) $$

最终结果 (Beltrami-Michell 协调方程):

在体力为常数或零的情况下，方程简化为：

$$ \nabla^2 \sigma_{ij} + \frac{1}{1+\nu} \frac{\partial^2 \Theta}{\partial x_i \partial x_j} = 0 $$

其中 $i, j = x, y, z$。

若考虑一般体力，完整形式为：

$$ \nabla^2 \sigma_{x} + \frac{1}{1+\nu} \frac{\partial^2 \Theta}{\partial x^2} = -\frac{\nu}{1-\nu} \nabla \cdot \mathbf{f} - 2\frac{\partial f_x}{\partial x} $$

(其他分量类似)

应力法求解时，必须同时满足：

1. 3个平衡方程 (在域内)

2. 6个 Beltrami-Michell 协调方程 (在域内)

3. 应力边界条件 (在边界上)

这使得直接用应力法求解三维问题非常困难（9个方程解6个未知数，存在冗余）。

有没有比较方便的寻找满足应力法方程应力分布的方法？

应力法虽然方程多，但在二维问题（平面应力/平面应变）中，可以通过引入艾里应力函数 (Airy Stress Function) $\phi(x,y)$ 来大大简化。通过将应力分量表示为 $\phi$ 的偏导数，可以自动满足平衡方程，此时只需解一个关于 $\phi$ 的双调和方程 ($\nabla^4 \phi = 0$) 即可。

在求解二维弹性力学问题（平面应力或平面应变问题）时，直接使用应力法求解仍然面临方程多于未知数的情况。为了简化求解过程，G.B. Airy 提出了一种通过引入一个标量函数将方程组降维的方法。

动机： 在二维问题中（不计体力），我们需要满足两个平衡方程： $$ \frac{\partial \sigma_x}{\partial x} + \frac{\partial \tau_{xy}}{\partial y} = 0, \quad \frac{\partial \tau_{yx}}{\partial x} + \frac{\partial \sigma_y}{\partial y} = 0 $$ 我们希望找到一种表示应力的方式，使其自动满足上述平衡方程，从而将注意力集中在相容方程上。

定义： 引入一个关于 $x, y$ 的函数 $\phi(x, y)$，称为艾里应力函数，定义如下（不计体力）：

$$ \sigma_x = \frac{\partial^2 \phi}{\partial y^2}, \quad \sigma_y = \frac{\partial^2 \phi}{\partial x^2}, \quad \tau_{xy} = -\frac{\partial^2 \phi}{\partial x \partial y} $$

验证平衡方程：

将上述定义代入第一个平衡方程：

$$ \frac{\partial}{\partial x}\left(\frac{\partial^2 \phi}{\partial y^2}\right) + \frac{\partial}{\partial y}\left(-\frac{\partial^2 \phi}{\partial x \partial y}\right) = \frac{\partial^3 \phi}{\partial x \partial y^2} - \frac{\partial^3 \phi}{\partial x \partial y^2} \equiv 0 $$

同理可验证第二个平衡方程也恒成立。

虽然平衡方程自动满足了，但应力分量还必须满足变形协调方程（Compatibility Equation）。

在二维问题中，Beltrami-Michell 协调方程简化为：

$$ \nabla^2 (\sigma_x + \sigma_y) = 0 $$

• (注：此式适用于不计体力且物体为单连通域的情况)*

推导过程：

1. 计算应力主不变量（第一不变量）：

$$ \sigma_x + \sigma_y = \frac{\partial^2 \phi}{\partial y^2} + \frac{\partial^2 \phi}{\partial x^2} = \left( \frac{\partial^2}{\partial x^2} + \frac{\partial^2}{\partial y^2} \right) \phi = \nabla^2 \phi $$

2. 代入协调方程：

$$ \nabla^2 (\nabla^2 \phi) = 0 $$

最终控制方程：

即双调和方程 (Biharmonic Equation)：

$$ \nabla^4 \phi = \frac{\partial^4 \phi}{\partial x^4} + 2\frac{\partial^4 \phi}{\partial x^2 \partial y^2} + \frac{\partial^4 \phi}{\partial y^4} = 0 $$

这意味着，只要找到一个函数 $\phi(x, y)$ 满足双调和方程，由它求导得出的应力分量就既满足平衡条件，又满足变形协调条件。此时，问题转化为根据边界条件确定 $\phi$ 的具体形式。

若物体受到有势体力（Conservative Body Forces）作用，设体力势为 $V$（即 $f_x = -\frac{\partial V}{\partial x}, f_y = -\frac{\partial V}{\partial y}$），则艾里应力函数的定义修正为：

$$ \sigma_x = \frac{\partial^2 \phi}{\partial y^2} + V, \quad \sigma_y = \frac{\partial^2 \phi}{\partial x^2} + V, \quad \tau_{xy} = -\frac{\partial^2 \phi}{\partial x \partial y} $$

此时的控制方程变为： $$ \nabla^4 \phi = -(1-\nu)\nabla^2 V $$ *(其中 $\nu$ 为泊松比)*

由于直接解四阶偏微分方程 $\nabla^4 \phi = 0$ 非常困难，弹性力学中常采用逆解法。

• 思路： 预先设定 $\phi$ 为某种形式的函数（通常是多项式），验证它是否满足双调和方程。
• 步骤：

1. 选取一个满足 $\nabla^4 \phi = 0$ 的多项式。
2. 求出对应的应力分量 $\sigma_x, \sigma_y, \tau_{xy}$。
3. 检查这些应力分量对应什么样的边界上的面力分布。
4. 如果这正是我们需要解决的物理问题，则该解即为真解（根据唯一性定理）。

常见的多项式尝试：

• 统一性： 在不计体力且单连通域的情况下，平面应力问题和平面应变问题的控制方程都是 $\nabla^4 \phi = 0$。这意味着对于应力边界值问题，两种状态下的应力分布是相同的（与材料常数 $E, \nu$ 无关）。
• 降维打击： 将含有3个未知函数（$\sigma_x, \sigma_y, \tau_{xy}$）的方程组，转化为求解1个标量函数 $\phi$ 的问题。