气动弹性颤振分析
理论与物理
颤振是由流体力和结构弹性力动态耦合产生的自激振动。当超过某个临界速度（颤振速度 $U_F$）时，气动阻尼变为负值，振动会发散性增长。这在飞机机翼、涡轮叶片、桥梁缆索等结构中会成为问题。
振动发散…意思是说会一直持续到结构破坏为止吗？
正是如此。1940年的塔科马海峡大桥倒塌就是一个著名的例子。在飞机设计中，必须在整个飞行包线内确保相对于颤振速度有足够的裕度。FAA规定（FAR 25.629）要求颤振裕度达到飞行速度的1.15倍以上。
控制方程
用数学公式描述颤振的话会是什么样的呢？
让我们从典型的2自由度翼型颤振模型开始。以弯曲位移 $h$ 和扭转角 $\alpha$ 为自由度，运动方程如下。
其中 $m$ 是单位展长质量，$S_\alpha$ 是静不平衡矩，$I_\alpha$ 是绕弹性轴的惯性矩，$K_h$ 和 $K_\alpha$ 分别是弯曲和扭转的弹簧常数。
右边的气动力项 $L$ 和 $M_{EA}$ 如何表达呢？
非定常气动力使用 Theodorsen 函数 $C(k)$ 来表达。对于折合频率 $k = \omega b / U$，
Theodorsen 函数 $C(k) = F(k) + iG(k)$ 定义为汉克尔函数的比值，表示环量气动力的振幅和相位滞后。当 $k \to 0$ 时退化为定常气动力，$k \to \infty$ 时退化为非环量力。
对于一般结构的情况如何扩展呢？
使用模态坐标 $\mathbf{q}$ 的广义形式如下。
对于动压 $q_\infty = \frac{1}{2}\rho U^2$，广义气动力矩阵 $\mathbf{Q}(k)$ 通过偶极子格网法 (DLM) 或 CFD 计算。
颤振解法
求解颤振速度的方法有哪些？
我们来整理一下代表性的解法。
| 解法 | 特征 | 适用场景 |
|---|---|---|
| V-g法（美国法） | 人为引入结构阻尼 $g$ 进行特征值分析 | 初步设计阶段 |
| p-k法 | 迭代更新折合频率，更符合物理 | 详细设计 |
| p法（状态空间法） | 有理函数近似将气动力时域化 | 控制系统耦合 |
| CFD-CSD耦合 | 直接计算非线性气动力 | 跨音速颤振 |
p-k法中求解以下特征值问题。
$p = \gamma \pm i\omega$ 的实部 $\gamma$ 穿过零点的速度即为颤振速度。
跨音速区域不能用线性理论了吗？
没错。在跨音速区域，激波与边界层的相互作用导致气动力非线性增强。这个区域需要采用 CFD（欧拉方程或 RANS）与结构 FEM 直接耦合的 CFD-CSD 方法。
商用工具中的实现
有哪些软件可以用于颤振分析？
| 工具 | 方法 | 特征 |
|---|---|---|
| MSC Nastran (SOL 145/146) | DLM + p-k/V-g法 | 航空航天行业标准 |
| Ansys Mechanical + Fluent | CFD-CSD耦合 | 跨音速对应 |
| ZAERO (ZONA Technology) | ZONA6/7 非定常面元法 | 高速颤振 |
| STAR-CCM+ + Abaqus | 协同仿真 | 非线性FSI |
Nastran的SOL 145是行业标准啊。我记住了！
“Vn图”与“颤振裕度”——让设计师夜不能寐的气动弹性壁垒
飞机设计中最紧张的试验之一就是“颤振飞行试验”。国际标准（FAR Part 25）强制要求拥有设计颤振速度（VF）15%以上的裕度，在实际飞行试验中，会分阶段将速度提升至设计最大速度（VD）以检测颤振的萌芽。问题在于颤振是“突然出现的不稳定现象”——即使发现萌芽，也可能已经无法返航，因此试飞员的风险极高。1990年代，欧洲下一代客机的认证试验中，曾发生过理论计算忽略了外挂油箱安装导致重心移动的影响，试验机在比设计VF低20%的速度下就显示出颤振征兆的案例。理论的“疏漏”在飞行试验现场显现——因此颤振理论必须始终保持谦逊和怀疑的态度。
各项的物理意义
- 结构-热耦合项：温度变化引起的热膨胀诱发结构变形，变形又影响温度场。$\sigma = D(\varepsilon - \alpha \Delta T)$。【日常例子】夏天铁轨伸长导致缝隙变窄——温度上升→热膨胀→产生应力的典型例子。电子电路板在焊接后翘曲也是由于不同材料热膨胀系数差异所致。发动机的气缸体因高温部和低温部的温差产生热应力，最坏情况下会导致裂纹。
- 流体-结构耦合（FSI）项：流体压力·剪切力使结构变形，结构变形又改变流体区域的双向相互作用。【日常例子】强风下悬索桥的缆索振动（涡激振动）——风力使结构摇晃，摇晃的结构改变气流，进而放大振动。心脏血流与血管壁的弹性变形、飞机机翼的颤振（气动弹性不稳定性）也是典型的FSI问题。有时仅需单向耦合，但变形较大时双向耦合是必须的。
- 电磁-热耦合项：焦耳发热 $Q = J^2/\sigma$ 引起温度上升，温度变化改变电阻的反馈循环。【日常例子】电暖炉的镍铬合金线通电后发热（焦耳热）变红——温度升高则电阻改变，电流分布也发生变化。IH电磁炉的涡流发热、输电线路温度升高导致的垂度增加也是这种耦合的例子。
- 数据传递项：通过插值解决不同物理场间网格不匹配的问题。【日常例子】天气预报中将“气温数据”和“风力数据”结合计算体感温度时，如果各自的观测地点不同就需要插值——在CAE的耦合分析中，结构网格和CFD网格通常也不一致，因此界面处的数据传递（插值）精度直接关系到结果的可信度。
假设条件与适用范围
- 弱耦合假设（单向耦合）：当一方物理场影响另一方而反向影响可忽略时有效
- 需要强耦合的情况：FSI中的大变形、电磁-热耦合中温度依赖性较强的情况
- 时间尺度的分离：各物理场特征时间差异较大时，可通过子循环提高效率
- 界面条件的协调性：确认耦合界面处的能量·动量守恒在数值上得到满足
- 不适用的情形：三个以上物理场同时强耦合时，有时需要整体式方法
量纲分析与单位制
| 变量 | SI单位 | 注意事项·换算备忘 |
|---|---|---|
| 热膨胀系数 $\alpha$ | 1/K | 钢: 约12×10⁻⁶、铝: 约23×10⁻⁶ |
| 耦合界面力 | N/m²（压力）或 N（集中力） | 确认流体侧与结构侧的力平衡 |
| 数据传递误差 | 无量纲（%） | 插值精度依赖于网格密度比。5%以下为目标 |
数值解法与实现
将翼面分割成梯形面元，在每个面元上布置加速度势的双重子（doublet）。求解满足面元3/4弦长点处下洗边界条件的联立方程，从而求得压力分布。
$A_{ij}$ 是包含核函数的气动影响系数，依赖于马赫数 $M$ 和折合频率 $k$。
DLM的精度极限在哪里？
它基于线性势流理论，因此在以下情况下精度会下降。
- 跨音速区域（$M \approx 0.8$〜1.2）：激波的影响
- 大迎角：伴随分离的流动
- 机翼厚度效应显著的情况
这些情况需要使用基于CFD的方法。
模态提取与气动样条
结构模态和气动面元的网格完全不同吧。如何将它们连接起来呢？
使用样条插值。代表性的方法如下。
| 样条 | 方法 | 用途 |
|---|---|---|
| Infinite Plate Spline (IPS) | 薄板弯曲的基本解 | 翼面内位移插值 |
| Thin Plate Spline (TPS) | IPS改良版 | 一般翼结构 |
| RBF (Radial Basis Function) | 径向基函数 | 3D形状 |
| Beam Spline | 1D插值 | 高展弦比机翼 |
在Nastran中使用 SPLINE1/SPLINE2 卡片定义。将结构节点的位移插值到气动面元的控制点，反之将气动载荷映射到结构节点。
使用转换矩阵 $\mathbf{G}$，
这种力的转换保证了虚功守恒。
p-k法的实现步骤
请告诉我p-k法具体的计算步骤。
1. 从结构模型进行模态分析（SOL 103）→ 得到固有模态 $\phi_i$、固有频率 $\omega_i$
2. 用DLM计算各马赫数·折合频率对应的广义气动力矩阵 $\mathbf{Q}(M, k)$
3. 对每个速度 $U$ 进行迭代求解：
- 用初始估计的 $k$ 对 $\mathbf{Q}$ 进行插值
- 求解特征值问题得到 $p = \gamma + i\omega$
- 用新的 $k = \omega b / U$ 更新 $\mathbf{Q}$
- 重复直到 $k$ 收敛
4. 根据 $\gamma(U)$ 的曲线图，$\gamma = 0$ 穿过的速度即为颤振速度
迭代大概多少次收敛？
通常3~5次收敛。但是，在模态接近的情况下追踪会变得困难，因此模态追踪算法（基于MAC值追踪等）很重要。
CFD相关主题
この記事の評価ご回答ありがとうございます！参考に
なったもっと
詳しく誤りを
報告
相关主题
なった
詳しく
報告