Random Vibration Fatigue: Method

Random Vibration Fatigue: PSD Integration, Spectral Moments, Miles, Dirlik, Narrow Band, Steinberg and Wirsching–Light

The problem

Equipment on a vehicle, a launcher or a shaker table sees a vibration that is random: no single frequency, no repeating waveform, only statistics. The standard description of a stationary Gaussian random vibration is its one-sided power spectral density W(f), the mean square acceleration per hertz, in g²/Hz against frequency in Hz. Its integral over frequency is the mean square acceleration, so the overall level is g_rms = √(∫W df). Test specifications (MIL-STD-810H Method 514.8, NAVMAT P-9492, NASA-HDBK-7005 and the GEVS) give W(f) as a few breakpoints joined by straight lines on log-log paper, with the overall g rms printed beside the curve. A structure mounted on that base amplifies the input near its natural frequencies, the stress at a critical point follows the response, and the question is how long the part lasts under a stress that is itself a Gaussian random process. The tool answers it in the frequency domain, without generating a time history, by the methods below.

Units and conventions

PSD levels are in g²/Hz and frequencies in Hz throughout; they are plain numbers and the unit switch leaves them alone. The stress per g (and the mass and area of the spring-mass form) go through the unit system (ksi or MPa, lbm or kg, in² or mm²) and are converted to SI base units for the calculation; stresses are reported in ksi (USCS) or MPa (SI). Durations are in hours, lives in hours (or minutes when short). The spectral moments are taken with frequency in Hz, m_n = ∫ fⁿ G(f) df, which gives rates directly per second; the dimensionless ratios (irregularity factor, Dirlik's x_m) are the same whether Hz or rad/s are used. The S-N curve is for fully reversed cycles: a zero-mean random stress has R = −1, and a static mean stress (gravity, preload) is not included.

Input PSD and its rms

Between two breakpoints (f₁, W₁) and (f₂, W₂) the specification is a straight line on log-log axes, W = W₁(f/f₁)^m with m = log(W₂/W₁)/log(f₂/f₁). The slope in dB per octave is 10 log₁₀(2^m) = 3.01 m: +3 dB/octave is m = 1 (W proportional to f), −6 dB/octave is m = −2. The area of a segment is exact, W₁f₁[(f₂/f₁)^m+1 − 1]/(m + 1), or W₁f₁ ln(f₂/f₁) when m = −1, and the overall level is the square root of the sum. The standard profiles on the dropdown reproduce the rms values printed in their sources: NAVMAT P-9492 Figure 5 (0.04 g²/Hz from 80 to 350 Hz with 3 dB/octave slopes to 20 and 2000 Hz; the document prints 6.0 g rms, the exact integral is 6.06 g), NASA-HDBK-7005 Table 8.9 (the NASA minimum workmanship level, 0.04 g²/Hz 80–500 Hz, 6.8 g rms) and Table 8.10 (USAF units, 0.04 g²/Hz 150–600 Hz with a −6 dB/octave roll-off, 6.1 g rms), MIL-STD-810H Figure 514.8E-1 (minimum integrity, 0.04 g²/Hz 20–1000 Hz then −6 dB/octave, 7.7 g rms, one hour per axis), Table 514.8C-I (common-carrier truck, vertical, 1.08 g rms for 1,000 miles in 60 minutes), Table 514.8C-II (the envelope for unknown orientation, 1.17 g rms) and Table 514.8C-IX (jet aircraft cargo general exposure, 4.02 g rms, one minute per take-off). The flat PSD of level W between f₁ and f₂ has rms √(W(f₂ − f₁)): 0.04 g²/Hz from 20 to 2000 Hz gives 8.9 g.

Response of a single-degree-of-freedom structure

A mass on a spring and damper fixed to the vibrating base has the absolute-acceleration transmissibility |H(f)|² = [1 + (2ζr)²]/[(1 − r²)² + (2ζr)²], r = f/f_n, Q = 1/(2ζ), and the response PSD is G_out(f) = |H|² W(f) (Bendat and Piersol). For a flat input the integral over all frequencies is exact, ∫|H|² df = (π/2) f_n Q (1 + 1/Q²), so g_rms = √[(π/2) f_n Q W (1 + 1/Q²)]; Miles (1954) dropped the 1/Q² term and his equation g_rms = √(π f_n Q W/2) is the everyday estimate of a resonant response, within 0.5 % for Q = 10. The tool integrates the actual PSD over its own band by Simpson's rule on a grid that clusters points across the resonance (about a thousand intervals), so the integrated rms differs from Miles' value when the PSD slopes near f_n or the band is cut: for the NAVMAT profile at f_n = 100 Hz, Q = 10 the integral gives 7.82 g rms, Miles 7.93 g. A measured or computed transmissibility can be pasted instead as |H| against frequency (interpolated linearly in |H| against log f, end values held), for example the response of a particular mode from a finite-element model.

Spectral moments, crossings and peaks

Rice (1944) showed that a stationary Gaussian process with one-sided PSD G(f) crosses zero upward at the expected rate E[0] = √(m₂/m₀) and has peaks at the rate E[P] = √(m₄/m₂), with m_n = ∫ fⁿ G df. Their ratio, the irregularity factor γ = E[0]/E[P] = m₂/√(m₀m₄), is 1 for a narrow-band process (every cycle crosses zero, so peaks and zero crossings are equally frequent) and falls toward 0 as high-frequency ripple adds peaks that do not cross zero; ε = √(1 − γ²) is the bandwidth parameter and α₁ = m₁/√(m₀m₂) the first bandwidth ratio; Dirlik's mean frequency is x_m = (m₁/m₀)√(m₂/m₄) = α₁γ. For a resonant response under a broadband input the fourth moment is dominated by the ripple above the resonance (G_out falls only as 1/r² there), so E[P] and γ depend on the upper frequency limit of the input: the NAVMAT response at 100 Hz, Q = 10, has E[0] = 101 Hz but E[P] = 143 Hz and γ = 0.71 when the input runs to 2000 Hz. The number of stress cycles in an exposure of T seconds is E[0]T, the number of peaks E[P]T.

Stress and the S-N curve

The stress at the critical point is proportional to the response acceleration, σ(t) = K_f σ_1g a(t)/g, where σ_1g is the stress for 1 g of response (from a static 1 g analysis or a unit load in a finite-element model; for a mass m on a member of stressed area A it is m g/A) and K_f the fatigue notch factor. The stress PSD is (K_fσ_1g)² G_out, so the stress rms is K_fσ_1g g_rms and its rates and irregularity are the response's. The S-N curve may be a Basquin law N = A/S_a^b for the amplitude S_a (with an optional endurance limit below which no damage is counted) or a MIL-HDBK-5J best-fit equation log N = A − B log(S_eq − C), S_eq = S_max(1 − R)^p, which for the fully reversed cycles of a zero-mean process (R = −1) becomes S_eq = 2^pS_a; a handbook curve with C = 0 is a Basquin law with b = B and A_bq = 10^A/2^pB, so the closed forms below apply to it exactly, while a curve with C > 0 is integrated numerically. Ten unnotched handbook curves (aluminum and steel) are on the dropdown, taken from the equations filed for the Spectrum Fatigue tool; their constants fill the fields and stay editable.

Damage by the narrow-band formula

For a narrow-band Gaussian process the peaks are Rayleigh distributed and each peak is followed by a trough of about the same size, so every peak is a cycle of range S = 2a: in the normalized range z = S/(2σ) the range density is p(z) = z e^−z²/2 and Miner's rule gives D_NB = E[P] T ∫ p(z)/N(σz) dz (Bendat 1964; Bishop 1988 eq. 2.57). With a Basquin curve the integral is Γ(1 + b/2)(√2σ)^b/A, the familiar closed form. Applied to a wide-band process the narrow-band formula pairs every peak with an equally large trough and so overestimates the large ranges: it is an upper bound on the rainflow damage, often by 30 to 50 % for a resonant response with a broadband tail (Bishop 1988, section 2.7).

Dirlik's rainflow-range density

Dirlik (1985) generated time histories from seventy spectra of different shapes, rainflow counted them and fitted a closed-form density to the counted ranges in terms of γ and x_m alone: p(z) = (D₁/Q) e^−z/Q + (D₂z/R²) e^−z²/2R² + D₃z e^−z²/2, an exponential term for the small ranges of the ripple, a Rayleigh term of reduced scale R and the full Rayleigh term, with D₁ = 2(x_m − γ²)/(1 + γ²), R = (γ − x_m − D₁²)/(1 − γ − D₁ + D₁²), D₂ = (1 − γ − D₁ + D₁²)/(1 − R), D₃ = 1 − D₁ − D₂ and Q = 1.25(γ − D₃ − D₂R)/D₁ (thesis pp. 128–129; the same relations are restated in Bishop 1988, section 3.2). The density integrates to one and its moments are closed, E[z^b] = D₁Q^bΓ(1 + b) + D₂(√2R)^bΓ(1 + b/2) + D₃2^b/2Γ(1 + b/2), which the test script uses to check the numerical integration. The damage is D = E[P] T ∫ p(z)/N(σz) dz with the same number of peaks. For a resonant response with a broadband tail (γ about 0.7 with α₁ near 1) the relations give a negative R; only R² enters the density, which remains a proper density, and the tool reports the value. Dirlik's method has been the reference frequency-domain method since Bishop and Sherratt's comparisons, usually within 10 to 20 % of a long rainflow count; the project's own seeded time-domain check (20,000 s of the default response, rainflow counted) is recorded in the References folder and in the test script.

Steinberg's three-band method

Steinberg (Vibration Analysis for Electronic Equipment) simplifies the Gaussian response to three stress levels: the response spends 68.27 % of its time within ±1σ, 27.18 % between 1σ and 2σ and 4.28 % between 2σ and 3σ (the normal probabilities; peaks beyond 3σ are ignored), and the E[0]T cycles of the exposure are counted at those amplitudes, D = n_1σ/N(σ) + n_2σ/N(2σ) + n_3σ/N(3σ). It is the method of the electronics industry for lead-wire and solder-joint fatigue and is easy to do by hand; its result lies between the narrow-band and Dirlik values for the usual S-N slopes because the 3σ truncation removes the rare large peaks that the Rayleigh tail carries.

Wirsching–Light correction

Wirsching and Light (1980) kept the narrow-band formula and corrected it with a factor fitted to rainflow counts of simulated wide-band histories: D = λ D_NB, λ = a + (1 − a)(1 − ε)^c, a = 0.926 − 0.033b, c = 1.587b − 2.323, ε = √(1 − γ²), for S-N slopes b from 3 to 10 and γ from 0.45 to 1 (Bishop 1988 eqs. 3.1–3.3). For the default case (γ = 0.71, b = 5.8) λ = 0.73. The factor depends on γ and b alone; Dirlik's density also uses the first moment, which Bishop shows reduces the scatter of the prediction.

Life, duration and life factor

Each method gives the damage of the entered exposure and the time to failure T/D; the two are proportional, so the time to failure does not depend on the duration entered. A life factor (scatter factor) multiplies the damage and divides the life: the handbook S-N curves are mean curves with a standard error of 0.2 to 0.5 decades on log life, and a factor of 4 is the usual aerospace allowance on a mean curve. Miner's rule takes the damage as linear and sequence-independent. D ≥ 1 means the exposure exceeds the mean life.

Fatigue damage spectrum

The fatigue damage spectrum (FDS) is the damage of the exposure for a family of single-degree-of-freedom structures of the same Q and stress per g, plotted against their natural frequency. It is the fatigue counterpart of the shock response spectrum and the basis of test tailoring (MIL-STD-810H Method 514.8 Annex A, Lalanne): two environments with the same FDS do the same damage to a structure of given f_n, which is how a long service exposure is compressed into a shorter, higher-level test. The response rms follows √(f_n W(f_n)) and the cycle count follows f_n, so across a flat PSD the damage rises with frequency as f_n^{1 + b/2}, and it peaks where the input PSD is high.

Reading the results

The input section tabulates every segment with its slope and area and prints the overall g rms beside the value the source prints. The response section gives Miles' estimate and the integrated rms, the moments table, E[0], E[P], γ and the number of cycles in the exposure, with the two PSDs drawn on log-log axes. The damage section lists the four methods with the damage, the time to failure and the ratio to the narrow-band result, in bold the headline method chosen on the form, and draws Dirlik's density against the Rayleigh density. The damage spectrum closes the page. A time to failure shorter than the specified test duration, or D ≥ 1, means the part will not pass the test on a mean S-N curve; the usual remedies are a lower response (a stiffer mount to move f_n where the PSD is low, or more damping), a lower stress per g (a thicker section, a shorter lever arm) or a notch-free detail.

Limits

Stationary Gaussian input and a linear structure (no clipping, no kurtosis correction, no sine-on-random); one mode at a time unless a full transmissibility is pasted; the S-N curve is used as entered with no mean-stress correction; the frequency-domain methods are approximations to a rainflow count, with the narrow-band formula as the upper bound and Dirlik the usual best estimate; short exposures where the number of cycles is small are outside the statistical description.