功率谱密度响应分析
理论与物理
PSD响应分析是什么
老师，什么是“PSD响应分析”？
随机振动没有确定的波形吗？
是的。地震、路面不平度、喷气发动机的声学噪声、火箭发射振动……这些不是确定性的时间历程，而是用统计特性（PSD）来描述的。
输入PSD与响应PSD
随机振动的基本公式：
- $S_{in}(f)$ — 输入PSD（加速度 $g^2$/Hz、力 $N^2$/Hz 等）
- $H(f)$ — 频率响应函数（FRF）
- $S_{out}(f)$ — 响应PSD（位移 $mm^2$/Hz、应力 $MPa^2$/Hz 等）
FRF的平方×输入PSD = 输出PSD。很简单呢。
PSD响应分析使用频率响应分析的结果。首先求出FRF $H(f)$（模态法或直接法），然后乘以输入PSD $S_{in}(f)$ 即可。几乎没有额外的计算成本。
RMS值
响应的RMS（均方根）值是响应PSD的积分：
PSD面积的平方根就是RMS。RMS是“振动大小”的指标呢。
在随机振动设计中，3σ（3×RMS）是最大响应的参考标准。假设为正态分布，则振幅有99.7%的概率在3σ以内。通常将3σ用作设计安全系数。
Nastran
```
SOL 111 $ 模态法频率响应
CEND
METHOD = 10
RANDOM = 20
BEGIN BULK
RANDOM, 20, , FREQ, TABRND1, 100
TABRND1, 100, , ,
, 20., 0.01, 200., 0.01, 500., 0.04, 2000., 0.04, ENDT
```
Abaqus
```
*STEP
*RANDOM RESPONSE
20., 2000., 100, 1.
*PSD DEFINITION, NAME=base_psd, TYPE=BASE
20., 0.01
200., 0.01
500., 0.04
2000., 0.04
*END STEP
```
输入PSD是用表格（频率 vs. PSD值）来定义的呢。
在航空航天领域，MIL-STD-810或NASA-STD-7001的随机振动谱是用PSD规定的。将这个PSD表格输入到FEM中。
总结
要点：
- $S_{out} = |H|^2 \cdot S_{in}$ — FRF的平方×输入PSD = 输出PSD
- RMS = $\sqrt{\int S_{out} df}$ — 响应的有效值
- 3σ是最大响应的参考标准 — 正态分布的99.7%
- 频率响应分析的PSD版本 — 额外成本几乎为零
- SOL 111 + RANDOM（Nastran）, *RANDOM RESPONSE（Abaqus）
PSD（功率谱密度）的定义
功率谱密度S(f)定义为每单位频率带宽的平均平方值，单位为G²/Hz或Pa²/Hz。根据维纳-辛钦定理（由Norbert Wiener和亚历山大·辛钦在20世纪30年代独立证明），PSD等于自相关函数的傅里叶变换。从实测数据估计时，Welch重叠平均法（1967年提出）是标准方法，通过窗函数宽度来控制频率分辨率与方差的权衡。
各项的物理意义
- 惯性项（质量项）：$\rho \ddot{u}$，即“质量×加速度”。您有过急刹车时身体被向前甩出的经历吗？那种“被带走的感觉”正是惯性力。物体越重越难启动，一旦启动也越难停止。建筑物在地震中摇晃，也是因为地面突然移动而建筑物的质量“被落下”。静力分析中此项设为零，这是“因为缓慢施力所以加速度可以忽略”的假设。在冲击载荷或振动问题中绝对不能省略。
- 刚度项（弹性恢复力）：$Ku$ 或 $\nabla \cdot \sigma$。拉弹簧时会感觉到“想要恢复的力”吧？那就是胡克定律 $F=kx$，也是刚度项的本质。那么提问——铁棒和橡皮筋，用相同的力拉，哪个伸得更长？当然是橡皮筋。这种“不易伸长性”就是杨氏模量 $E$，它决定了刚度。常见的误解是：“刚度高=强度高”。刚度是“不易变形的程度”，强度是“不易破坏的程度”，是不同的概念。
- 外力项（载荷项）：体积力 $f_b$（重力等）和表面力 $f_s$（压力、接触力等）。可以这样想——桥上卡车的重量是“作用在整个内容物上的力”（体积力），轮胎压路面的力是“只作用在表面的力”（表面力）。风压、水压、螺栓紧固力……这些都是外力。这里容易犯的错误是：弄错载荷的方向。本想“拉伸”却变成了“压缩”——听起来像笑话，但在3D空间中坐标系发生旋转时确实会发生。
- 阻尼项：瑞利阻尼 $C\dot{u} = (\alpha M + \beta K)\dot{u}$。试着弹一下吉他的弦。声音会一直持续吗？不，会逐渐变小。这是因为振动能量通过空气阻力或弦的内部摩擦变成了热。汽车的减震器也是同样的原理——故意吸收振动能量来提高乘坐舒适性。如果阻尼为零会怎样？建筑物在地震后会一直摇晃下去。实际上不会这样，所以设置适当的阻尼很重要。
假设条件与适用范围
量纲分析与单位制
| 变量 | SI单位 | 注意事项·换算备忘 |
|---|---|---|
| 位移 $u$ | m（米） | 输入为mm时，载荷·弹性模量也需统一为MPa/N系 |
| 应力 $\sigma$ | Pa（帕斯卡）= N/m² | MPa = 10⁶ Pa。与屈服应力比较时注意单位制不一致 |
| 应变 $\varepsilon$ | 无量纲（m/m） | 注意工程应变与对数应变的区别（大变形时） |
| 弹性模量 $E$ | Pa | 钢：约210 GPa，铝：约70 GPa。注意温度依赖性 |
| 密度 $\rho$ | kg/m³ | mm制时为tonne/mm³（钢约为 10⁻⁹ tonne/mm³） |
| 力 $F$ | N（牛顿） | mm制时为N，m制时也统一为N |
数值解法与实现
PSD分析的计算
PSD响应是怎么计算的？
1. 固有频率分析（SOL 103 / *FREQUENCY）
2. 频率响应分析（SOL 111 / *SSD）→ 获取FRF $H(f)$
3. PSD响应计算 — $S_{out}(f) = |H(f)|^2 \cdot S_{in}(f)$
4. RMS值计算 — $x_{rms} = \sqrt{\int S_{out} df}$
步骤3和4是FRF的后处理，不需要额外的分析。
多自由度系统的PSD响应
在多自由度系统中，还需考虑模态间的相关性（交叉模态项）：
相关项反映了“模态的重叠”。在密集模态中，相关项很重要。
输入PSD的典型示例
| 标准 | 频率范围 | PSD值 |
|---|---|---|
| MIL-STD-810 Method 514 | 20〜2000 Hz | 0.04 g²/Hz（平坦部分） |
| NASA-STD-7001 | 20〜2000 Hz | 取决于等级 |
| IEC 60068-2-64 | 10〜500 Hz | 取决于产品类别 |
总结
基于Welch法的PSD估计步骤
Welch法的PSD估计步骤为：①将时间序列数据以50%重叠率分割（典型为512〜4096点），②对每个分段加Hann窗并进行FFT，③用采样间隔和分段数对FFT绝对值的平方进行归一化并平均。在Python中，用scipy.signal.welch(x, fs, nperseg=1024)一行代码即可计算。振动试验标准ASTM E1049建议设置记录时间，使统计自由度至少达到120以上（约60个分段以上）。
线性单元（一阶单元）
节点间线性插值。计算成本低，但应力精度低。注意剪切锁定（可通过减缩积分或B-bar法缓解）。
二阶单元（带中间节点）
可以表现曲线变形。应力精度大幅提高，但自由度约增加2〜3倍。推荐：应力评估很重要的情况。
完全积分 vs 减缩积分
完全积分：有过约束（锁定）风险。减缩积分：沙漏模式（零能量模式）风险。根据情况选择。
自适应网格
基于误差指标（ZZ估计量等）的自动细化。有效提高应力集中部位的精度。有h法（单元分割）和p法（增加阶次）。
牛顿-拉弗森法
非线性分析的标准方法。每次迭代更新切线刚度矩阵。在收敛半径内具有二次收敛性，但计算成本高。
修正牛顿-拉弗森法
切线刚度矩阵使用初始值或每隔几次迭代更新。每次迭代成本低，但收敛速度为线性。
收敛判定标准
力残差范数: $||R|| / ||F_{ext}|| < \epsilon$（通常 $\epsilon = 10^{-3}$〜$10^{-6}$）。位移增量范数: $||\Delta u|| / ||u|| < \epsilon$。能量范数: $\Delta u \cdot R < \epsilon$
载荷增量法
不一次性施加全部载荷，而是分小步增加。弧长法（Riks法）可以越过载荷-位移关系的极值点进行追踪。
なった
詳しく
報告