Abstract
Gaussian signal is produced by ordinary random vibration controllers to test the products in the laboratory, while the field data usually is nonGaussian. To synthesize nonGaussian random vibration, both the probability density function (PDF) and the damage effects must be considered. A new method is presented in this paper to synthesize nonGaussian random vibration that is characterized by running RMS (root mean square). The essential idea is to model the nonGaussian signal by a Gaussian signal multiplied by an amplitude modulation function (AMF). A twoparameter Weibull distribution is used to model the PDF of the running RMS and to create the AMF. The shock response spectrum (SRS) is used to detect significant shocks within the nonGaussian signal. A case study is presented to show that the synthesized nonGaussian signal has the same power spectral density (PSD), kurtosis, PDF and fatigue damage spectrum (FDS) as the field data.
1. Introduction
Random vibration testing is usually used to bring a test item to failure to identify weaknesses in the product or to verify if the product can survive a particular random vibration environment. Random vibration controllers accomplish this goal by producing a power spectral density (PSD) that would expose the test item to the type of vibratory environment that the test item would experience in a realworld setting. Since PSD does not contain phase information, vibration controllers always assume the phase of each frequency component follows a 0 to 2$\pi $ uniform distribution when producing time history using the PSD. This leads to data with Gaussian distribution. NonGaussian vibration is usually encountered, especially in the road transportation. W. H. Connon and D. Charles points out that when vehicles are running on irregular roads, the random excitation placed on vehicles is often nonGaussian. The tail of measured acceleration PDF is wider and the middle is sharper when compared with Gaussian distribution [1, 2]. This difference may lead to totally different accumulated fatigue damage, because most fatigue damage is caused by 2 to 4$\sigma $ [3]. When peak value exceeds 3$\sigma $, products fail quickly. MILSTD810 points out that we should always check if the field tested data is nonGaussian and the testing hardware and software is appropriate [4]. The time waveform replication (TWR) is frequently referred to as a methodology for nonGaussian testing. The basic idea of TWR is to reproduce a sequence of instantaneous values of the vibration process. Such a test may be nonGaussian, however this is only a replication of one particular measured record, not a simulation of a specified road type. Besides, MILSTD810F says that TWR is usually used for controlling transient signals or short time lasting random vibration [4]. With these drawbacks of TWR, many researchers presented some other methods to simulate nonGaussian vibration. Smallwood used three zeromemory nonlinear (ZMNL) functions to transform Gaussian vibration to nonGaussian vibration [5, 6]. These functions all give similar, but not identical results. Each covers similar but slightly different ranges of skewness and kurtosis. Each will result in slightly different distributions. Initially, this method was recommended in the form where nonGaussian time histories are first prepared offline before running a test and, then, these time histories are reproduced on the shaker in the Time Waveform Replication mode. It was only later when this method was extended and modified into a method that can be used in the closedloop frequency domain control mode [8, 9]. This method is easy to implement but it significantly reduces the controller’s dynamic range compared to the traditional Gaussian random testing and may disrupt the PSD shape. Stein wolf developed a special phase selection method [79]. The essence of this method is that the amplitudes are still responsible for the PSD, whereas the phases are not all random as they were in the classic IFFT generation. Instead, some of the phases are now used to adjust kurtosis. Thus, both characteristics, the PSD and the kurtosis, are controlled independently of each other. Rouillard et al. presented a novel technique by which nonGaussian vibrations are synthesized by generating a sequence of random Gaussian processes of varying RMS levels and durations [1015]. John et al. studied the relationship between kurtosis and fatigue damage spectrum (FDS) [16]. He declares that the fatigue damage spectrum shows the increasing fatigue damage across all frequencies that accompany the increased kurtosis levels that cause the faster product failure.
To synthesize the nonGaussian random vibration, both the PDF and the damage effects must be considered. A new method is presented in this paper to synthesize the nonGaussian random vibration that is characterized by the running RMS (root mean square). The essential idea is to model the nonGaussian signal by a Gaussian signal multiplied by a slowly varying amplitude modulation function (AMF). A twoparameter Weibull distribution is used to model the PDF of the running RMS and to create the AMF. The shock response spectrum (SRS) is used to detect significant shocks within the nonGaussian field signal. No sequence of random Gaussian processes of varying RMS levels and durations needs to be generated in this method, which makes it much easier to be implemented than that Rouillard et al. used [1015]. A case study is presented to show that the synthesized nonGaussian signal has the same PSD, kurtosis, PDF and FDS as the field data.
2. Theory
2.1. NonGaussian random signal
It is well known that the PSD is sufficient for describing the characteristics of stationary random signals only if the process under consideration is Gaussian. To check if a random process is Gaussian, the PDF of the process must be compared with the Gaussian model:
where ${\mu}_{x}$ is the mean value, ${\sigma}_{x}$ is the standard deviation.
Two principal characters describing nonGaussian PDF features are skewness and kurtosis:
where ${m}_{n}$ is the central moments of the PDF:
Eq. (3) is the theoretical definitions of the central moments. For a given digitized data record, it can be calculated by time averaging:
For a truly Gaussian process, the skewness is 0 and the kurtosis is 3. Any deviation from these indicates that the process is nonGaussian. From the perspective of simulating nonGaussian vibration, kurtosis is a more important parameter than skewness, because it represents the probability of peak values in time history. The PDF of data produced by vibration testing can be more realistic using this parameter. The kurtosis describes the sharpness or flatness in the middle range of PDF and the narrowness or wideness in the tail of PDF. Kurtosis larger than 3 indicating the middle range of PDF of field tested signal is sharper than a Gaussian signal with the same standard deviation, and the tail is wider. Wider tail means that the probability of peak values is larger, which leads to larger fatigue damage and quicker fail.
2.2. Fatigue damage analysis
The Rain Flow Cycle Counting (RFCC) and Dirlik methods are frequently referred to as two methodologies for fatigue damage analysis. The RFCC method was initially proposed by M. Matsuiski [17] and T. Endo [18] to count the cycles or the halfcycles of straintime signals. For fatigue damage analysis under vibration, the cyclic numbers ${n}_{RFCC}^{i}$ at $i\text{th}$ stress level ($i=$ 1, 2,..., $p$) is counted directly from the stress time history. The number of stress levels $p$ depends on how many bins are used to perform the RFCC. See the counting algorithm in [19] for detail.
T. Dirlik [20] established empirical expressions of the PDF of the ordinary halfranges and those counted with the RFCC method using a digital simulation. Let $x\left(t\right)$ be the response of a SDOF system to a Gaussian signal, ${x}_{rms}$ the rms value of $x\left(t\right)$, $\delta $ a range defined starting from RFCC and $u$ the reduced variable corresponding to a halfcycle:
The PDF of $u$ is expressed as:
where:
Using the given PSD, the mean number of positive maxima per second can be expressed as [21]:
With Eqs. (6) and (10), the cyclic numbers ${n}_{Dirlik}^{i}$ at $i\text{th}$ stress level can be expressed as:
where $T$ is the time duration of the vibration environment that the product experienced, in seconds.
At each stress level, a comparion between ${n}_{Dirlik}^{i}$ and ${n}_{RFCC}^{i}$ for both Gaussian and nonGaussian vibrations will be made later, to show the limitation of the Dirlik method when dealing with nonGaussian signal.
2.3. Fatigue damage spectrum and shock response spectrum
All structures experience fatigue as they are repeatedly exposed to adequate levels of stress. The total damage a product experiences in a particular time period can be calculated from field data and plotted for a specific range of frequencies. The FDS is essentially a plot that shows the responses of a number of single degree of freedom (SDOF) systems to a base input acceleration time history. Many SDOF systems tuned to a range of natural frequencies are assessed using the same input. The final plot, the FDS, shows the fatigue damage encountered for a particular SDOF system anywhere within the analyzed time, as shown in Fig. 1 and Fig. 2.
Fig. 1SDOF systems for calculation of FDS
Fig. 2Calculation procedure of the FDS
For a SDOF system with a natural frequency ${f}_{n}$ and a damping ratio $\xi $, the output relative displacement ${x}_{d}$ to an input acceleration ${x}_{a}$ can be calculated using a ramp invariant digital filter method [22]. With the output ${x}_{d}$, the cumulative damage can be calculated in both time domain and frequency domain. In time domain, RFCC method is usually used to count cyclic numbers at different stress levels. Then SN curve [23] and Miner’s rule [24] are combined to calculate fatigue damage:
where ${n}_{i}$ is the number of cycles exposure at ${S}_{i}$ ($i=\text{1, 2,...,}\text{}p$). $p$ is the number of stress levels considered. ${N}_{i}$ is fatigue life at stress ${S}_{i}$, $b$ is the fatigue exponent, ${D}_{t}$ is the total damage index calculated in time domain.
The stress is assumed to be proportional to relative displacement [25]:
The total damage can be calculated as:
In frequency domain, Rayleigh distribution of response stress maxima is assumed and used to calculate the stress cycles [26]:
where $S$ is the stress value of peaks, ${\sigma}_{S}$ is the RMS of the stress time history.
The total damage can be calculated as:
where ${n}_{p}^{+}$ can be replaced by the natural frequency ${f}_{n}$ for a lightly damped SDOF system response, $T$ is the total time of exposure to the stress environment.
Solving Eq. (17) with Eq. (16) leads to:
where $\mathrm{\Gamma}$ is the Gamma function, ${\sigma}_{d}$ is the RMS of relative displacement.
The RMS of relative displacement can be calculated using:
where $H\left(f\right)$ is the transmissibility of a SDOF system with a natural frequency of ${f}_{n}$ and a damping ratio of $\xi $.
Since different kurtosis leads to different number of peak accelerations in the test, and these peak accelerations are primarily responsible for damage that a product experiences in the field, the FDS should show the difference when calculated for Gaussian and nonGaussian vibrations.
The effects of the shock loads are generally depicted in a SRS. Like the FDS, the SRS is essentially a plot that shows the responses of a number of SDOF systems to a base input acceleration time history, as shown in Fig. 1. The only difference is that the SRS is generated by calculating the maximum response of a SDOF system to the input. The final plot, the SRS, shows the largest response encountered for a particular SDOF system anywhere within the analyzed time.
Since the essential idea of the new synthesis method is to model the nonGaussian signal by a Gaussian signal multiplied by an AMF:
significant shocks in the field data must be first identified and separated. To detect if there are significant shocks in the field data, the response of a number of SDOF systems to the field data is firstly calculated in time domain. Then the response is divided into $M$ consecutive segments and the maximum value is calculated for each segment. Finally the ratio between this maximum value and the SRS computed from the PSD, which is calculated from the whole time history using Welch method, is used to show where the significant shock is:
where ${y}_{n,i}$ indicates the maximum value calculated for the $i$th segment and the $n$th SDOF system, $SR{S}_{n}$ indicates the SRS calculated form PSD for the $n$th SDOF system.
3. Synthesis method
A new method is presented to synthesize the nonGaussian random vibration characterized by the running RMS. The procedure uses three phases to streamline the process:
Phase 1: Check the characteristics of measured vibration signal.
1. Obtain the environment acceleration time history, resample if needed.
2. Calculate the kurtosis of measured signal to see if it is NonGaussian or not.
3. Calculate the PSD of measured signal using Welch’s method.
4. Use SRS to detect significant shocks (see Eq. (21)).
Phase 2: synthesize the NonGaussian random vibration.
1. Create the stationary Gaussian signal using the PSD calculated in step 3, Phase 1.
2. Calculate the probability distribution function ${F}_{RMS}$ of the running RMS of the field data.
3. Use a two parameter Weibull distribution to model the probability distribution of the running RMS.
The cumulative distribution function of Weibull distribution can be expressed as:
where $x$ is the random variable, $L$ is the scale parameter, $K$ is the shape parameter.
To find the correct $L$ and $K$, the probability distribution function of the assumed Weibull distribution ${F}_{Weibull}$ and the running RMS ${F}_{RMS}$ are used to form an objective function using least squares method:
where $x$ are the values of the running RMS, depending on how many bins are used to calculated the probability distribution function of the running RMS, $i$ indicates the $i$th iteration, $std$ indicates standard deviation. Note that both ${F}_{RMS}$ and ${F}_{Weibull}$ are vectors.
The MATLAB function “fminsearch” is used find the right $L$ and $K$, so that the standard deviation of error between these two probability distribution functions is minimized.
1. Use the determined $L$ and $K$ to create the AMF and to get the required NonGaussian signal.
2. To create the appropriate AMF, the kurtosis of the field data is used as an objective value. The AMF is created from the Weibull distribution function and multiplied with the Gaussian signal created in step 1, Phase 2, to get the synthesized signal. Compare the kurtosis of the synthesized signal with the objective value. Do iterations until the objective value is reached. The synthesized nonGaussian signal is scaled to have the same RMS as the field data.
Phase 3: compare the characteristics of the synthesized NonGaussian signal and the field data.
1. Compare the PSD.
2. Compare the PDF.
3. Compare the fatigue, e.g., cyclic numbers and FDS.
4. Case study
In this case study, a new method for synthesizing the nonGaussian random vibration is developed and examined by the field data. The field data provided by Prof. Kjell Ahlin is from an Ericsson Mast Project and is shown in Fig. 3. The test item and setup are shown in Fig. 4.
Fig. 3Field data
Fig. 4Test item and setup
This is a nonGaussian signal with kurtosis equals to 9.4262. Eq. (21) is used to detect significant shocks (see Fig. 5). The ratio is shown in the third dimension by the color. As we can see from Fig. 5, there are only some insignificant shocks superimposed on the underlying random vibration signal.
The PSD of the field data is calculated using Welch method with 4096 points hanging window and 50 % overlapping. Gaussian signal is generated using the same but smoothed PSD and is shown in Fig. 6.
Fig. 5Detection of shocks
Fig. 6Synthesized Gaussian signal using the same PSD calculated from field data
The PSD and PDF of both the field data and the synthesized Gaussian signal are shown in Fig. 7 and Fig. 8. As we can see, the PDF of the field data is obviously different from that of the synthesized Gaussian signal, although the PSD are basically the same.
Fig. 7PSD of field data and synthesized Gaussian signal
Fig. 8PDF of the field data and the synthesized Gaussian signal
Fig. 9The running RMS
Fig. 10PDF of the running RMS using different number of bins and Rayleigh distribution
Fig. 11The iteration to determine L and K
The running RMS and its PDF using different number of bins are calculated and shown in Fig. 9 and Fig. 10. The running RMS is calculated using 10 s frames with 50 % overlapping. From Fig. 10 we can see that, despite of the numbers of bins, the PDF of the running RMS is similar to the PDF of a given Rayleigh distribution. To be more general, a twoparameter Weibull distribution is used to model the probability distribution function of the running RMS Eq. (22).
Iteration method described in Section 3 (step 3, phase 2) is used to determine $L$ and $K$. The iteration process is shown in Fig. 11. The probability distribution function of resulted Weibull distribution is shown in Fig. 12, together with the probability distribution function of the running RMS.
As the parameters of the Weibull distribution are determined, the AMF can be created and the nonGaussian signal can be synthesized using the iteration method presented in Section 3 (step 4, phase 2). The synthesized nonGaussian signal is shown in Fig. 13, together with the field data (the time duration of the synthesized nonGaussian signal is adjustable). Note that some care must be taken to have the time scale for the AMF selected properly, since inappropriate selection of time scale may lead to inaccurate PDF and/or PSD.
The PSD and PDF of the field data and the synthesized nonGaussian signal are very similar (see Fig. 14 and 15), which proves the effectiveness of the new method.
For comparison, the cycle counting with the RFCC and Dirlik method can be performed directly on the field data (acceleration signal), the synthesized Gaussian signal and nonGaussian signal. The results are shown in Fig. 16. As we can see, the Dirlik method is equivalent to the RFCC for Gaussian signal. Big error can be introduced by the Dirlik method for nonGaussian signal. The RFCC results for the field data and synthesized nonGaussian signal are very similar.
Fig. 12The probability distribution function of the running RMS and Weibull distribution
Fig. 13Field data and synthesized nonGaussian signal
Fig. 14Mast vibration signal power spectral density. The PSD of the field data and the synthesized nonGaussian signal
Fig. 17 shows the FDS calculated for the field data, the synthesized Gaussian and nonGaussian signal in the time domain using Eq. (15). For comparison, $k$ and $c$ are normalized to 1. From Fig. 17 we can see that the FDS of the field data and the synthesized nonGaussian signal are very similar. This again proves the effectiveness of the new method. Besides, we can also see that larger kurtosis leads to larger fatigue damage.
Fig. 15Amplitude probability density. PDF for field data and synthesized Gaussian and nonGaussian signal
Fig. 16Number of cycles versus data ranges for RFCC and Dirlik
Fig. 17Fatigue damage spectrum. The FDS of the field data, synthesized Gaussian and nonGaussian signal
5. Conclusion
A new method for synthesizing the running RMSinduced nonGaussian Random Vibration is developed and examined by field data. Conclusions are as below:
1. This method is limited to the case of the running RMSinduced nonGaussian signals, which contain no significant shocks.
2. The synthesized nonGaussian signal has the same PSD, kurtosis, PDF and FDS with the field data, which proves the effectiveness of the new method.
3. The Dirlik methods are equivalent to the RFCC method for Gaussian signal and can introduce big error for nonGaussian signal.
4. Larger kurtosis leads to larger fatigue damage.
5. Weibull distribution can be used to simulate the statistical distribution of the running RMS and to synthesize the running RMSinduced nonGaussian signal.
References

Connon W. H. Comments on kurtosis of military vehicle vibration data. Journal of the IES, Vol. 34, Issue 6, 1991, p. 3841.

Charles D. Derivation of environment descriptions and test severities from measured road transportation data. Journal of the IES, Vol. 36, Issue 1, 1993, p. 3742.

Lambert R. G. Fatigue life prediction for various random stress peak distributions. Shock and Vibration Inform. Center The Shock and Vibration Bulletin, Vol. 52, Issue 4, 1982, p. 110.

MILSTD810F. Department of Defense Test Method Standard for Environmental Engineering Considerations and Laboratory Tests, USA, 2000.

Smallwood D. O. Generation of stationary nonGaussian time histories with a specified crossspectral density. Shock and Vibration, Vol. 4, Issue 5, 1997, p. 361377.

Smallwood D. O. Generating nonGaussian vibration for testing purposes. Sound and Vibration, Vol. 39, Issue 10, 2005, p. 1823.

Steinwolf A. Shaker simulation of random vibration with a high kurtosis value. Journal of the Institute of Environmental Sciences, Vol. 40, Issue 3, 1997, p. 3343.

Steinwolf A. Closedloop shaker simulation of nonGaussian random vibrations. Part 1. Discussion and methods. Test Engineering and Management, Vol. 68, Issue 3, 2006, p. 1013. Part 2. Numerical and experimental results. Test Engineering and Management, Vol. 68, Issue 5, 2006, p. 1419.

Steinwolf A. Random vibration testing beyond PSD limitations. Sound and Vibration, Vol. 40, Issue 9, 2006, p. 1221.

Rouillard V. On the Synthesis of NonGaussian Road Vehicle Vibrations. Monash University, 2007.

Rouillard V., Lamb M. On the effects of sampling parameters when surveying distribution vibrations. Packaging Technology and Science, Vol. 431, Issue 8, 2008, p. 467477.

Rouillard V., Sek M. A. Synthesizing nonstationary, nonGaussian random vibrations. Packaging Technology and Science, Vol. 23, Issue 8, 2010, p. 423439.

GarciaRomeu M. A., Rouillard V. On the statistical distribution of road vehicle vibrations. Journal of Packaging Technology and Science, Vol. 24, Issue 8, 2011, p. 451467.

Rouillard V., Sek M. Creating transport vibration simulation profiles from vehicle and road characteristics. Packaging Technology and Science, Vol. 26, Issue 2, 2013, p. 8295.

Rouillard V. Quantifying the nonstationarity of vehicle vibrations with the run test. Packaging Technology and Science, Vol. 27, Issue 3, 2014, p. 203219.

Van Baren J., Van Baren P., Jenison M. I. The fatigue damage spectrum and kurtosis control. Sound and Vibration, Vol. 46, Issue 10, 2012, p. 1013.

Matsuishi M., Endo T. Fatigue of metals subjected to varying stress. Japan Society of Mechanical Engineers, Fukuoka, Japan, 1968, p. 3740.

Endo T. Damage evaluation of metals for random on varying loadingthree aspects of rain flow method. Proceedings of Symposium on MBM, Vol. 1, 1974, p. 374.

Lalanne C. Mechanical Vibration and Shock, Fatigue Damage. Wiley, 2010, p. 111120.

Lalanne C. Mechanical Vibration and Shock, Fatigue Damage. Wiley, 2010, p. 172174.

Lalanne C. Mechanical Vibration and Shock, Random Vibration. Wiley, 2010, p. 281.

Ahlin K. Comparison of test specifications and measured field data. Sound and Vibration, Vol. 40, Issue 9, 2006, p. 2225.

Lalanne C. Mechanical Vibration and Shock, Fatigue Damage. Wiley, 2010, p. 1013.

Lalanne C. Mechanical Vibration and Shock, Fatigue Damage. Wiley, 2010, p. 4950.

Lalanne C. Mechanical Vibration and Shock, Specification Development. Wiley, 2010, p. 76.

Papoulis A. Narrowband systems and Gaussianity. IEEE Transactions on Information Theory, Vol. 18, Issue 1, 1972, p. 2027.