结构力学基础 — 应力、应变、有限元法与非线性分析完整解析
老师，我在学FEA软件，发现里面有好多概念——应力、应变、刚度矩阵……感觉很混乱，应该从哪里开始理解？
哈，这个问题太常见了。简单来说，结构力学就是研究"物体受力后会发生什么"。核心就两件事：应力（内部的力）和应变（形变程度）。搞清楚这两个，后面所有东西都会清晰很多。
1. 应力与应变——结构分析的两大主角
应力和应变……字面意思我懂，但感觉在软件里看到的值不知道代表什么意义。
好，举个例子。你把一根钢棒两端拉开，截面上单位面积承受的力就是正应力。比如截面积1cm²的棒，拉力1000N，应力就是10 MPa。应变呢，就是棒拉长了多少比例——原来100mm，拉长后变101mm，应变就是1%。
1.1 正应力与剪应力
应力 $\sigma$ 的基本定义：
其中 $F$ 是垂直于截面的法向力（单位：N），$A$ 是截面积（单位：m²），$\sigma$ 的单位是 Pa（帕斯卡）= N/m²。工程中常用 MPa = 10⁶ Pa。
除了拉压的正应力，截面上还有平行于截面的剪应力 $\tau$：
这里 $V$ 是切向力（剪力）。
1.2 正应变与剪应变
正应变（线应变）$\varepsilon$ 描述长度变化：
剪应变 $\gamma$ 描述角度变化（弧度）：
应变是无量纲量（m/m），实际工程中常用 $\mu\varepsilon$（微应变，10⁻⁶）或百分比表示。
1.3 三维应力状态——应力张量
三维空间中，一个点的完整应力状态需要用 应力张量（3×3矩阵）表示：
由力矩平衡条件，应力张量是对称的：$\tau_{xy} = \tau_{yx}$，$\tau_{xz} = \tau_{zx}$，$\tau_{yz} = \tau_{zy}$。因此独立分量只有6个，通常写成向量形式：
哦，这就是为什么软件里有σx, σy, σz, τxy这些分量！那通常我们关心的是哪个应力值呢？
实际工程里最常用的是冯米塞斯等效应力（Von Mises），它把三维应力状态综合成一个数，方便跟材料屈服强度比较。汽车车身分析、压力容器设计，基本都是看这个值来判断是否安全。
2. 胡克定律与弹性常数
那应力和应变之间是什么关系？我在材料设置里输入了杨氏模量E，这是什么？
这就是著名的胡克定律！简单说就是：弹性范围内，应力和应变成正比，比例系数就是杨氏模量 E，单位是 GPa。钢的 E 约210 GPa，铝约70 GPa，橡胶约0.001 GPa——差了好几个数量级。
2.1 单轴胡克定律
E 越大，材料越"硬"，同样应力下变形越小。
2.2 泊松效应与泊松比 ν
拉伸一根棒时，除了轴向伸长，横截面也会收缩。这个现象叫泊松效应，比例系数是泊松比 $\nu$：
金属材料 $\nu$ 约为0.25～0.35，不可压缩材料（如橡胶）$\nu \approx 0.5$。$\nu$ 必须满足 $-1 < \nu < 0.5$（热力学稳定性要求）。
2.3 剪切模量 G
G 描述材料抵抗剪切变形的能力。对于各向同性材料，$E$、$\nu$、$G$ 三者只有两个是独立的，因为它们通过上式相关联。
2.4 三维广义胡克定律
三维情况下，应力向量与应变向量通过弹性刚度矩阵 $\mathbf{D}$（6×6）关联：
对于各向同性材料：
原来如此。那软件里还有个"体积模量"K，它和E有什么关系？
体积模量 $K = E / [3(1-2\nu)]$，描述材料在均匀压力下的体积压缩能力。当 $\nu \to 0.5$ 时 $K \to \infty$，意味着材料不可压缩。这在做橡胶或生物软组织仿真时很重要——$\nu$ 取0.499而不是0.5，就是为了避免数值上"锁死"（体积锁定问题）。
| 材料 | 杨氏模量 E (GPa) | 泊松比 ν | 密度 (kg/m³) | 屈服强度 (MPa) |
|---|---|---|---|---|
| 结构钢 | 200–210 | 0.27–0.30 | 7,800 | 250–350 |
| 铝合金（6061） | 68–70 | 0.33 | 2,700 | 270 |
| 钛合金（Ti-6Al-4V） | 110–115 | 0.34 | 4,430 | 880 |
| 碳纤维复合材料（纵向） | 130–200 | 0.25–0.35 | 1,600 | —（各向异性） |
| 混凝土（压） | 20–40 | 0.15–0.20 | 2,400 | —（张拉强度低） |
3. 多轴应力状态与屈服准则
我理解了应力状态，但怎么判断材料要屈服了？材料手册给的是单轴拉伸强度，但实际是三维应力状态啊。
对，这就是屈服准则要解决的问题——把三维应力状态"折算"成一个等效应力，再跟单轴屈服强度比较。金属材料最常用的是冯米塞斯准则，它基于应变能理论；另一个是特雷斯卡准则，基于最大剪应力。
3.1 冯米塞斯屈服准则（Von Mises）
冯米塞斯等效应力（也叫 $\sigma_{vM}$ 或 $\sigma_{eq}$）定义为：
屈服条件：$\sigma_{vM} \geq \sigma_y$（单轴屈服强度）
冯米塞斯准则对延性金属（钢、铝等）非常准确，是 FEA 软件默认的屈服准则。
3.2 特雷斯卡屈服准则（Tresca）
基于最大剪应力理论：
其中 $\sigma_1 \geq \sigma_2 \geq \sigma_3$ 是主应力。特雷斯卡准则比冯米塞斯更保守（偏安全），在压力容器设计（ASME规范）中有应用。
3.3 莫尔应力圆
二维应力状态下，通过莫尔圆可以图形化求出任意截面上的应力，以及主应力和最大剪应力：
那什么时候用冯米塞斯，什么时候用特雷斯卡呢？
一般结构分析、汽车零件、航空结构件——都用冯米塞斯，因为实验验证更准。压力容器、核反应堆压力边界这类安全级别高的——用特雷斯卡，因为它更保守，留的安全余量更多。脆性材料（铸铁、陶瓷）两个都不适用，要用莫尔-库仑准则或德鲁克-普拉格准则。
4. 梁理论——欧拉-伯努利与铁木辛柯
我做梁分析的时候，选单元类型有时候有两个选项，一个是"欧拉梁"一个是"铁木辛柯梁"，有什么区别？
关键区别是：欧拉-伯努利梁忽略剪切变形，假设截面始终垂直于梁轴线；铁木辛柯梁考虑剪切变形，截面可以转一个角度。细长梁（跨高比L/h > 10）两者差别小，但短梁或夹层结构一定要用铁木辛柯，不然误差很大。
4.1 欧拉-伯努利梁方程
挠度 $w(x)$ 满足：
其中 $I$ 是截面二阶矩（惯性矩），$q(x)$ 是分布载荷。简支梁在均布载荷 $q_0$ 下的最大挠度：
4.2 弯曲应力公式
$M$ 是弯矩，$y$ 是到中性轴的距离。最大弯曲应力发生在截面顶部和底部（$y = \pm h/2$）：
4.3 铁木辛柯梁的修正
铁木辛柯梁引入独立的截面转角 $\psi$，挠度方程变为两个耦合方程：
其中 $\kappa$ 是剪切修正系数（矩形截面约0.833，圆截面约0.886），$G$ 是剪切模量，$A$ 是截面积。
| 特性 | 欧拉-伯努利梁 | 铁木辛柯梁 |
|---|---|---|
| 剪切变形 | 忽略 | 考虑 |
| 适用条件 | 细长梁（L/h > 10） | 短梁、深梁、夹层结构 |
| 自由度/节点 | 2（w, θ） | 2（w, θ，独立） |
| 数值问题 | 无锁死 | 需注意剪切锁死（需缩减积分） |
| 计算复杂度 | 低 | 略高 |
5. 有限元法（FEM）数学基础
FEA软件里核心就是有限元法，我知道最终是解一个 Ku=f 的方程组，但这个K从哪里来？为什么网格细了结果就更精确？
好问题！有限元法的核心思想是：把连续体（梁、壳、体）分成有限个小单元，在每个单元里用简单的多项式函数近似位移场，再通过能量原理组装出整体方程组。网格越细，近似越接近精确解，这就是收敛性。
5.1 弱形式与伽辽金法
平衡方程（强形式）$\nabla \cdot \boldsymbol{\sigma} + \mathbf{b} = 0$ 需要解精确满足，在实际中很困难。有限元法的出发点是弱形式（虚功原理）：
左边是内虚功（应变能），右边是外力（体力 $\mathbf{b}$ + 面力 $\mathbf{t}$）做的虚功。
5.2 形函数（插值函数）
在每个单元内，位移场用节点位移 $\mathbf{u}_e$ 通过形函数矩阵 $\mathbf{N}$ 插值：
应变由位移的梯度计算，通过B矩阵（应变-位移矩阵）：
线性四面体单元（TET4）：每节点3个位移自由度，形函数是线性的；二次四面体单元（TET10）：10节点，形函数是二次的，精度高得多。
5.3 单元刚度矩阵的组装
将弱形式代入并利用胡克定律 $\boldsymbol{\sigma} = \mathbf{D}\boldsymbol{\varepsilon}$，得到单元刚度矩阵 $\mathbf{K}_e$：
将所有单元的 $\mathbf{K}_e$ 按节点自由度编号"拼装"成全局刚度矩阵 $\mathbf{K}$，载荷向量 $\mathbf{f}$ 同理组装，最终得到：
5.4 高斯积分
单元刚度矩阵的积分通常用高斯数值积分计算。在自然坐标 $(-1,1)$ 上，$n$ 点高斯积分对最高次数 $2n-1$ 的多项式精确：
三维单元（六面体）用 $2\times2\times2 = 8$ 个积分点；四面体用4个（对应二次单元）。积分点数不足会导致"欠积分"，过多会增加计算量。
那边界条件是怎么加进去的？约束是固定某些位移为零吗？
对。本质边界条件（Dirichlet条件，如固定支撑 $u=0$）通过修改刚度矩阵或用罚函数法施加。自然边界条件（Neumann条件，如面力 $\mathbf{t}$）已经自动包含在弱形式的右端项里了。没有施加足够约束的话，$\mathbf{K}$ 是奇异的，方程组无解——这就是"刚体运动未约束"错误的来源。
5.5 常见单元类型与选择
| 单元类型 | 节点数 | 精度 | 适用场景 | 备注 |
|---|---|---|---|---|
| TET4（线性四面体） | 4 | 低 | 快速网格划分 | 过刚，避免用于精度要求高的分析 |
| TET10（二次四面体） | 10 | 高 | 复杂几何体，主流选择 | 自动网格划分首选 |
| HEX8（线性六面体） | 8 | 中 | 规则几何体 | 需注意剪切锁死（缩减积分） |
| HEX20（二次六面体） | 20 | 最高 | 高精度分析 | 手动网格成本高 |
| SHELL（壳单元） | 4/8 | 高 | 薄板、薄壳结构 | 需正确定义厚度和积分点 |
| BEAM（梁单元） | 2/3 | 高 | 细长杆、框架结构 | 需选择欧拉梁或铁木辛柯梁 |
6. 非线性分析概述
什么时候必须做非线性分析？我总是用线性分析，会不会有问题？
线性分析有三个假设：①小变形、②线弹性材料、③无接触。任何一个不满足就要用非线性。比如汽车碰撞——大变形+材料塑性+接触，三个非线性都有。板材冲压成形，应变超过20%，必须用几何+材料双非线性。螺栓连接中的接触分析，也是非线性的。
6.1 几何非线性
当变形大到"不能用变形前的几何来表征变形后的状态"，就需要几何非线性（大变形分析）。需要区分：
- 格林-拉格朗日应变：参考构形（变形前），适合固体大变形
- 欧拉-阿尔曼西应变：当前构形（变形后），流体力学中常用
几何非线性通过更新参考构形来处理（Total Lagrangian 或 Updated Lagrangian 格式）。
6.2 材料非线性
超过弹性极限后，材料进入塑性区域。弹塑性本构关系需要：
- 屈服准则（冯米塞斯等）
- 流动法则（相关或非相关流动）
- 强化规则（各向同性强化、随动强化或混合强化）
塑性应变 $\varepsilon^p$ 不可恢复，需在每个增量步结束后更新内部状态变量（应力返回映射算法）。
6.3 接触非线性
零件之间的接触是强非线性——接触面积会随载荷变化。常用算法：
- 罚函数法：在接触面引入弹簧刚度，实现简单，但穿透量与罚刚度有关
- 拉格朗日乘子法：精确满足不穿透条件，但自由度增加
- 增广拉格朗日法：两者的折中，大多数商业软件采用
6.4 牛顿-拉弗森迭代法
非线性问题 $\mathbf{R}(\mathbf{u}) = \mathbf{f} - \mathbf{K}_{int}(\mathbf{u}) = \mathbf{0}$ 用牛顿迭代求解：
$\mathbf{K}_T$ 是切线刚度矩阵，每次迭代更新。收敛准则通常是残差范数 $\|\mathbf{R}\| / \|\mathbf{f}\| < 10^{-3}$ 或位移增量范数 $\|\Delta\mathbf{u}\| < 10^{-6}$。
非线性分析经常遇到"不收敛"的问题，怎么处理？
不收敛原因很多，但最常见的几个：第一，载荷步太大——把载荷分成更多小步（自动时间步长）；第二，初始猜测离解太远——检查边界条件和初始状态；第三，模型有刚体运动——检查约束是否完整；第四，接触摩擦系数设得太大——从0开始逐步增加。遇到不收敛，先看残差向量的位置，通常就在问题所在的节点附近。
7. 动力分析基础
我的产品要做振动测试，CAE要怎么对应？模态分析是什么？
结构动力学最基本的方程是 $\mathbf{M}\ddot{\mathbf{u}} + \mathbf{C}\dot{\mathbf{u}} + \mathbf{K}\mathbf{u} = \mathbf{f}(t)$，三项分别是惯性力、阻尼力、弹性力。模态分析就是求这个方程的固有频率（自然频率）和振型——相当于找结构"最容易振动的方式"。
7.1 运动方程
$\mathbf{M}$ 是质量矩阵（由密度积分得到），$\mathbf{C}$ 是阻尼矩阵，$\mathbf{K}$ 是刚度矩阵。
7.2 特征值问题（模态分析）
自由振动（$\mathbf{f}=0$，$\mathbf{C}=0$）假设解为 $\mathbf{u} = \boldsymbol{\phi} e^{i\omega t}$，代入得广义特征值问题：
求解得 $n$ 对特征值-特征向量 $(\omega_i^2, \boldsymbol{\phi}_i)$，即固有频率 $f_i = \omega_i / (2\pi)$ 和振型 $\boldsymbol{\phi}_i$。振型关于质量矩阵和刚度矩阵正交：
7.3 数值时间积分——Newmark-β法
瞬态响应分析用数值积分求解运动方程。Newmark-β法（隐式）是最常用的一种：
常用参数：$\gamma = 0.5$，$\beta = 0.25$（平均加速度法），无条件稳定。$\gamma = 0.5$，$\beta = 0$ 为中心差分法（显式），有条件稳定（时间步需满足 $\Delta t < \Delta t_{cr} = 2/\omega_{max}$）。
7.4 频率响应分析（谐波分析）
对于频率为 $\Omega$ 的谐波激励 $\mathbf{f}(t) = \hat{\mathbf{f}} e^{i\Omega t}$，稳态响应为 $\mathbf{u}(t) = \hat{\mathbf{u}} e^{i\Omega t}$，代入运动方程：
扫频分析通过遍历频率范围求解，得到响应的频响函数（FRF）。当激励频率接近固有频率时，响应急剧放大——这就是共振。
那阻尼怎么设定？软件里的"瑞利阻尼"是什么？
瑞利阻尼是最常用的阻尼模型：$\mathbf{C} = \alpha\mathbf{M} + \beta\mathbf{K}$，用两个系数 $\alpha$ 和 $\beta$ 来控制阻尼比随频率的变化。实际工程中，金属结构阻尼比 $\zeta$ 通常在1%～5%之间。先从振动试验获得某两阶模态的阻尼比，反推 $\alpha$ 和 $\beta$。注意：瑞利阻尼在频率范围两端会过阻尼或欠阻尼，要小心。
8. 实际注意事项与常见错误
有没有新手最容易犯的错误汇总？我想提前避开这些坑。
太多了，我列几个最典型的。第一个：单位不一致——比如用mm做几何，Pa做应力，密度用kg/m³，结果质量矩阵差了9个数量级，动力结果完全错误。第二个：边界条件过约束——两个零件在接触面都施加了相同固定约束，相当于把刚度人为增加了。第三个：网格不够细——应力集中处（孔边、焊缝）网格粗的话，应力严重低估。
常见错误清单
- 单位系统混用：选择一套单位并全程一致（推荐：N, mm, MPa, t）
- 忽略网格收敛性检查：至少做3套不同密度的网格，确认结果收敛
- 弹性分析超出适用范围：应变超过1%就要考虑非线性
- 模态分析漏掉重要频率范围：提取的模态数应覆盖载荷频率的1.5～2倍以上
- 接触设置错误：忘记定义摩擦系数、接触刚度不合理
- 不检查反力：求解后第一件事是核对支反力是否与外力平衡
- 材料属性方向定义错误：复合材料、各向异性材料的坐标系要仔细确认
8.1 网格质量指标
| 质量指标 | 理想值 | 可接受范围 | 差网格警示 |
|---|---|---|---|
| 纵横比（Aspect Ratio） | 1.0（正四面体） | < 5 | > 20 → 结果不可信 |
| 雅可比比值（Jacobian Ratio） | 1.0 | > 0.6 | < 0.1 → 单元反转 |
| 偏斜度（Skewness） | 0（正形） | < 0.7 | > 0.9 → 严重问题 |
| 翘曲角（Warpage） | 0° | < 15° | > 45° → 重新网格 |
8.2 结构分析精度目标
实际工程中合理的精度预期：
- 线性静力分析：与解析解误差 <5%（网格已收敛时）
- 非线性分析：与实验数据误差 <10%（复杂接触可到15%）
- 模态分析：固有频率误差 <3%，振型相关系数（MAC）> 0.9
- 疲劳寿命预测：对数坐标下±1个数量级（疲劳本身离散性大）