模态法过渡响应分析（模态瞬态）
模态法过渡响应分析（模态瞬态）的理论基础
什么是模态法过渡响应
老师，"模态法的过渡响应分析"和频率响应的模态法有什么区别？名字很像，但做的事情是不是不一样？
问得很好。共同点是"通过固有模态展开来缩约问题"。区别在输入和输出的领域。频率响应分析是求"恒定频率持续加振时的稳态响应"，过渡响应分析是对随时间变化的任意外力——冲击、地震波、阶跃荷载——求解时间历程响应。
具体在什么场景使用呢？
典型例子有3个。(1) 火箭发射时的卫星振动——发动机推力变化在几秒内驱动结构。(2) 汽车门关闭冲击——关门瞬间的冲击力传遍车体。(3) 建筑地震响应——在弹性范围内，数十个模态可预测免震建筑的行为。这些都是"荷载随时间变化的线性问题"，模态法特别有效。
支配方程与模态展开
起点是多自由度系的运动方程：
$[M]$为质量矩阵，$[C]$为阻尼矩阵，$[K]$为刚性矩阵，$\{F(t)\}$为随时间变化的外力向量。若有$N$自由度，每个时间步都要解$N \times N$的联立方程——这就是直接法。
100万自由度、1万个时间步的话，每次都解100万×100万的矩阵？那太可怕了……。
这就是模态叠加法（Mode Superposition Method）的威力所在。我们用固有模态展开位移：
$\{\phi_i\}$是第$i$个固有模态，$q_i(t)$是模态坐标（一般化坐标），$m$是采用的模态数（$m \ll N$）。将这个展开式代入运动方程，从左边乘以$[\Phi]^T$——
模态的正交性会让联立方程解耦，对吧？
完全正确。对于质量正规化的模态，$[\Phi]^T[M][\Phi] = [I]$、$[\Phi]^T[K][\Phi] = \mathrm{diag}(\omega_i^2)$成立。若阻尼是比例阻尼（Rayleigh型），$[\Phi]^T[C][\Phi] = \mathrm{diag}(2\zeta_i\omega_i)$也对角化。结果是各模态成为独立的单自由度系：
$f_i(t) = \{\phi_i\}^T\{F(t)\}$叫模态荷载（一般化力）。100万自由度简化为30个单自由度系统的计算——计算时间剧烈缩短。
模态参加系数与有效质量
"需要多少个模态"怎么判断呢？
关键指标是模态参加系数和有效质量。对于某个方向$\{d\}$（如X方向=[1,0,0,…]），第$i$模态的参加系数是：
有效质量为$M_{\mathrm{eff},i} = \Gamma_i^2$。全模态有效质量总和等于全质量：$\sum_{i=1}^{N} M_{\mathrm{eff},i} = M_{\mathrm{total}}$。实务中有效质量累计达到全质量的90%以上时采用该模态数——这是基本规则。原子能法规（NQA-1）和宇宙规格（ECSS-E-ST-32C）都明确规定90%的标准。
达到90%需要多少模态呢？
取决于结构。建筑一般10～30个模态，机器设备或卫星50～200个，管道系统可能需要300个以上。管道的低阶模态占有效质量的比例很分散。但即使增加模态数，时间积分的计算量也只是线性增长，远低于直接法。
与直接法的比较
| 项目 | 模态法（Modal Transient） | 直接法（Direct Transient） |
|---|---|---|
| 支配方程 | $m$个独立单自由度系 | 每个时间步解$N \times N$联立方程 |
| 计算速度（100万DOF） | 极快（与模态数成正比） | 慢（强烈依赖于$N$） |
| 非线性 | 不可用（接触、塑性、大变形） | 可以 |
| 频率范围 | 至采用模态的最高固有频率 | 受时间步长限制的全频带 |
| 输出灵活性 | 可单独评估各模态的贡献 | 直接输出物理量 |
| 阻尼处理 | 模态阻尼（比例阻尼前提） | 任意阻尼矩阵 |
模态法比直接法快多少呢？具体数字想知道。
100万自由度、缩约到100个模态的话，时间积分成本下降至1/10,000以下。加上前期的特征值解析开销，总体也能达到直接法的1/50～1/100。但要注意——冲击荷载这样的高频成分占主导的问题，所需模态数会膨胀，速度优势就消退。那时应考虑切换到直接法或显式算法。
NASA孕育的模态叠加法
模态叠加法在60年代NASA的阿波罗计划和土星五号火箭开发中得到大规模应用。那时代的计算机无法应对100万自由度的直接法，但通过缩约到数十～百个模态，使得振动解析成为实用化。Nastran的SOL 112（模态法过渡响应）正是那个时代诞生的求解器架构，至今60多年后仍是航空航天业界的标准手法。
模态法过渡响应分析（模态瞬态）的数值计算手法
Duhamel积分与解析解
把运动方程分离后，那个单自由度系怎么解呢？
解析上可用Duhamel积分（卷积积分）来表示：
这里$\omega_{di} = \omega_i\sqrt{1-\zeta_i^2}$是阻尼固有圆频率。物理上讲，"过去每一时刻$\tau$加的力$f_i(\tau)$，以指数衰减的形式影响到现在$t$"。
这个积分能解析求解吗？地震波那样的复杂波形……。
阶跃荷载或三角形波这样的简单波形可以解析求解，但实际的荷载波形几乎都需要数值积分。
数值时间积分的实现
单自由度系的时间积分运用Newmark-$\beta$法或中央差分法进行标量运算。Newmark法（$\beta=1/4, \gamma=1/2$，梯形法则）如下：
单自由度就是标量运算——不需要矩阵分解。
各模态能并行计算吧？GPU加速可以吗？
原则上完全并行化。各模态互相独立，所以可以。Nastran 2023以后支持模态法的GPU加速。实务上瓶颈通常在特征值解析和输出的物理坐标变换，而非时间积分。
残余向量（RESVEC）的作用
"残余向量"是什么？求解器手册里经常提到，但没明白。
非常重要的概念。模态法只用$m$个模态，荷载向量$\{F\}$中不能用采用模态表现的分量就丧失。残余向量（RESVEC）就是把这个"丧失的分量"作为伪模态补充。
具体怎么算呢？
从荷载向量减去采用模态的贡献，把余项用刚性矩阵的逆（静位移）来求解：
把这个$\{r\}$正规化后加为额外模态。效果明显——特别是荷载点附近的局部响应大幅改善。Nastran的PARAM,RESVEC,YES、Abaqus的RESIDUAL MODES选项对应这个。冲击解析时必须启用。
输出的物理坐标变换
模态坐标$q_i(t)$求出后，实际位移和应力怎么算呢？
物理坐标的变换：
应力是各模态应力分量乘以模态坐标后叠加：$\{\sigma(t)\} = \sum_{i=1}^{m} \{\sigma_i\}\, q_i(t)$。这里有实务技巧——全节点、全时间步做物理坐标变换会产生巨大数据，所以通常只在关注点（应力评估点、传感器位置）做变换。Nastran用SET卡限定输出节点。
输出数据爆炸的战争
100万自由度 × 1万时间步 × 6个应力分量 = 600亿个浮点数。双精度的话480GB。模态法让计算速度快了，但输出也会成瓶颈。某宇宙机构工程师说过"解析3分钟，后处理30分钟，文件复制2小时"。聪明的做法是保存模态坐标$q_i(t)$，后处理时按需对特定位置做物理坐标变换。
模态法过渡响应分析（模态瞬态）的实务应用
模态数的决定方法
有效质量90%是最低线，但实务上具体怎么决定模态数？
有3个判断标准。(1) 有效质量90%以上——这是最低底线。(2) 荷载的频率带宽——对输入荷载做FFT，找出包含99%能量的最高频率$f_{\max}$，采用达到$1.5～2$倍$f_{\max}$的固有振动数为止的模态。(3) 与直接法的对比验证——用小规模模型先与直接法对比，事先预估必要模态数。
比如火箭卫星整流罩分离冲击解析的话？
整流罩分离需要2000Hz左右的带宽。卫星结构的低阶固有振动数通常5～50Hz，所以2000Hz所需的全模态——有时需要数千个模态。这个量级模态法就失去速度优势，应考虑用RESVEC或转向直接法。而正弦波振动试验（5～100Hz）就50～100模态足够了。
时间步长的设置
时间步长$\Delta t$的设置基准有吗？
基本规则是最高关注振动数$f_{\max}$的周期的1/20以下：
另外重要的是，荷载波形的变化也要被充分采样。如果阶跃荷载的上升时间是0.001秒，那么$\Delta t = 0.0001$秒的量级必需。模态法通常用Newmark法（无条件稳定），不受稳定性约束，但精度受影响。
阻尼的设置与注意点
模态法的阻尼怎么设置？用Rayleigh阻尼吧。
模态法的大优点之一就是能单独给每个模态设定阻尼比$\zeta_i$。Rayleigh阻尼的话$\zeta_i = \frac{\alpha}{2\omega_i} + \frac{\beta\omega_i}{2}$，低阶和高阶的阻尼比会扭曲。模态法则可直接指定"1阶2%、5阶5%"这样物理合理的值。
没有实测数据的时候，阻尼比怎么决定呢？
没有实测就用经验值。钢结构：0.5～2%、混凝土：3～5%、螺栓连接结构：2～4%、橡胶、树脂：5～15%。这些都是参考值，振幅依赖性存在——小振幅时低，大振幅时高。NASA手册（NASA-HDBK-7005）有材料别推荐值列表。
SRS（冲击响应谱）评估
SRS（冲击响应谱）和模态法有关系吗？
关系很密切。SRS就是"对某一冲击输入，固有振动数$f_n$的单自由度系所示最大响应"作为$f_n$的函数的曲线。模态法过渡响应的各模态最大响应$|q_i|_{\max}$直接就是SRS的构成元素。宇宙机器的环境试验条件（NASA-STD-7003、ECSS-E-ST-10-03C）都用SRS规定，所以模态法过渡响应→SRS评估是宇宙业界的标准工作流。
实务检查表
模态法过渡响应解析的实务检查项目总结：
- 特征值解析正常完成 — 确认没有负特征值，剛体模式无
- 有效质量3方向都≥90% — 在求解器输出表中确认
- RESVEC（残余向量）启用 — 冲击荷载必须启用
- 时间步长$\Delta t \leq T_{\min}/20$ — $T_{\min}$是关注的最小周期
- 荷载波形采样充分 — 立上/立下的急峻部分要能捕捉
- 模态阻尼$\zeta_i$设置完整 — 未设置的话默认$\zeta=0
- 初始条件正确 — 通常为零，但前级阶段传递时注意
- 与直接法对比验证 — 小规模模型确认一致性
H-IIA火箭的卫星振动解析
JAXA H-IIA火箭发射时，卫星结构受到20～2000Hz的振动环境。卫星侧缩约到1000～5000自由度的Craig-Bampton缩约模型，与火箭侧的发射架模型联结，采用上位200个模态的模态叠加过渡响应来评估发射荷载。这个解析得出的SRS成为卫星试验等级（认证/验收）的根据。
模态法过渡响应分析（模态瞬态）的软件比较
Nastran SOL 112
Nastran SOL 112是模态法过渡响应的业界标准。输入的骨架：
SOL 112 $ 模态法过渡响应
CEND
METHOD = 10 $ 特征值解析法参照
TSTEP = 100 $ 时间步定义
DLOAD = 200 $ 动荷载定义
PARAM,RESVEC,YES $ 残余向量启用
BEGIN BULK
EIGRL, 10, , , 50 $ 抽取50模态
TSTEP, 100, 1000, 0.001 $ 1000步，Δt=0.001s
TLOAD1, 200, 300, , 0, 400 $ 时间荷载
TABLED1, 400, ... $ 荷载波形表重点是PARAM,RESVEC,YES。没这个的话冲击解析结果与直接法会偏离。
Abaqus *MODAL DYNAMIC
Abaqus的做法是连续执行2个阶段：
** 阶段1：特征值解析
*STEP
*FREQUENCY
50, , $ 抽取50模态
*RESIDUAL MODES $ 残余模态启用
*END STEP
** 阶段2：模态法过渡响应
*STEP
*MODAL DYNAMIC
0.001, 1.0 $ Δt=0.001s, 总时间=1.0s
*MODAL DAMPING
1, 50, 0.02 $ 模态1～50，阻尼比2%
*CLOAD
... $ 时间荷载
*END STEPAbaqus的特点是*MODAL DAMPING可按模态范围设不同阻尼比，很灵活。
Ansys MSUP过渡响应
Ansys Mechanical怎么设置呢？
Ansys Mechanical的GUI里依次设定"Modal"解析→"Transient (Mode Superposition)"。APDL级别是ANTYPE,MODAL→ANTYPE,TRANS流程。APDL命令骨架：
/SOLU
ANTYPE,MODAL ! 特征值解析
MODOPT,LANB,50 ! Block Lanczos抽取50模态
MXPAND,50 ! 50模态展开
SOLVE
FINISH
/SOLU
ANTYPE,TRANS ! 过渡响应
TRNOPT,MSUP ! MSUP（模态叠加）法
DELTIM,0.001 ! Δt=0.001s
MDAMP,ALL,0.02 ! 全模态阻尼比2%
F,NODE,FY,LOAD ! 荷载定义
TIME,1.0
SOLVETRNOPT,MSUP是选择模态法而非直接法的关键字。
功能比较表
| 功能 | Nastran SOL 112 | Abaqus *MODAL DYNAMIC | Ansys MSUP |
|---|---|---|---|
| 残余向量 | PARAM,RESVEC,YES | *RESIDUAL MODES | RESVEC,ON |
| 模态别阻尼 | TABDMP1 | *MODAL DAMPING | MDAMP |
| SRS输出 | PARAM,SRS（内置） | Python脚本 | POST26宏 |
| GPU加速 | 2023以后（特征值） | 受限 | 2024R1以后 |
| Craig-Bampton结合 | SEAM+SOL 112 | *SUBSTRUCTURE | CMS |
| 优势领域 | 航宇、防卫 | 通用、科研 | 汽车、通用产业 |
结局什么时候用什么软件？
这样分。航宇SRS评估、卫星振动→Nastran一択。SRS计算内置，与业界标准的CLA（耦合荷载解析）流程整合。通用模态法过渡+直接法频繁切换→Abaqus。同一模型轻松在模态法与直接法间切换。GUI工作流、汽车→Ansys Mechanical。Workbench的工作流管理便捷。
从模态法向直接法的切换判断
某项目用模态法算出的结果与直接法相差两成。原因是模型中混进了非线性弹簧单元（GAP单元）。GAP单元接触/分离时刚性变化，破坏固有模态的正交性。最后只好改为直接法。教训是："模态法能用吗"的判断得放在最前面。模型内哪怕混1个非线性单元，模态法也不能用。
模态法过渡响应分析（模态瞬态）的故障排除
直接法结果不一致
相同模型、相同荷载，模态法和直接法的结果峰值相差30%。什么原因？
"模态法直接法不一致"是模态法最常见的问题。优先顺序原因：
- 模态数不足 — 确认有效质量累计是否3方向都≥90%。特别Z方向（重力方向）只有80%的情况常被遗漏
- RESVEC未启用 — Nastran默认RESVEC=YES，但旧数据继承可能是RESVEC=NO
- $\Delta t$过大 — 模态法的$\Delta t$比直接法粗就捕捉不到荷载峰值
- 阻尼设置差异 — 直接法的Rayleigh阻尼与模态法的模态阻尼，各模态实效阻尼比会不同
- 混入非线性单元 — GAP单元、CBUSH非线性、MPC从属单元在模型中没有
5番的非线性单元最恐怖啊……求解器不报错吗？
Nastran在SOL 112中用非线性单元会出WARNING，但不会停止计算。把它线性化继续算。所以有输出，但与直接法不一致。"计算出了结果=结果正确"绝不能这么想。
冲击荷载初期响应不准确
冲击加载后的0.001秒间，模态法和直接法完全不合。但过一会儿会一致。
典型的RESVEC不足症状。冲击荷载的傅里叶变换是宽带的——数kHz成分都有。采用模态的最高固有振动数500Hz的话，500Hz以上完全消失。RESVEC补偿这个高频的静态分量，但动态分量（冲击后的高频振动）补不了。对策3种：
- 启用RESVEC（还没的话）
- 大幅增加模态数（最高固有振动数至少荷载带宽的2倍以上）
- 转为直接法（冲击初期精度要求高的话）
响应不衰减
荷载除去后响应不衰减，永远振动。
阻尼没设置。模态法里不明示就默认$\zeta=0$（很多求解器）。各求解器的设置方法不同：
- Nastran：TABDMP1卡定义模态阻尼表，CASE里设SDamp=ID
- Abaqus：*MODAL DAMPING行指定模态范围和阻尼比
- Ansys：MDMP命令单独或全模态一括设置
先全模态一律2%试，如有实测数据再调。
出现负特征值
特征值解析出现负的特征值。这样还能做模态法过渡吗？
绝对不能直接用。负特征值表示刚性矩阵正定性崩坏——模型哪里不稳定。常见原因：
- 约束不足 — 还有剛体模态（6个剛体模態的特征值≒0是正常，其他负值异常）
- 单元变形/翘曲 — 雅可比值降低到负程度的严重扭曲单元存在
- 初应力导入错误 — 预应力符号反或过大
负特征值存在就模态坐标$q_i(t)$发散。特征值解析阶段消除原因。
故障解决的铁则：最小重现模型
模态法过渡响应出故障，首先造个最小重现模型。100万单元的大模型跑1小时后出问题很难调试。不如造一个能复现问题的最小构成（10个单元的悬臂梁+同样荷载条件），几秒就解完。如此构成出现相同症状，就用它来逐个改变参数找根源。大模型里没症状的话，就把小模型逐步加大元素，找到哪个要素引起问题。专业调试是"减法"。