Power Spectral Density Response Analysis
Theory and Physics
What is PSD Response Analysis?
Professor, what is "PSD Response Analysis"?
It's a method to calculate a structure's response to Random Vibration (stochastic vibration input) in the frequency domain. PSD (Power Spectral Density) represents the frequency distribution of vibration energy.
Does random vibration not have a deterministic waveform?
Correct. Earthquakes, road surface roughness, jet engine acoustics, rocket launch vibrations... these are described not by deterministic time histories, but by statistical characteristics (PSD).
Input PSD and Response PSD
Basic formula for random vibration:
- $S_{in}(f)$ — Input PSD (acceleration $g^2$/Hz, force $N^2$/Hz, etc.)
- $H(f)$ — Frequency Response Function (FRF)
- $S_{out}(f)$ — Response PSD (displacement $mm^2$/Hz, stress $MPa^2$/Hz, etc.)
FRF squared × input PSD = output PSD. That's simple.
PSD response analysis uses the results of frequency response analysis. First, obtain the FRF $H(f)$ (modal method or direct method), then just multiply by the input PSD $S_{in}(f)$. Additional computational cost is almost zero.
RMS Value
The RMS (Root Mean Square) value of the response is the integral of the response PSD:
The square root of the PSD area is the RMS. So RMS is an indicator of the "magnitude of vibration," right?
In random vibration design, 3σ (3×RMS) is used as an estimate for the maximum response. Assuming a normal distribution, the amplitude is within 3σ with 99.7% probability. 3σ is often used as a design safety factor.
Nastran
```
SOL 111 $ Modal frequency response
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
```
So the input PSD is defined by a table (frequency vs. PSD value).
In aerospace, random vibration spectra from MIL-STD-810 or NASA-STD-7001 are specified as PSD. This PSD table is input into the FEM.
Summary
Key points:
- $S_{out} = |H|^2 \cdot S_{in}$ — FRF squared × input PSD = output PSD
- RMS = $\sqrt{\int S_{out} df}$ — Effective value of response
- 3σ is an estimate for maximum response — 99.7% of normal distribution
- PSD version of frequency response analysis — Almost zero additional cost
- SOL 111 + RANDOM (Nastran), *RANDOM RESPONSE (Abaqus)
Definition of PSD (Power Spectral Density)
Power spectral density S(f) is defined as the mean square value per unit frequency bandwidth, with units like G²/Hz or Pa²/Hz. According to the Wiener-Khinchin theorem (proven independently by Norbert Wiener and Aleksandr Khinchin in the 1930s), PSD is equal to the Fourier transform of the autocorrelation function. For estimation from measured data, Welch's overlapped segment averaging method (proposed in 1967) is standard, controlling the trade-off between frequency resolution and variance via window width.
Physical Meaning of Each Term
- Inertia term (mass term): $\rho \ddot{u}$, i.e., "mass × acceleration". Haven't you experienced your body being thrown forward during sudden braking? That "feeling of being pulled" is precisely the inertial force. Heavier objects are harder to set in motion and harder to stop once moving. Buildings shake during earthquakes because the ground moves suddenly while the building's mass "gets left behind". In static analysis, this term is set to zero, which assumes "forces are applied slowly so acceleration can be ignored". It absolutely cannot be omitted for impact loads or vibration problems.
- Stiffness term (elastic restoring force): $Ku$ or $\nabla \cdot \sigma$. When you stretch a spring, you feel a "force trying to return it", right? That's Hooke's law $F=kx$, the essence of the stiffness term. So here's a question—an iron rod and a rubber band, which stretches more under the same pulling force? Obviously the rubber. This "resistance to stretching" is Young's modulus $E$, which determines stiffness. A common misconception: "high stiffness = strong". Stiffness is "resistance to deformation", strength is "resistance to failure"—they are different concepts.
- External force term (load term): Body force $f_b$ (gravity, etc.) and surface force $f_s$ (pressure, contact force, etc.). Think of it this way—the weight of a truck on a bridge is a "force acting on the entire contents" (body force), while the force of the tires pushing on the road surface is a "force acting only on the surface" (surface force). Wind pressure, water pressure, bolt tightening force... all are external forces. A typical mistake here: getting the load direction wrong. Intending "tension" but it becomes "compression"—sounds like a joke, but it actually happens when coordinate systems are rotated in 3D space.
- Damping term: Rayleigh damping $C\dot{u} = (\alpha M + \beta K)\dot{u}$. Try plucking a guitar string. Does the sound continue forever? No, it gradually fades away. That's because vibration energy is converted to heat by air resistance and internal friction in the string. Car shock absorbers work on the same principle—they intentionally absorb vibration energy to improve ride comfort. What if damping were zero? Buildings would continue shaking forever after an earthquake. Since that doesn't happen in reality, setting appropriate damping is crucial.
Assumptions and Applicability Limits
- Continuum assumption: Treats material as a continuous medium, ignoring microscopic heterogeneity.
- Small deformation assumption (for linear analysis): Deformation is sufficiently small compared to initial dimensions, and stress-strain relationship is linear.
- Isotropic material (unless otherwise specified): Material properties are independent of direction (anisotropic materials require separate tensor definition).
- Quasi-static assumption (for static analysis): Ignores inertial and damping forces, considering only equilibrium between external and internal forces.
- Non-applicable cases: For large deformation/large rotation problems, geometric nonlinearity is required. For plasticity, creep, and other nonlinear material behaviors, constitutive law extension is needed.
Dimensional Analysis and Unit Systems
| Variable | SI Unit | Notes / Conversion Memo |
|---|---|---|
| Displacement $u$ | m (meter) | When inputting in mm, unify loads and elastic modulus to MPa/N system |
| Stress $\sigma$ | Pa (Pascal) = N/m² | MPa = 10⁶ Pa. Be careful of unit system inconsistency when comparing with yield stress |
| Strain $\varepsilon$ | Dimensionless (m/m) | Note the distinction between engineering strain and logarithmic strain (for large deformation) |
| Elastic modulus $E$ | Pa | Steel: ~210 GPa, Aluminum: ~70 GPa. Note temperature dependence |
| Density $\rho$ | kg/m³ | In mm system: tonne/mm³ (= 10⁻⁹ tonne/mm³ for steel) |
| Force $F$ | N (Newton) | Unify as N in mm system, N in m system |
Numerical Methods and Implementation
PSD Analysis Calculation
How is the PSD response calculated?
1. Natural Frequency Analysis (SOL 103 / *FREQUENCY)
2. Frequency Response Analysis (SOL 111 / *SSD) → Obtain FRF $H(f)$
3. PSD Response Calculation — $S_{out}(f) = |H(f)|^2 \cdot S_{in}(f)$
4. RMS Value Calculation — $x_{rms} = \sqrt{\int S_{out} df}$
Steps 3 and 4 are post-processing of the FRF; no additional analysis is required.
Multi-Degree-of-Freedom System PSD Response
For multi-degree-of-freedom systems, cross-modal correlation terms are also considered:
The correlation terms reflect "mode overlap". They are important for closely spaced modes.
Typical Examples of Input PSD
| Standard | Frequency Range | PSD Value |
|---|---|---|
| MIL-STD-810 Method 514 | 20 to 2000 Hz | 0.04 g²/Hz (flat portion) |
| NASA-STD-7001 | 20 to 2000 Hz | Level dependent |
| IEC 60068-2-64 | 10 to 500 Hz | Product category dependent |
Summary
PSD Estimation Procedure via Welch's Method
PSD estimation via Welch's method follows these steps: ① Split time series data with 50% overlap (typically 512~4096 points per segment), ② Apply Hann window to each segment and perform FFT, ③ Normalize the squared absolute value of the FFT by sampling interval and number of segments, then average. In Python, this can be done in one line with scipy.signal.welch(x, fs, nperseg=1024). Vibration testing standard ASTM E1049 recommends setting recording time so that statistical degrees of freedom are at least 120 (approximately 60 segments or more).
Linear Elements (1st Order Elements)
Linear interpolation between nodes. Computational cost is low, but stress accuracy is low. Beware of shear locking (mitigated with reduced integration or B-bar method).
Quadratic Elements (with Mid-Side Nodes)
Can represent curved deformation. Stress accuracy improves significantly, but degrees of freedom increase by about 2~3 times. Recommended: when stress evaluation is important.
Full Integration vs Reduced Integration
Full Integration: Risk of over-constraint (locking). Reduced Integration: Risk of hourglass modes (zero-energy modes). Choose appropriately for the situation.
Adaptive Mesh
Automatic refinement based on error indicators (e.g., ZZ estimator). Efficiently improves accuracy in stress concentration areas. Includes h-method (element subdivision) and p-method (order increase).
Newton-Raphson Method
Standard method for nonlinear analysis. Updates tangent stiffness matrix every iteration. Achieves quadratic convergence within convergence radius, but computational cost is high.
Modified Newton-Raphson Method
Updates tangent stiffness matrix at initial value or every few iterations. Cost per iteration is low, but convergence speed is linear.
Convergence Criteria
Force residual norm: $||R|| / ||F_{ext}|| < \epsilon$ (typically $\epsilon = 10^{-3}$ to $10^{-6}$). Displacement increment norm: $||\Delta u|| / ||u|| < \epsilon$. Energy norm: $\Delta u \cdot R < \epsilon$
Load Increment Method
Instead of applying full load at once, apply in small increments. The arc-length method (Riks method) can trace beyond extremum points on the load-displacement relationship.
Related Topics
なった
詳しく
報告