====== 弹性力学求解体系 (Solving System of Elasticity) ======

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

===== 1. 基本架构：15个未知数与15个方程 =====

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

**未知数 (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个未知量。

^ 方程类型 ^ 个数 ^ 描述 ^
| **平衡方程** (Equilibrium) | 3 | 描述应力与体力之间的平衡关系 |
| **几何方程** (Geometric) | 6 | 描述应变与位移之间的微分关系 (柯西方程) |
| **本构方程** (Constitutive) | 6 | 描述应力与应变之间的物理关系 (广义胡克定律) |

===== 2. 边界条件 (Boundary Conditions) =====

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

=== 2.1 应力边界条件 (在 $S_\sigma$ 上) ===
当物体边界上的面力（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$) 是边界外法线的方向余弦。

=== 2.2 位移边界条件 (在 $S_u$ 上) ===
当物体边界上的位移已知时：
$$ u = \bar{u}, \quad v = \bar{v}, \quad w = \bar{w} $$

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

===== 3. 求解方法与推导 =====

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

==== 3.1 位移解法 (Displacement Method) ====

**核心思想：** 选取位移 $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.2 应力解法 (Stress Method) ====

**核心思想：** 选取应力为基本未知量。除了满足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个未知数，存在冗余）。


===== 4. 思考 =====

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

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

===== 艾里应力函数求解法 (Airy Stress Function Method) =====

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

==== 1. 引入动机与定义 ====

**动机：**
在二维问题中（不计体力），我们需要满足两个平衡方程：
$$
\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 $$

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


==== 2. 求解方程：双调和方程 ====

虽然平衡方程自动满足了，但应力分量还必须满足**变形协调方程**（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$ 的具体形式。

==== 3. 考虑体力的情况 ====

若物体受到有势体力（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$ 为泊松比)*

==== 4. 常用求解策略：逆解法 (Inverse Method) ====

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

  *   **思路：** 预先设定 $\phi$ 为某种形式的函数（通常是多项式），验证它是否满足双调和方程。
  *   **步骤：**
  -  选取一个满足 $\nabla^4 \phi = 0$ 的多项式。
  -  求出对应的应力分量 $\sigma_x, \sigma_y, \tau_{xy}$。
  -  检查这些应力分量对应什么样的边界上的面力分布。
  -  如果这正是我们需要解决的物理问题，则该解即为真解（根据唯一性定理）。

**常见的多项式尝试：**

^ 多项式次数 ^ $\phi(x,y)$ 形式 ^ 对应物理问题举例 ^
| **二次** | $A x^2 + B xy + C y^2$ | 均匀拉伸、纯剪切 |
| **三次** | $A x^3 + B x^2 y + C x y^2 + D y^3$ | 纯弯曲、受重力坝体 |
| **更高次** | 四次及以上 | 简支梁受均布载荷等复杂情况 |

==== 5. 物理意义总结 ====

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