# Review and assessment of frequency-based fatigue damage models.

ABSTRACTSeveral popular frequency-based fatigue damage models (Wirsching and Light, Ortiz and Chen, Larsen and Lutes, Benascuitti and Tovo, Benascuitti and Tovo with [alpha].75, Dirlik, Zhao and Baker, and Lalanne) are reviewed and assessed. Seventy power spectrum densities with varied amplitude, shape, and irregularity factors from Dirlik's dissertation are used to study the accuracies of these methods. Recommendations on how to set up the inverse fast Fourier transform to synthesize load data and obtain accurate rainflow cycle counts are given. Since Dirlik's method is the most commonly used one in industry, a comprehensive investigation of parameter setups for Dirlik's method is presented. The mean error and standard deviation of the error between the frequency-based model and the rainflow cycle counting method was computed for fatigue slope exponent m ranging from 3 to 12. The results showed Ortiz and Chen, Benascuitti and Tovo with [alpha].75, Larsen and Lutes, Dirlik, and Benascuitti and Tovo to be significantly more accurate than Lalanne and Zhao and Baker. These five models have a tendency for the error variation to increase as the fatigue exponent m increases. For this study using Dirlik's seventy spectra, Ortiz and Chen, Benascuitti and Tovo with [alpha].75, and Larsen and Lutes Single Moment methods had lower mean error than Dirlik's method.

CITATION: Quigley, J., Lee, Y., and Wang, L., "Review and Assessment of Frequency-Based Fatigue Damage Models," SAE Int. J. Mater. Manf. 9(3):2016.

INTRODUCTION

Fatigue and durability of automotive components mainly focuses on two types of mechanical loadings: deterministic cyclic loading from the periodical engine firing events and random loading from the road excitations. Fatigue for engine related components/structures typically adopt an infinite life approach based on a sampled time series loading from engine idling to fuel cut off [1]. While for frame/body mounted components, due to the random nature of the loading, fatigue analysis is generally accomplished in two ways: the time series approach or the spectral method based on the power spectral density of the load inputs. For the time series approach, an adequate duration of the time-domain mechanical loading covering all random loading conditions is required. Rainflow cycle counting [2, 3] and the linear damage rule [4] are then used to estimate the damage and fatigue life. There are drawbacks for applying the time series approach on random loading. It is costly to record a sufficiently long time history of load data to capture all frequency contents especially the less frequent but highly damaging events. Based on the study by Bishop [5], it is also computationally expensive to evaluate stresses and the fatigue damage from time series data of long duration.

The most common way to represent the dynamic response of the random loading input is by spectral methods in the frequency domain using a power spectral density with the assumption that the loading is random, stationary and Gaussian. Following Rice's [6] description of the probability density function for load peaks, Bendat [7] showed how to calculate the fatigue damage for a narrowband random loading input. However, it was found that the narrow band approach is too conservative for estimating the damage because the random loadings on structures are more or less wide band in nature. Approaches were developed by Wirsching and Light [8], Ortiz and Chen [9], and Larsen and Lutes [10] to account for the difference between narrow-band and wide-band fatigue damages, by introducing an empirical correction factor which is a function of irregularity factor and fatigue slope exponent (a negative reciprocal of the fatigue strength exponent in Basquin's equation). A refined method which claimed having sound theoretical framework was introduced recently by Benascuitti and Tovo [11]. Both numerical and experimental verifications of the method have been reported elsewhere [12, 13, 14]. The other approach was originally introduced by Dirlik [15] to empirically describe the probability density function for approximating the rainflow amplitude distribution, which uses a combination of an exponential and two Rayleigh probability densities. Following a similar concept, Zhao-Baker [16] estimated the amplitude distribution using a combination of one Weibull and one Rayleigh density, and Lalanne [17] proposed a method using Rice's formula which Halfpenny [18] found gave excellent results.

The purpose of this study is three-fold: (1) to present the frequency-to-time domain process for random loading, (2) to establish the guidance for estimating the parameters in Dirlik's approach, and (3) to review all the mentioned frequency-based fatigue damage models and recommend the best one for use.

Applications of the spectral methods have been presented in numerous papers [5, 10, 11, 18, 19, 20, 21]. It was concluded that the Dirlik method is considered by far the most accurate and superior method in estimating rainflow fatigue damage. However, the papers do not address in detail the parameter setups for the spectral moment and damage calculation processes. A significant error could result if inaccurate parameters are used. Therefore, in this paper a thorough study of all the parameters in order to achieve accurate damage results from Dirlik's and other methods is presented. Seventy power spectral densities (PSD's) from Dirlik's dissertation were used for the calculations. Explanations are given on how to numerically integrate Dirlik's method and get results with high precision. A method to calculate the upper bound for integration is given and the effect on accuracy of integration step size is discussed.

The way to evaluate any frequency-based models is to compare the estimated damage value from the frequency domain method with that calculated from rainflow cycle counting histogram of the original time domain loading. Since no detailed study on the subject of converting a power spectral density to a time-domain loading has been reported, it is another purpose of the this paper to present the frequency-to-time domain conversion process including the details such as sampling rate, duration, resolution.

The third objective of this paper is to review the above mentioned frequency-based fatigue damage models and to present our findings based on the damage error assessment process and results. Fatigue damage versus the fatigue slope exponent is also studied for each method. The results indicate the accuracy of spectral methods from previously published sources may not hold if a different fatigue slope exponent is used.

PROCESSES FOR CALCULATING FATIGUE DAMAGE FROM RANDOM LOADS

Fatigue damage calculations for random loading are initially conducted in the time domain followed by rainflow cycle counting and the linear damage rule. By relating a random load's power spectral density in the frequency domain to the load's amplitude distribution in the time domain, fatigue damage calculations are then extended to the frequency domain. In order to verify the accuracy of any proposed frequency spectral methods, a synthesized time domain approach is developed. The general calculation process and methods for all three approaches are illustrated in Figure 1. In this paper, each approach will be considered briefly but the synthesized time series approach and primarily the frequency spectral approaches are discussed in detail.

Case 1: Time History Fatigue Damage Calculations

During data collection, obeying Nyquist's criteria of using sampling rate higher than two times the highest frequency of interest will reveal the maximum frequency of the loading. For accurate rainflow cycle counting one needs to identify the peaks and troughs and use a sampling rate that is at least ten to twenty times higher than the highest frequency of interest [22]. Time history tests are the most costly since they require real time data recording on actual hardware.

Case 2: Synthesized Time History Fatigue Damage Calculations

Taking a power spectral density of random loading specification, a time history can be synthesized on a computer by performing an inverse fast Fourier transform for a given time duration. As with time history calculations, a rainflow cycle counting of the peaks and troughs is extracted from the synthesized load data, followed by a damage calculation using Miner's linear damage rule. Again, for accurate rainflow cycle counting results, one needs to synthesize a time-domain loading history with a sampling rate that is at least ten to twenty times higher than the highest frequency of interest [22]. With the assumption that the data is stationary and Gaussian, this approach can simulate the results obtained from costly data collection. If the test duration is long enough, the damage rate will converge towards a mean value. The process of performing inverse fast Fourier transforms and rainflow cycle counting on long synthesized time histories can be computationally expensive. However, this approach has commonly been used in the development of frequency domain methods for comparisons with rainflow cycle counting histograms.

There are various ways to convert the load PSD input to time history load data. In this paper a process as shown in Figure 2 is used. A custom Fourier Filter is used which multiplies the digitized PSD curve and the Fourier transform of the input white noise. Then, an inverse Fourier transform is used to obtain the time series which is Gaussian and has the same PSD shape and amplitude as the original input. A confirmation is conducted by comparing the original PSD and the synthesized PSD for each load case.

Several parameters such as the sampling rate, duration, and custom filter buffer size are the keys for generating a representative time history load data from a PSD input. A pseudo damage analysis was conducted for various generated time series loads with different sampling rates. The analysis result is shown in Figure 3 where the X axis is the ratio between sampling rate and the highest frequency in the original PSD input and the Y axis is the pseudo damage of the time series load assuming a constant S-N curve. As shown, the sampling rate should be higher than ten times the peak frequency content to achieve a converged pseudo damage value. In this paper, the sampling rate was set as 5000 Hz for all the PSDs cases, which is more than 20 times higher than the highest loading frequency.

Duration of the generated time history load also affects the accuracy. Figure 4 shows the original versus synthesized PSD comparisons with varied durations for the same rectangular PSD input. As shown, the PSD of the generated time series load is much closer to the actual PSD input for a longer duration data than for a short data. To better represent the actual PSD shape, time series data with long enough duration is required. In this paper, the duration of the generated time series load is 10,000 seconds.

The accuracy of the Fourier transforms and inverse Fourier transforms is also affected by the buffer size of the analysis, which determines resolution when digitizing the PSD input and FFT of the white noise, and also the IFFT process. As shown in Figure 5, the higher the buffer size of the analysis, the better the correlation turns out to be. All analysis used in this paper has a buffer size of 8,192 for all PSD inputs.

Case 3: Frequency Domain Fatigue Damage Calculations

Frequency methods utilize parameters calculated from the spectral moments about the origin and integrated over the power spectral density. One advantage is that all frequency domain methods find a single mean damage rate regardless of the test duration. The usual assumption is that the power spectral density is stationary and Gaussian and is representative of the expected field loading condition. Some fundamental vibration and fatigue theories can be found in References [17, 23, 24, 25].

The general form of Basquin's damage equation for frequency domain calculations is

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (1)

E[freq] is the expected or average frequency, T is the test duration, E[Sa^m] is the expected load or stress raised to the fatigue exponent and Ca is the Basquin material constant. The subscript "a" in the above equation stands for "amplitude" which is the convention used throughout this paper, however range can also be used. Also, stress can be replaced by strain or displacement with a suitable change in the material constant, [C.sub.a].

The spectral moments about the origin of any order are found by

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (2)

The term [W.sub.x] is the PSDs with units (load^2/Hz) or (stress ^2/Hz) and [f.sup.j] is the frequency of order j.

The PSD is often a simplified envelope of the fast Fourier transform of time history data but with a reduced number of break points. Typically a PSD is given in the form of a discontinuous function shown as a table of values, or breakpoints. It is necessary to numerically integrate a discontinuous PSD to obtain the spectral moments. Lalanne [17] provides a number of equations for calculating spectral moments for various bases of the PSD axes. Most frequently, PSD's are defined on a log-log plot scale. The equation for integrating a segment on log-log scale from one break point, A, to the next one, break point B is

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (3)

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (4)

In the above equations, [W.sub.A] is the PSD magnitude at the first break point, k is the slope between breakpoints on the PSD plot, and j is the order of the spectral moment. Note that if k = -(j + 1) then

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (5)

The moments of any order are found by summing the individual segments between break points found using equation 3.

Two spectral parameters that are derived from the spectral moment of considerable note are the [zero.sup.th] order moment, [M.sub.0], which is the variance, and [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII], which is the standard deviation and also the area enveloped by the PSD. In addition, the expected rate of positive zero crossings E[[0.sup.+]] (the rate that load or stress changes sign from negative to positive) and the average rate of peaks E[P] (the average rate that load or stress reaches a peak value) are given by

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (6)

and

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (7)

Another important function is the irregularity factor which is given by

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (8)

A narrowband process has an irregularity factor of close to one and a for wideband process it is close to zero.

If the mean is not zero, then the mean stress effect must be accounted for also by the Morrow, Walker, modified Walker, or Smith-Watson-Topper equations [26]. It will be assumed for purpose of this paper that the mean value is zero, the structure is linear, and the loading is stationary and Gaussian.

In his often cited paper, S.O. Rice [6] described the probability density function of peak values for a narrowband process follows a Rayleigh distribution.

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (9)

[[sigma].sub.Sa] is the standard deviation of the random process.

Rice stated further that for any random stationary Gaussian process, narrowband or wideband, the probability density function is the weighted sum of a Rayleigh and normal distribution. This function is known as the Rice distribution and is given in equation 10. The function is plotted in Figure 6 versus normalized load for various irregularity factors. Note that for an irregularity factor of 1, the shape is a Rayleigh distribution and for an irregularity factor of zero, the shape is Gaussian.

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (10)

[PHI] is the standard normal distribution function.

A variety of methods for calculating fatigue damage from random loading have been created and many of the most popular and simple to calculate methods are listed in the following.

Narrowband Method (NB) [7]

By using the probability density function of narrowband peaks, Bendat [7] showed that the fatigue damage for a narrow band process can be found by

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (11)

Where [GAMMA]( ) is the gamma function.

The narrowband method tends to overestimate damage for wideband random loading. A number of methods apply a correction factor, [zeta], to the Narrowband method to obtain a more correct result so that

D = [zeta][D.sub.NB] (12)

In the following, several methods with varied correction factors are introduced.

Wirsching and Light Method (WL) [8]

D = [[zeta].sub.WL][D.sub.NB] (13)

where [[zeta].sub.WL] = a + [(1 - a)(1 - [lambda]).sup.b] (14)

a = .926 - .033m (15)

b = 1.587m - 2.323 (16)

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (17)

Ortiz and Chen Method (OC) [9]

D = [[zeta].sub.K][D.sub.NB] (18)

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (19)

The moments [M.sub.k] and [M.sub.k+2] are found using equation 3 and the order term j is replaced by k.

k =2.0/m (20)

Larsen and Lutes' Single Moment (SM) Method [10]

The Larsen and Lutes' single moment method was developed after extensive examination of simulation data and rainflow analysis, particularly for bimodal PSD's.

D = [[zeta].sub.SM][D.sub.NB] (21)

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (22)

The moment [M.sub.k] is found using equation 3 and the order term j is replaced by k.

k = 2.0/m (23)

The damage equation reduces to require the calculation of only a single moment, [M.sub.k], and hence its name.

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (24)

Benascuitti and Tovo's (BT) Method [11, 14]

Benascuitti and Tovo [11, 14] derived an expression for finding the expected rainflow fatigue damage for Gaussian loading. Rychlik [27] found that the expected damages for rainflow cycle counting damages will be less than or equal to the damage estimate found by the narrow-band approximation and greater than or equal to the damage estimate found by the range counting method. Recall the range counting method defines a cycle as the difference between a local peak and the preceding local minimum. Benascuitti and Tovo therefore sought to find a solution to the weighted linear combination between these two bounds as follows

[D.sub.RFC] - b[D.sub.NB] + (1 - b)[D.sub.RC] (25)

The range counting approximation was provided by Madsen [28],

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (26)

The damage is found by solving

[D.sub.BT] = [[zeta].sub.BT][D.sub.NB] (27)

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (28)

It is required to determine the following parameters

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (29)

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (30)

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (31)

The evaluation of the b coefficient was achieved by sampling 286 PSDs with [[alpha].sub.1] values ranging from 0.1 to 0.9 and [[alpha].sub.2] (irregularity factors) varying from 0.1 to 0.7. Rainflow fatigue damages were computed for these synthesized time histories. Note that the b coefficient requires only four spectral moments: [M.sub.0], [M.sub.1], [M.sub.2], and [M.sub.4], which in combination form two parameters, [[alpha].sub.1] and [[alpha].sub.2]. The approximation formula was found by a least squares best fit method that minimized the error between the simulated and approximated results.

Benascuitti and Tovo's [[alpha].sub..75] Method ([[alpha].sub..75]) [12]

Benascuitti and Tovo developed this method based on observations that the damage correction factor is a function of the bandwidth parameter, [[alpha].sub..75].

D = [[zeta].sub.[alpha]][D.sub.NB] (32)

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (33)

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (34)

For example, [M.sub.(.75)] is the moment when j=0.75 based on the moment calculation as in Equation 5.

The following methods solve for the probability density function for the load.

Dirlik's Method [15]

Dirlik's method solves for the probability density function for the stress amplitudes using the equation below.

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (35)

Dirlik's probability density function for rainflow stress amplitudes is

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (36)

Dirlik's stress probability density function parameters are

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (37)

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (38)

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (39)

[D.sub.3] = 1 - [D.sub.1] - [D.sub.2] (40)

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (41)

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (42)

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (43)

Dirlik's probability density function for stress cannot be integrated directly so discrete numerical integration must be used. To do that, upper bound for the load or stress value, [S.sub.a], must be estimated which will be discussed later.

Zhao and Baker's Method (ZB) [16]

The probability density function for the normalized stress amplitudes is expressed as a linear combination of a Rayleigh and Weibull PDF and functions of the irregularity factor [gamma] and a bandwidth parameter, [beta].

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (44)

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (45)

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (46)

[infinity] = 8 - 7[gamma] (47)

[beta] = 1.1 for [gamma] < 0.9 and [beta] = 1.1 + 9([gamma] - 0.9) for [gamma] [greater than or equal to] 0.9 (48)

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (49)

Lalanne's Method (LaL) [17]

Lalanne used Rice's formula as follows

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (50)

where the probability density function for stress amplitude is obtained from Rice's formula, Eq. 10.

Upper Bound Estimation

For methods requiring integration (Dirlik, Zhao and Baker, Lalanne), an upper band of integration needs to be estimated. Rice's formula and Figure 6 show that the probability density for the highest magnitude stress peaks is greatest for a narrowband process. It follows that the probability density for the highest magnitude stress peaks for a wideband process will be less. Put another way, the probability density for maximum peak stress values of a narrowband process is equal or greater than any wideband process. Therefore the narrowband process can safely be used for estimating the upper bound of high level stress peaks for any random process. The probability density of peak values of a narrowband process was given in equation 9.

The fatigue slope exponent, m, should also be considered since the fatigue slope exponent amplifies the damage for higher stress values. Therefore, a damage probability density function for a normalized stress term [S.sub.a.sup.m] can be written as follows:

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (51)

For probability density functions, the probability P of an event occurring between two bounds is a and b is

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (52)

Given a probability of stress values for a narrowband process, equation 51 can be used to solve for the upper bound of the stress amplitude. Filling in the appropriate variables gives

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (53)

Integrating from the lower bound 0 to the upper bound [UB.sub.Sa] gives

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (54)

From equation 1 and Bendat's formula equation 11,

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (55)

By noting that [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] and making the further substitution

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (56)

The probability reduces to

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (57)

The term in the denominator is the complete gamma function and term in the numerator is known as the upper incomplete gamma function. Factoring out Z and also [UB.sub.Sa] requires computing the inverse gamma function. These are tabulated below for various probabilities and fatigue slope exponent, m.

The simple steps to find the upper bound for stress amplitude are:

1. Determine or select a probability value P

2. Obtain a value for Z under the appropriate fatigue slope exponent m

3. Determine the upper bound for stress with the relation

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (58)

where [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] is the standard deviation of the load for the random process in terms of stress amplitude.

For a component with an expected life of 10 million stress cycles and a fatigue slope exponent of 3, the upper bound for normalized stress is [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII]. The range of m for most materials is 3<=m<=10.

For rainflow cycle counting, it is normal industry practice to sort the range and mean counts into at least 64 range bins and 64 mean bins. Estimating the upper bound and bin sizes for the stress amplitude can be done in the manner just given.

The upper bound mean values typically have a positive and negative extreme of similar magnitude. The maximum normalized mean values for a number of spectra with similar numbers of cycles and varying irregularity factor are plotted below in Figure 7. The data suggests that the maximum mean stress, positive or negative, with respect to irregularity factor appears to be related to the spectral width parameter (Equation 17).

To estimate the upper or lower bound for the means of the rainflow cycle counting (RFC) histogram, use the following formula.

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (59)

The same value for Z selected from Equation 56 and Table 1 can be used in Equation 59.

To test different numerical integration methods, needed for Dirlik's method, the Midpoint Rule, Trapezoidal Rule and Simpson's Rule were compared by integrating over a Rayleigh distribution for various probabilities from .001 to .999999999. The numbers of steps were varied and compared with the exact result using a spreadsheet (Excel). The Trapezoidal Rule converges approximately 2 times faster than the Midpoint Rule and Simpson's Rule about 8 times faster. However, when integrating up to high probability levels the Trapezoidal and Simpson's Rule methods become similar.

To integrate a Rayleigh distribution and get accuracy better than .01%, the minimum number of steps needed using the Midpoint of Trapezoidal Rule is 256 and when using Simpson's Rule was 32 steps. Simpson's Rule is more difficult to set up than the other methods but the potential accuracy is orders of magnitude better. Later it was found the number of steps needed was higher.

EXPERIMENTS

Dirlik's seventy random loading spectra were used for this experiment. Fifty six are rectangular bimodal spectra and 14 spectra are smooth, continuous functions. The spectral moments were calculated using equations 3, 4 and 5 between each break point and summed. For the continuous spectra, the PSD function was subdivided into very small .01 Hz segments, and the moment for each segment summed in the same fashion. The frequency range is 2 Hz to 242 Hz. Dirlik's spectra all have an expected number of peaks E[P]=108 peaks/block, or cycles/second, and an rms of 176.78.

Moments were tabulated for each PSD's. The process used was:

1. Synthesize a time history of each spectrum using an inverse Fourier transform for the specified time duration.

2. Perform rainflow cycle counting on the time series data.

3. Compute damage from the RFC results based on specified material properties.

4. Using the specified wideband method, compute the damage accumulated during the specified time duration.

5. Compute the damage ratio DR, damage found by wideband method divided by damage found by rainflow cycle counting. (If DR =1, there is no error)

6. Compare the methods by analyzing the error defined by

error = DR - 1. (60)

7. Compare the means and standard deviations of the error for all seventy spectra and various values of m.

For synthesizing the time history data the inverse fast Fourier transform had a sampling frequency of 5000 Hz. To get a sufficiently high number of samples, the duration for each test chosen was 10,000 seconds making the expected number of peaks was 1,080,000 peaks per test.

For the rainflow cycle counting data, the maximum stress amplitude was estimated to be 1200 based on a chosen values of Z=6.8 from Equation 57 and Table 1. The maximum and minimum histogram size for the mean stress was found by the method described earlier. The rainflow cycle count data was grouped into 1200 bins for stress amplitude and 60 bins for mean stress. The mean stress for each test was verified to be zero and all of the counts for mean stress were summed into a single value for each stress amplitude bin.

For all runs, the Basquin material constant used was Ca= 5E08 MPa and the fatigue strength coefficient used was Sf'=1000 MPa. These constants and the absolute damage results were not important to the relative comparisons sought. The fatigue slope exponent was also varied to investigate this variable.

The upper bound for integrating Dirlik's probability density function was found assuming Z=8, which was found to be sufficiently high based on table 1 and the expected number of peaks. A conversion was found to be needed for Dirlik's spectral values in order to obtain the required root mean square value of 176.78 as well as the same relative mean [X.sub.m] or irregularity factor, [gamma]. Dirlik's spectral parameters, [A.sub.1], [A.sub.2], [f.sub.1], [f.sub.2], [p.sub.1], and [p.sub.2] were used to create new PSD amplitudes [B.sub.1], and [B.sub.2], that matched Dirlik's rms value and other spectral moments, the following conversion was made:

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (61)

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (62)

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] (63)

This resulted in the approximate conversion of [MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII]

DISCUSSION

During preliminary tests integrating the Rayleigh probability density function, the trapezoidal rule method was found to converge faster than the midpoint rule and to be easier than Simpson's rule to set up the formulas. The trapezoidal rule also gave the same results as Simpson's rule for 64 steps and higher. However, when the trapezoidal rule was applied to this experiment, it took far more steps to converge Dirlik's density function than the Rayleigh density function. To converge to a point where a further doubling of steps produced less than 0.5% change, 2048 steps was needed for Dirlik's method and this value was also used for Lalanne's and Zhao and Baker's methods.

The effect of fatigue slope exponent was studied from m=3 to m=12. The error was calculated using Equation 60 for all seventy spectra. The standard deviation and mean error of the error for each method were then calculated. Results for the standard deviation and average error for varied fatigue slope exponent m are tabulated in Table 2 and Table 3. For each value of m in Table 2 and Table 3, the best value for standard deviation and mean for all the methods is shown in bold.

The standard deviations for all methods for a fatigue slope exponent from m=3 to m=12 are plotted in Figure 8. The methods with the lowest standard deviations are BT, Dirlik, [[alpha].sub..75], SM and OC with WL also fitting in this category except for m=3 to m=8. Standard deviation of the errors for all methods, with the exception of Lalanne, NB and WL, generally increases as fatigue slope exponent m increases.

The absolute values of the mean error for all methods for fatigue slope exponents from m=3 to m=12 are plotted in Figure 9. As the exponent increases, the mean error increases for all methods except for LaLanne and WL. Note WL reverses slope around m=8 in this plot of absolute values but in fact the mean error changes from positive to negative above m=8. The least mean error was for OC, [[alpha].sub..75], SM, Dirlik and BT having 6% or less at m=3.

Each method is plotted with mean error and standard deviation error in Figures 10 to 18.

LaL, NB and the ZB methods have a large positive mean error and a large variation in error. WL shows a negatively sloping error and the least variation for high fatigue exponents. In terms of mean error and standard deviation of the error, OC, [[alpha].sub..75], SM, Dirlik and BT are the best methods overall. Larsen and Lutes [10] did report that OC was also less accurate for certain bimodal spectra with modes of greater dissimilarity than are found in Dirlik's spectra. This experiment was done using Dirlik's spectra so it might be expected that Dirlik's method would do better than the others since it was originally derived from these. However, SM, [[alpha].sub..75], and OC methods generally were found to have more accurate mean errors than Dirlik's. Only BT had a better standard deviation of the error than Dirlik's.

Only Dirlik's, ZB and LaL methods solve for the probability density function for the load by numerical integration. Dirlik's method was found to have less than half the mean error of the NB and LaL methods. There are three potential sources of error when numerically integrating the probability density function for the stress amplitude.

These are 1) the choice of integration method, 2) the upper bound for load selected, and 3) the number of integration steps. Because these are unnecessary for calculating by OC, [[alpha].sub..75], SM, and BT methods, it can be concluded that these methods have an advantage over integration methods by having less risk of computational error.

These last four methods utilize the Narrowband method with a correction factor. SM, [[alpha].sub..75,,], OC, and BT methods require the calculation of 1, 3, 4, and 4 spectral moments respectively and this may be the only significant difference to the user. However, by utilizing equations 3, 4 and 5, the number of spectral moments to calculate is simple and not a deciding factor.

It is difficult to rank methods overall based on mean errors, standard deviations of the error, fatigue slope exponents and ease of application to the user. We offer this ranking of methods for the reader's consideration, in order of preference: [[alpha].sub..75], OC, Dirlik, BT and SM. We leave the final decision to the reader.

CONCLUSIONS

The Narrowband method with a correction factor continues to be a cornerstone for calculating fatigue damage for wideband random loading. The key to accuracy is calculating the spectral moments for the random loading power spectral density function and the equations provided here will work for continuous and discontinuous functions. In setting up the inverse fast Fourier transform to synthesize load data and obtain accurate synthetic rainflow cycle counts, a high sampling rate, high buffer size and long test durations are needed. Recommendations for these settings were given.

The upper bound for stress is related to the standard deviation of the power spectral density. An equation and table was provided for estimating the upper bound utilizing the Narrowband probability density function as a function of the probability of interest. Estimating the upper bound for stress is important for setting the bounds of integration for Dirlik's, Zhao and Baker and Lalanne's methods and is also required for binning data outputted from rainflow cycle counting. An equation for estimating the maximum mean values for binning rainflow counts was also given as a function of the spectral width parameter.

An attempt was made to estimate the best step size and numerical integration method to use for Dirlik's method by checking convergence using a Rayleigh probability density for which an exact solution can be found. However, in practice, Dirlik's method required 2048 steps to reach convergence of less than 0.5%.

Nine fatigue damage calculation methods for random loading were tested using the seventy spectra from Dirlik's thesis. The mean error and standard deviation of the error was computed for fatigue slope exponent ranging from 3 to 12. The results showed the Ortiz and Chen, Benascuitti and Tovo [[alpha].sub..75], Larsen and Lutes Single Moment, Dirlik and Benascuitti and Tovo methods to be significantly more accurate than the Lalanne, the Narrowband and Zhao and Baker methods. These five methods have a tendency for the error variation to increase as the fatigue slope exponent m increases. For this study using Dirlik's seventy spectra, Ortiz and Chen, Benascuitti and Tovo [[alpha].sub..75], and Larsen and Lutes Single Moment methods had a lower mean error than Dirlik's method.

REFERENCES

[1.] Wang, L., Lee, Y.L., Burger, R., and Li, K. "Multiple Sinusoidal Vibration Test Development for Engine Mounted Components," Journal of Failure Analysis and Prevention, Vol. 13, No. 2, 2013, doi: 10.1520/JTE20120035.

[2.] Matsuishi M. and Endo T., Fatigue of Metals Subjected to Varying Stress," Presented at Japan Society of Mechanical Engineers, Fukuoka, Japan, 1968.

[3.] Lee, Y.L., Pan, J., Hathaway, R., Barkey, M., Fatigue Testing and Analysis: Theory and Practice, Elsevier Butterworth-Heinemann, Oxford, UK, 2005.

[4.] Miner, M. A., "Cumulative Damage in Fatigue," Journal of Applied Mechanics, Vol. 67, A159-A164, 1945.

[5.] Bishop, N.W.M., The Use of Frequency Domain Parameters to Predict Structural Fatigue, PhD Thesis, University of Warwick, 1988.

[6.] Rice, S.O., "Mathematical Analysis of Random Noise, " Selected Papers on Noise and Stochastic Process, Dover, New York, 1954.

[7.] Bendat, J.S., Probability Functions for Random Responses, NASA Report on Contract NAS-5-4590, 1964.

[8.] Wirsching, P.H. and Light, M.C., "Fatigue under Wide Band Random Stresses," Journal of Structural Engineering, ASCE, 106(7), p. 1593-1607, 1980.

[9.] Ortiz, K. and Chen, N.K. "Fatigue Damage Prediction for Stationary Wide-Band Stresses," presented at the 5th International Conference on the Applications of Statistics and Probability in Civil Engineering, Vancouver, Canada, 1987.

[10.] Larsen, C.E. and Lutes L. D., "Predicting the Fatigue Life of Offshore Structures by the Single-Moment Spectral Method," Probabilistic Engineering Mechanics, Vol. 6, N0. 2, 96-108, 1991.

[11.] Benasciutti, D. and Tovo, R., "Spectral Methods for Lifetime Prediction under Wide-band Stationary Random Processes," International Journal of Fatigue, Vol. 27, 867-877, 2005.

[12.] Mrsik, M., Slavic, J., Boltezar, M., "Frequency-domain Methods for a Vibration-Fatigue-Life Estimation Application to Real Data," International Journal of Fatigue, Vol. 47, 8-17, 2013.

[13.] Bouyssy, V., Nboishikov, S.M. and Rackwitz R., "Comparison of Analytical Counting Methods for Gaussian Processes," Structural Safety, Vol. 12, 35-57, 1993.

[14.] Benascuitti, D., Fatigue Analysis of Random Loadings, PhD Thesis, University of Ferrara, 2004.

[15.] Dirlik, T., Application of Computers in Fatigue Analysis, PhD Thesis, University of Warwick. 1985.

[16.] Zhao, W. and Baker M.J. "On the Probability Density Function of Rainflow Stress Range for Stationary Gaussian Processes." International Journal of Fatigue, 14(2), 121-135, 1992.

[17.] Lalanne, C., Mechanical Vibration and Shock, Vol. 3, 4 and 5, Hermes Penton Science, London, 2009.

[18.] Halfpenny, A. and Kihm, F., "Rainflow Cycle Counting and Acoustic Fatigue Analysis Techniques for Random Loading," 10th International Conference RASD, 12-14 July 2010 Southampton.

[19.] Halfpenny, A., "Frequency Domain Approach for Fatigue Life Estimation from Finite Element Analysis," International Conference on Damage Assessment of Structures (DAMAS 99), Dublin, Ireland, June 28-30, 1999

[20.] Nieslony, A., Ruzicka, M., Papuga, J., Hodr, A., Balda, M., and Svoboda, J., "Fatigue Life Prediction for Broad Band Multiaxial Loading with Various PSD Curve Shapes," International Journal of Fatigue, Vol. 44, 74-88, 2012.

[21.] Gao, Z. and Moan, T., "Frequency-domain Fatigue Analysis of Wide-band Stationary Gaussian Processes Using a Trimodal Spectral Formulation," International Journal of Fatigue, Vol. 30, 1944-1955, 2008.

[22.] Donaldson, K., "Field Data Classification and Analysis Techniques," SAE Technical Paper 820685, 1982, doi:10.4271/820685.

[23.] Wirsching, P.H., Paez, T.L. and Ortiz, H., Random Vibrations Theory and Practice, Mineola: Dover, 1995.

[24.] Lee, Y.L., Barkey, M., Kang, H.K., Metal Fatigue Analysis in CAE: Practical Problem-Solving Techniques for Computer-Aided Engineering, Elsevier Butterworth-Heinemann, Oxford, UK, 2012.

[25.] Lutes, L.D. and Sarkani, S., Stochastic Analysis of Structural and Mechanical Vibrations, Englewood Cliffs: Prentice-Hall; 1997.

[26.] Dowling, N., "Mean Stress Effects in Stress-Life and Strain-Life Fatigue," SAE Technical Paper 2004-01-2227, 2004, doi:10.4271/2004-01-2227.

[27.] Rychlik, I., "Note on cycle counts in irregular loads," Fatigue and Fracture of Engineering Materials and Structures, 16:377-90, 1993.

[28.] Madsen, H.O., Krenk, S. and Lind, N.C., Methods of Structural Safety. Englewood Cliffs: Prentice-Hall; 1986.

CONTACT INFORMATION

John Quigley

General Dynamics Land Systems

38500 Mound Road, M/Z 435-01-40

Sterling Heights, MI 48310

quiglejp@gdls.commailto

quiglejp@gdls.com

(Approved for public release, LogNo.-2015-55, Distribution Unlimited, January 11, 2016)

Yung-Li Lee

FCA US LLC

800 Chrysler Dr.

Auburn Hills, MI 48326

yung.lee@fcagroup.com

Liang Wang

FCA US LLC

800 Chrysler Dr.

Auburn Hills, MI 48326

liang.wang@fcagroup.com

ACKNOWLEDGMENTS

The authors would like to thank HBM nCode, Southfield, MI for their generous support.

DEFINITIONS/ABBREVIATIONS

[A.sub.1] - Dirlik spectral parameter

[A.sub.2] - Dirlik spectral parameter

[A.sub.1-2] - area under PSD between points 1 and 2

a - WL parameter

[B.sub.1] - Dirlik spectral conversion parameter

[B.sub.2] - Dirlik spectral conversion parameter

b - WL parameter, BT parameter

BT - Benascuitti and Tovo method

[BT[alpha].sub..75] - Benascuitti and Tovo [[alpha].sub..75] method

[C.sub.a] - Basquin material constant in terms of amplitude

D, [D.sub.i] - damage

[D.sub.1] - Dirlik parameter

[D.sub.2] - Dirlik parameter

[D.sub.3] - Dirlik parameter

[D.sub.WB] - wideband damage

[D.sub.Dirlik] - damage by Dirlik method

[D.sub.RFC] - damage by rainflow cycle counting method

DR - damage ratio

E[freq] - expected, or average, frequency

E[[0.sup.+]] - expected number of positive zero crossings

E[P] - average rate of peaks

FFT - Fast Fourier Transform

[f.sub.1] - Dirlik spectral parameter

[f.sub.2] - Dirlik spectral parameter

f, [f.sub.i] - frequency

[f.sup.j] - frequency to the jth power

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] - Dirlik probability density function

j - spectral moment order

k - slope on PSD plot, special moment order for OC and SM

LaL - Lalanne method

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] - lower bound of mean stress or load

m - fatigue slope exponent

[M.sub.j] - spectral moment of order j

Max - maximum

Min - minimum

N - number of load bins

NB - narrowband

[N.sub.f] - cycles to failure at a given load

[n.sub.i] - number of cycles for given load

OC - Ortiz and Chen method

[p.sub.1] - Dirlik spectral parameter

[p.sub.2]] - Dirlik spectral parameter

P - probability of an event

PSD - power spectral density

Q - Dirlik parameter

R - Dirlik parameter

r - Dirlik spectral conversion ratio

Rect - rectangular

RFC - rainflow cycle counting

[S.sub.a] - load or stress amplitude

[S'.sub.f] - fatigue strength coefficient

SM - Larsen and Lutes Single Moment method

Std Dev - standard deviation

T - time

[UB.sub.Sa] - upper bound of stress in terms of amplitude

[UB.sub.S.sub.mean] - upper bound of mean stress

w - Zhao and Baker parameter

[W.sub.A] - load^2/hz at breakpoint "A"

[W.sub.x] - load^2/hz

WB - wideband

WL - Wirsching and Light method

[X.sub.m] - relative mean frequency

Z - load normalized by standard deviation, Dirlik parameter

ZB - Zhao and Baker method

[alpha] - Zhao and Baker parameter

[[alpha].sub.1] - BT bandwidth parameter

[[alpha].sub.2] - BT bandwidth parameter, irregularity factor

[[alpha].sub..75] - BT bandwidth parameter

[beta] - Zhao and Baker Weibull bandwidth parameter

[zeta], [[zeta].sub.i] - narrowband correction factor

[GAMMA]( ) - Gamma function

[lambda] - spectral width parameter

[PHI] - standard normal distribution function

[MATHEMATICAL EXPRESSION NOT REPRODUCIBLE IN ASCII] - standard deviation of load with respect to amplitude

[gamma] - irregularity factor

John P. Quigley

General Dynamics Land Systems

Yung-Li Lee and Liang Wang

FCA US LLC

Table 1. Tabulated results of the upper bound factor Z Values for upper fatigue slope factor m bound Z 3 4 5 6 7 8 9 Probability P 0.9 3.039 3.263 3.467 3.655 3.832 3.998 4.156 0.99 3.884 4.100 4.298 4.482 4.655 4.818 4.972 0.999 4.529 4.739 4.932 5.111 5.280 5.440 5.591 0.9999 5.074 5.278 5.466 5.642 5.807 5.964 6.113 0.99999 5.555 5.754 5.938 6.110 6.272 6.426 6.573 0.999999 5.991 6.185 6.366 6.535 6.694 6.846 6.990 0.9999999 6.392 6.583 6.760 6.926 7.083 7.233 7.375 0.99999999 6.767 6.954 7.128 7.292 7.446 7.594 7.735 0.999999999 7.120 7.304 7.475 7.636 7.788 7.934 8.073 0.9999999999 7.454 7.635 7.804 7.962 8.113 8.256 8.394 0.99999999999 7.772 7.950 8.117 8.273 8.422 8.564 8.700 0.999999999999 8.077 8.253 8.417 8.571 8.718 8.858 8.993 0.9999999999999 8.370 8.543 8.705 8.857 9.003 9.141 9.275 0.99999999999999 8.652 8.822 8.982 9.133 9.277 9.415 9.546 Values for upper fatigue slope factor m Frequency bound Z 10 11 12 1 in Probability P 0.9 4.307 4.451 4.590 10 0.99 5.120 5.262 5.398 100 0.999 5.737 5.876 6.010 1,000 0.9999 6.256 6.393 6.525 10,000 0.99999 6.714 6.849 6.980 100,000 0.999999 7.129 7.263 7.392 1,000,000 0.9999999 7.512 7.644 7.771 10,000,000 0.99999999 7.870 8.000 8.126 100,000,000 0.999999999 8.207 8.336 8.460 1,000,000,000 0.9999999999 8.526 8.654 8.777 10,000,000,000 0.99999999999 8.831 8.957 9.079 100,000,000,000 0.999999999999 9.122 9.247 9.368 1,000,000,000,000 0.9999999999999 9.403 9.527 9.647 10,000,000,000,000 0.99999999999999 9.673 9.796 9.915 100,000,000,000,000 Table 2. Standard deviation versus fatigue slope exponent m m NB WL OC BT SM [[alpha].sub..75] ZB LaL Dirlik 3 0.27 0.21 0.06 0.02 0.05 0.04 0.21 0.41 0.04 4 0.32 0.24 0.08 0.03 0.08 0.06 0.26 0.41 0.06 5 0.33 0.25 0.11 0.06 0.10 0.09 0.28 0.39 0.08 6 0.33 0.24 0.14 0.08 0.12 0.12 0.29 0.37 0.11 7 0.32 0.22 0.16 0.09 0.14 0.14 0.29 0.36 0.12 8 0.32 0.21 0.18 0.11 0.16 0.16 0.30 0.34 0.14 9 0.31 0.20 0.20 0.12 0.18 0.18 0.31 0.33 0.15 10 0.31 0.18 0.22 0.14 0.19 0.19 0.31 0.33 0.17 11 0.31 0.17 0.24 0.15 0.20 0.20 0.33 0.32 0.18 12 0.30 0.16 0.25 0.16 0.22 0.22 0.35 0.32 0.19 Table 3. Average error versus fatigue slope exponent m m NB WL OC BT SM [[alpha].sub..75] ZB LaL 3 0.36 0.14 0.00 -0.04 -0.05 -0.03 0.12 0.47 4 0.43 0.15 0.00 -0.09 -0.07 0.01 0.16 0.51 5 0.47 0.12 -0.01 -0.13 -0.09 0.04 0.19 0.53 6 0.48 0.08 -0.03 -0.17 -0.10 0.05 0.20 0.53 7 0.49 0.04 -0.04 -0.20 -0.11 0.06 0.20 0.53 8 0.49 -0.01 -0.04 -0.22 -0.12 0.06 0.21 0.53 9 0.49 -0.06 -0.05 -0.24 -0.13 0.07 0.21 0.52 10 0.49 -0.11 -0.05 -0.25 -0.13 0.07 0.22 0.52 11 0.49 -0.16 -0.05 -0.27 -0.13 0.08 0.23 0.52 12 0.50 -0.21 -0.05 -0.27 -0.14 0.08 0.25 0.52 m Dirlik 3 -0.06 4 -0.09 5 -0.12 6 -0.14 7 -0.16 8 -0.17 9 -0.18 10 -0.19 11 -0.19 12 -0.19

Printer friendly Cite/link Email Feedback | |

Author: | Quigley, John P.; Lee, Yung-Li; Wang, Liang |
---|---|

Publication: | SAE International Journal of Materials and Manufacturing |

Article Type: | Report |

Date: | Aug 1, 2016 |

Words: | 8039 |

Previous Article: | Experimental study of the plasticity responses of TRIP780 steel subjected to strain-path changes. |

Next Article: | Effect of humidity on the very high cycle fatigue behavior of a cast aluminum alloy. |

Topics: |