An Analytical Solution for Predicting the Vibration-Fatigue-Life in Bimodal Random Processes

Abstract

Predicting the vibration-fatigue-life of engineering structures subjected to random loading is a critical issue for. Frequency methods are generally adopted to deal with this problem. This paper focuses on bimodal spectra methods, including Jiao-Moan method, Fu-Cebon method, and Modified Fu-Cebon method. It has been proven that these three methods can give acceptable fatigue damage results. However, these three bimodal methods do not have analytical solutions. Jiao-Moan method uses an approximate solution, Fu-Cebon method, and Modified Fu-Cebon method needed to be calculated by numerical integration which is obviously not convenient in engineering application. Thus, an analytical solution for predicting the vibration-fatigue-life in bimodal spectra is developed. The accuracy of the analytical solution is compared with numerical integration. The results show that a very good agreement between an analytical solution and numerical integration can be obtained. Finally, case study in offshore structures is conducted and a bandwidth correction factor is computed through using the proposed analytical solution.

Chaoshuai Han 1 and Yongliang Ma 1 and Xianqiang Qu 1 and Mindong Yang 2

Academic Editor:Nuno M. Maia

1, College of Shipbuilding Engineering, Harbin Engineering University, Harbin, Heilongjiang 150001, China

2, Guangdong Electric Power Design Institute Co. Ltd. of China Energy Engineering Group, Guangzhou, China

Received 6 September 2016; Accepted 26 December 2016; 31 January 2017

This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

1. Introduction

Engineering structures from different fields (e.g., aircrafts, wind energy utilizations, and automobiles) are commonly subjected to random vibration loading. These loads often cause structural fatigue failure. Thus, it is significant to carry out a study on assessing the vibration-fatigue-life [1, 2].

Vibration fatigue analysis commonly consists of two part processes: structural dynamic analysis and results postprocessing. Structural dynamic analysis provides an accurate prediction of the stress responses of fatigue hot-spots. Once the stress responses are obtained, vibration fatigue can be successfully performed. Existing technologies such as operational modal analysis [3], finite element modeling (FEM), and accelerated-vibration-tests are mature and applicable to obtain the stress of structures [4, 5]. Therefore, the crucial part of a vibration fatigue analysis focuses on results postprocessing.

The postprocessing is usually used to calculate fatigue damage based on known stress responses. When the stress responses are time series, fatigue can be evaluated using a traditional time domain method. However, the stress responses of real structures are mostly characterized by the power spectral density (PSD) function. Thus, frequency domain method becomes popular in vibration fatigue analysis [6, 7].

A bimodal spectrum is a particular PSD in the random vibration stress response of a structure. For some simple structures, the stress response of structures will show explicit characterization of two peaks. One peak of the bimodal spectrum is governed by the first-order natural frequency of the structure; another is dominated by the main frequency of applied loads. Therefore, some bimodal methods for fatigue analysis can be adopted [8-10]. Moreover, several experiments (e.g., vibration tests on mechanical components) and numerical studies (e.g., virtual simulation of dynamic using FEM) also obtain have shown that the stress PSD is a typical bimodal spectrum [5, 7, 11]. However, for some complex flexible structures, the PSD of the stress response of structures usually is a multimodal and wide-band spectrum. For this situation, existing general wide-band spectral methods such as Dirlik method [12], Benasciutti-Tovo method [10], and Park method [13, 14] can be used to evaluate the vibration-fatigue-life. Recently, Braccesi et al. [15, 16] proposed a bands method to estimate the fatigue damage of a wide-band random process in the frequency domain. In order to speed up the frequency domain method, Braccesi et al. [17, 18] developed a modal approach for fatigue damage evaluation of flexible components by FEM.

For fatigue evaluation in bimodal processes, some specific formulae have been proposed. Jiao and Moan [8] provided a bandwidth correction factor from a probabilistic point of view, and the factor is an approximate solution derived by the original model. The approximation inevitably leads to some errors in certain cases. Based on an similar idea, Fu and Cebon [9] developed a formula for predicting the fatigue life in bimodal random processes. In the formula, there is a convolution integral. The author claimed that there is no analytical solution for the convolution integral which has to be calculated by numerical integration. Benasciutti and Tovo [10] compared the above two methods and established a Modified Fu-Cebon method. The new formula improves the damage estimation, but it still needs to calculate numerical integration.

In engineering application, the designers prefer an analytical solution rather than numerical methods. Therefore, the purpose of this paper is to develop an analytical solution to predict the vibration-fatigue-life in bimodal spectra.

2. Theory of Fatigue Analysis

2.1. Fatigue Analysis

The basic S-N curve for fatigue analysis can be given as [figure omitted; refer to PDF] where S represents the stress amplitude; N is the number of cycles to fatigue failure, and K and m are the fatigue strength coefficient and fatigue strength exponent, respectively.

The total fatigue damage can then be calculated as a linear accumulation rule after Miner [19] [figure omitted; refer to PDF] where _ni is the number of cycles in the stress amplitude _Si , resulting from rainflow counting (RFC) [20], and _Ni is the number of cycles corresponding to fatigue failure at the same stress amplitude.

When the stress amplitude is a continuum function and its probability density function (PDF) is _fS (S), the fatigue damage in the duration T can be denoted as follows: [figure omitted; refer to PDF] where _νc is the frequency of rainflow cycles.

For an ideal narrowband process, _fS (S) can be approximated by the Rayleigh distribution [21]; the analytical expression is given as [figure omitted; refer to PDF] Furthermore, the frequency of rainflow cycles _νc can be replaced by rate of mean zero upcrossing _ν0 .

According to (3), an analytical solution of fatigue damage [22] for an ideal narrowband process can be written as [figure omitted; refer to PDF] where Γ(·) is the Gamma function.

For general wide-band stress process, fatigue damage can be calculated by a narrowband approximation (i.e., (5)) first, and bandwidth correction is made based on the following model [23]: [figure omitted; refer to PDF]

In general, bimodal process is a wide-band process; thus, the fatigue damage in bimodal process can be calculated through (6).

2.2. Basic Principle of Bimodal Spectrum Process

Assume that a bimodal stress process X(t) is composed of a low frequency process (LF) _XL (t) and a high frequency process (HF) _XH (t). [figure omitted; refer to PDF] where _XL (t) and _XH (t) are independent and narrow Gaussian process.

The one-sided spectral density function of X(t) can be summed from the PSD of LF and HF process. [figure omitted; refer to PDF]

The ith-order spectral moments of S(ω) are defined as [figure omitted; refer to PDF]

The rate of mean zero upcrossing corresponding to _XL (t) and _XH (t) is [figure omitted; refer to PDF]

The rate of mean zero upcrossing of X(t) can be expressed as [figure omitted; refer to PDF]

According to (5), (9), and (11), narrowband approximation of bimodal stress process X(t) can be given as [figure omitted; refer to PDF] Equation (12) is known as the combined spectrum method in API specifications [24].

The existing bimodal methods proposed by Jiao and Moan, Fu and Cebon, and Benasciutti and Tovo are based on the idea: two types of cycles can be extracted from the rainflow counting, one is the large stress cycle, and the other is the small cycle [8-10]. The fatigue damage due to X(t) can be approximated with the sum of two individual contributions. [figure omitted; refer to PDF] where _Dl represents the damage due to the large stress cycle and _Ds denotes the damage due to the small stress cycle.

3. A Review of Bimodal Methods

3.1. Jiao-Moan (JM) Method

To simplify the study, X(t), _XL (t), and _XH (t) are normalized as ^{X[low *]} (t), ^{XL[low *]} (t), and ^{XH[low *]} (t) through the following transformation: [figure omitted; refer to PDF] and then [figure omitted; refer to PDF] where [figure omitted; refer to PDF]

Jiao-Moan points out that the small stress cycles are produced by the envelope of the HF process, which follows the Rayleigh distribution. The fatigue damage due to the small stress cycles can be obtained according to (5).

While the large stress cycles are from the envelop process, P(t) (see Figure 1), the amplitude of P(t) is equal to [figure omitted; refer to PDF] where _RL (t) and _RH (t) are the envelopes of ^{XL[low *]} (t) and ^{XH[low *]} (t), respectively.

Figure 1: Bimodal process ^{X[low *]} (t), the envelope process P(t), and the amplitude process Q(t).

[figure omitted; refer to PDF]

The distribution of Q(t) can be written as a form of a convolution integral [figure omitted; refer to PDF] _RL (t) and _RH (t) obey the Rayleigh distribution; therefore, (18) has an analytical solution which is given [8] [figure omitted; refer to PDF]

The rate of mean zero upcrossing due to P(t) can be calculated as [figure omitted; refer to PDF] where [figure omitted; refer to PDF]

An approximation was made by Jiao and Moan for (19) as follows [8]: [figure omitted; refer to PDF]

After the approximation, a closed-form solution of the bandwidth correction factor can be then derived [8] [figure omitted; refer to PDF]

Finally, the fatigue damage can be obtained as (6) and (12).

3.2. Fu-Cebon (FC) Method

Similarly to JM method, Fu and Cebon also considered that the total damage is produced by a large cycle (_SH +_SL ) and a small cycle (_SL ), as depicted in Figure 2. The small cycles are from the HF process, and the distribution of the amplitude _PSs (S) is a Rayleigh distribution, as shown in (4). However, the number of cycles associated with the small cycles _ns is different from JM method and equals (_ν0,H -_ν0,L )·T. According to (5), the damage due to the small cycles is [figure omitted; refer to PDF]

Figure 2: The large cycles and small cycles for a random stress process.

[figure omitted; refer to PDF]

The amplitude of the large cycles _Sl can be approximated as the sum of amplitude of the LF and HF processes, the distribution of which can be expressed by a convolution of two Rayleigh distributions [9]. [figure omitted; refer to PDF] where U=1/2_λ0,L +1/2_λ0,H and V=1/_λ0,H .

The number of cycles of the large cycles is _nl =_ν0,L ·T. Thus, the fatigue damage due to the large stress cycles can be expressed by [figure omitted; refer to PDF]

Equation (26) can be calculated with numerical integration [9, 10]. Therefore, the total damage can be obtained according to (13).

3.3. Modify Fu-Cebon (MFC) Method

Benasciutti and Tovo made a comparison between JM method and FC method and concluded that using the envelop process is more suitable [10]. Thus, a hybrid technique is adopted to modify the FC method. More specifically, the large cycles and small cycles are produced according to the idea of FC method. The number of cycles associated with the large cycles is defined similarly to JM method. That is, _nl =_ν0,P ·T, while the number of cycles corresponding to the small cycles is _ns =(_ν0,H -_ν0,P )·T. The total damage for MFC method can be then written according to (13).

Although the accuracy of the MFC method is improved, the fatigue damage still has to be calculated with numerical integral.

3.4. Comparison of Three Bimodal Methods

Detailed comparison of the aforementioned three bimodal methods can be found in Table 1. In all methods, the amplitude of the small cycle obeys Rayleigh distribution, and the corresponding fatigue damage has an analytical expression as in (5); the distribution of amplitude of the large cycle is convolution integration of two Rayleigh distributions, and the relevant fatigue damage can be calculated by (26).

Table 1: Comparison of large cycles and small cycles for different spectral methods.

Method	Large cycles	Small cycles
n l	PDF of amplitude	n s	PDF of amplitude
JM	ν 0 , P [sm middot] T	Eq. (22)	ν 0 , H [sm middot] T	Rayleigh
FC	ν 0 , L [sm middot] T	Eq. (25)	( ν 0 , H - ν 0 , L ) [sm middot] T	Rayleigh
MFC	ν 0 , P [sm middot] T	Eq. (25)	( ν 0 , H - ν 0 , P ) [sm middot] T	Rayleigh

Because of complexity of the convolution integration, several researches assert that (26) has no analytical solution [9, 10]. To solve this problem, Jiao and Moan used an approximate model (i.e., (22)) to obtain a closed-form solution [8]. However, the approximate model may lead to errors in some cases as in Figure 3 which illustrates the divergence of (19) and (22) for different values of ^{λ0,L[low *]} and ^{λ0,H[low *]} . It is found that (22) becomes closer to (19) with the increase of ^{λ0,L[low *]} .

Figure 3: The divergence of (19) and (22) for JM method in case of (a) ^{λ0,L[low *]} =0.01 and (b) ^{λ0,L[low *]} =0.99.

(a) [figure omitted; refer to PDF]

(b) [figure omitted; refer to PDF]

For FC and MFC methods, (26) was calculated by numerical technique. Although the numerical technique can give a fatigue damage prediction, it is complex and not convenient when applied in real engineering. In addition, the solutions in some cases are not reasonable. In Section 4, an analytical solution of (26) will be derived to evaluate the fatigue damage, and the derivation of the analytical solution focuses on the fatigue damage of the large cycles.

4. Derivation of an Analytical Solution

4.1. Derivation of an Analytical Solution for (25)

Equation (25) can be rewritten as [figure omitted; refer to PDF]

Equation (27) will be divided into two items.

( 1) The First Item. It is as follows: [figure omitted; refer to PDF]

( 2) The Second Item. It is as follows: [figure omitted; refer to PDF]

The analytical solution of (25) can be then obtained [figure omitted; refer to PDF]

Note that when _λ0,L +_λ0,H =1, (30) is just equal to (19) derived by Jiao and Moan [8]. Therefore, (19) is a special case of (30).

4.2. Derivation of an Analytical Solution for (26) Based on (30)

The derivation of an analytical solution for (26) is on the basis of (30), as [figure omitted; refer to PDF]

Equation (31) will be divided into five parts.

( 1) The First Part. It is as follows: [figure omitted; refer to PDF]

( 2) The Second Part. It is as follows: [figure omitted; refer to PDF]

( 3) The Third Part. It is as follows: [figure omitted; refer to PDF]

( 4) The Fourth Part. It is as follows: [figure omitted; refer to PDF]

Equation (35) contains a standard Normal cumulative distribution function Φ(·). It is difficult to get an exact solution directly. Thus, a new variable is introduced, as follows: [figure omitted; refer to PDF] With a method of variable substitution, (35) can be simplified: [figure omitted; refer to PDF] By defining [figure omitted; refer to PDF] (37) becomes [figure omitted; refer to PDF] and using integration by parts, (38) reduces to [figure omitted; refer to PDF]

Equation (40) is a recurrence formula; when m is an odd number, it becomes [figure omitted; refer to PDF] where (·)!! is a double factorial function and H(_λ0,L ,_λ0,H ,1) has an analytical expression which can be derived conveniently [figure omitted; refer to PDF] When m is an even number, H(_λ0,L ,_λ0,H ,m) is [figure omitted; refer to PDF] where [figure omitted; refer to PDF] Specific derivation of (44) can be seen in Appendix A.

In addition, a Matlab program has been written to calculate H(_λ0,L ,_λ0,H ,m) in Appendix B.

( 5) The Fifth Part. It is as follows: [figure omitted; refer to PDF]

Similarly to the fourth part, the analytical solution of the fifth part can be derived as [figure omitted; refer to PDF] The final solution is [figure omitted; refer to PDF]

5. Numerical Validation

In this part, the accuracy of FC method and JM method and the derived analytical solution will be validated with numerical integration. Transformation of (26) will be carried out first.

5.1. Treatment of Double Integral Based on FC Method

As pointed out by Fu and Cebon and Benasciutti and Tovo [9, 10], FC's numerical integration can be calculated as the following processes.

Equation (26) contains a double integral, in which variables S and y are in the range of (0,∞) and (0,S), respectively. Apparently, the latter is not compatible with the integration range of Gauss-Legendre quadrature formula. Therefore, by using a integration transformation [figure omitted; refer to PDF]

the integral part of (26) can be simplified to [figure omitted; refer to PDF]

Equation (49) can be thus calculated with Gauss-Legendre and Gauss-Laguerre quadrature formula.

5.2. Treatment of Numerical Integral for (31)

Direct calculation of (31) may lead to some mathematical accumulative errors. To obtain a precise integral solution, (31) has to be handled with a variable substitution, and the result is [figure omitted; refer to PDF] Note that _Z1 and _Z2 for (31) have analytical solution; therefore, only _Z3 , _Z4 , and _Z5 are dealt with in (50).

The solution of (50) can be obtained through Gauss-Laguerre quadrature formula.

The accuracy of the numerical integral in (49) and (50) depends on the order of nodes and weights which can be obtained from handbook of mathematics [25]. The accuracy increases with the increasing orders. However, from the engineering point of view, too many orders of nodes and weights will lead to difficulty in calculation. The integral results of (50) are convergent when the orders of nodes and weights are equal to 30. Therefore, the present study takes the order of 30.

5.3. Discussion of Results

It is very convenient to use JM method in the case of _λ0,L +_λ0,H =1, while FC method and MFC method can be used for any case regardless of _λ0,L and _λ0,H . Therefore, the analytical results are divided into two: _λ0,L +_λ0,H =1 and _λ0,L +_λ0,H ≠1.

The solutions calculated by different mathematical methods for m=3 and m=4 in the case of _λ0,L +_λ0,H =1 are plotted in Figures 4 and 5. It turns out that the proposed analytical solution gives the same results with the exact numerical integral solution for any value of _λ0,L and _λ0,H . FC's numerical integral solution matches the exact numerical integral solution very well in most cases. However, for relatively low value of _λ0,L or when _λ0,L tends to 1, it may not be right. JM's approximate solution is close to the exact numerical integral solution only in a small range.

Figure 4: Comparison of different solutions for m=3.

[figure omitted; refer to PDF]

Figure 5: Comparison of different solutions for m=4.

[figure omitted; refer to PDF]

FC's numerical integral solutions, the exact numerical integral solutions, and the analytical solutions for m=3 and m=4 in the case of _λ0,L +_λ0,H ≠1 are shown in Tables 2 and 3. The results indicate that the latter two solutions are always approximately equal for any value of _λ0,L and _λ0,H , while FC's numerical integral solutions show good agreement with the exact integral solutions only in a few cases. As like the case of _λ0,L +_λ0,H =1, for relatively high or low values of _λ0,L and _λ0,H (i.e., _λ0,L =100, _λ0,H =1×1^0-5 ), FC's integral solutions may give incorrect results. This phenomenon is in accord with the analysis of Benasciutti and Tovo (e.g., FC's numerical integration may be impossible for too low values of _λ0,L and _λ0,H ) [10].

Table 2: Comparison of different solutions for m=3 in the case of _λ0,L +_λ0,H ≠1.

λ 0 , L	λ 0 , H	FC's numerical integral solution for (49)	The exact numerical integral solution for (50)	The analytical solution for (47)
1.000 E - 05	9.000 E - 05	3.657 E - 12	6.183 E - 06	6.183 E - 06
3.000 E - 03	7.000 E - 03	3.292 E - 03	7.591 E - 03	7.591 E - 03
5.000 E - 02	5.000 E - 02	2.591 E - 01	2.522 E - 01	2.522 E - 01
2.000 E - 01	6.000 E - 01	5.262 E + 00	5.267 E + 00	5.267 E + 00
5.000 E + 00	4.000 E + 00	2.146 E + 02	2.146 E + 02	2.146 E + 02
8.000 E + 01	2.000 E + 01	7.062 E + 03	7.062 E + 03	7.062 E + 03
4.000 E + 02	6.000 E + 02	1.144 E + 05	2.493 E + 05	2.493 E + 05
5.000 E + 03	3.000 E + 03	1.238 E + 04	5.602 E + 06	5.602 E + 06
6.000 E + 07	5.000 E + 07	8.468 E - 05	9.180 E + 12	9.180 E + 12
7.000 E + 09	4.000 E + 09	9.073 E - 09	8.999 E + 15	8.999 E + 15
1.000 E + 02	1.000 E - 05	3.084 E - 01	3.762 E + 03	3.762 E + 03

Table 3: Comparison of different solutions for m=4 in the case of _λ0,L +_λ0,H ≠1.

λ 0 , L	λ 0 , H	FC's numerical integral solution for (49)	The exact numerical integral solution for (50)	The analytical solution for (47)
1.000 E - 05	9.000 E - 05	2.580 E - 13	1.437 E - 07	1.437 E - 07
3.000 E - 03	7.000 E - 03	1.187 E - 03	1.832 E - 03	1.832 E - 03
5.000 E - 02	5.000 E - 02	2.198 E - 01	1.942 E - 01	1.942 E - 01
2.000 E - 01	6.000 E - 01	1.131 E + 01	1.130 E + 01	1.130 E + 01
5.000 E + 00	4.000 E + 00	1.567 E + 03	1.567 E + 03	1.567 E + 03
8.000 E + 01	2.000 E + 01	1.682 E + 05	1.682 E + 05	1.682 E + 05
4.000 E + 02	6.000 E + 02	6.632 E + 06	1.915 E + 07	1.915 E + 07
5.000 E + 03	3.000 E + 03	7.706 E + 05	1.216 E + 09	1.216 E + 09
6.000 E + 07	5.000 E + 07	5.308 E - 03	2.344 E + 17	2.344 E + 17
7.000 E + 09	4.000 E + 09	5.688 E - 07	2.289 E + 21	2.289 E + 21
1.000 E + 02	1.000 E - 05	7.230 E - 01	8.006 E + 04	8.006 E + 04

In a word, for any value of _λ0,L and _λ0,H , the analytical solution derived in this paper always gives an accurate result. Furthermore, it can be solved very conveniently and quickly with the aid of a personal computer through a program given in Appendix B.

6. Case Study

In this section, the bandwidth correction factor is used to compare different bimodal spectral methods.

A general analytical solution (GAS) of fatigue damage for JM, FC, and MFC method can be written as [figure omitted; refer to PDF]

Z can be obtained according to (47), and _nl and _ns can be chosen as defined in Table 1, which represent different bimodal spectra methods, that is, JM, FC, and MFC method.

According to (6), the analytical solution of the bandwidth correction factor is [figure omitted; refer to PDF]

Likewise, a general integration solution (GIS) as given in FC and MFC method can be written as [figure omitted; refer to PDF] J can be obtained from (49).

The integration solution of the damage correction factor is [figure omitted; refer to PDF]

6.1. Ideal Bimodal Spectra

Different spectral shapes have been investigated in the literatures [5, 26]. Bimodal PSDs with two rectangular blocks are used to carry out numerical simulations, as shown in Figure 6. Two blocks are characterized by the amplitude levels _SωL and _SωH , as well as the frequency ranges _ωb -_ωa and _ωd -_ωc .

Figure 6: Ideal bimodal spectra.

[figure omitted; refer to PDF]

_{ω L} and _ωH are the central frequencies, as defined in [figure omitted; refer to PDF]

_{A 1} and _A2 represent the areas of block 1 and block 2 and are equal to the zero-order moment, respectively. For convenience, the sum of the two spectral moments is normalized to unity; that is, _A1 +_A2 =1.

To ensure that these two spectra are approximately narrow band processes, _ωb /_ωa =_ωd /_ωc =1.1.

Herein, two new parameters are introduced [figure omitted; refer to PDF]

The parameter values in the numerical simulations are conducted as follows: B=0.1, 0.4, 1, 2, and 9; R=6 and 10 which ensure two spectra are well-separated; m=3, 4, and 5 which are widely used in welded steel structures. _ωa does not affect the simulated bandwidth correction factor ρ [10, 27]. Thus, _ωa can be taken as the arbitrary value in the theory as long as _ωa >0; herein, _ωa =5 rad/s.

In the process of time domain simulation, the time series generated by IDFT [28] contains 20 million sample points and 150-450 thousand rainflow cycles, which is a sufficiently long time history, so that the sampling errors can be neglected.

Figure 7 is the result of JM method. The bandwidth correction factor calculated by (23) is in good agreement with the results obtained from (52). For m=3, JM method can provide a reasonable damage prediction. However, for m=4 and m=5, this method tends to underestimate the fatigue damage.

Figure 7: The bandwidth correction factor of JM method for (a) R=6 and (b) R=10.

(a) [figure omitted; refer to PDF]

(b) [figure omitted; refer to PDF]

Figures 8 and 9 display the results of FC method and MFC method, which are both calculated by the proposed analytical solution (see (52)) and FC's numerical integration solution (see (54)). The bandwidth correction factor calculated by (52) is very close to the results obtained from (54). Besides, the FC method always provides conservative results compared with RFC. The MFC method improves the accuracy of the original FC method to some extent. However, it may underestimate the fatigue damage in some cases.

Figure 8: The bandwidth correction factor of FC method for (a) R=6 and (b) R=10.

(a) [figure omitted; refer to PDF]

(b) [figure omitted; refer to PDF]

Figure 9: The bandwidth correction factor of MFC method for (a) R=6 and (b) R=10.

(a) [figure omitted; refer to PDF]

(b) [figure omitted; refer to PDF]

As has been discussed above, the bandwidth correction factor calculated by the analytical solution has divergence with that computed from RFC. The disagreement is not because of the error from the analytical solution but arises from the fact that the physical models of the original bimodal methods and the rainflow counting method are different.

6.2. Real Bimodal Spectra

In practice, the rectangular bimodal spectra in previous simulations cannot represent real spectra encountered in the structures. Therefore, more realistic bimodal stress spectra will be chosen to predict the vibration-fatigue-life by using the analytical solution proposed in this paper. For offshore structures subjected to random wave loading, the PSDs of fatigue stress of joints always exhibit two predominant peaks of frequency. Wirsching gave a general expression to characterize the PSD, and the model has been widely used in a few surveys [10, 23, 29]. The analytical form is [figure omitted; refer to PDF] where A is a scale factor, _Hs is the significant wave height, _TD is the dominant wave period, _ωn is the natural angular frequency of structure, and ζ is the damping coefficient.

A and _Hs do not affect the value of the bandwidth correction factor ρ [10, 27]. Thus, they are chosen to be equal to unity for simplicity. The value of other parameters can be seen in Table 4.

Table 4: Parameters for real bimodal spectra.

	A	H s (m)	T D (s)	ω N (rad/s)	ζ
Group 1	1	1	15	2.2	1.5%
Group 2	1	1	5	10	1%

Real stress spectra corresponding to group 1 and group 2 can be seen in Figure 10. S(ω) is a double peak spectrum, the first peak is produced by the peak of random wave spectrum, and the second one is excitated by the first mode response of structures.

Figure 10: Real bimodal spectra corresponding to (a) group 1 and (b) group 2.

(a) [figure omitted; refer to PDF]

(b) [figure omitted; refer to PDF]

The results of different bimodal methods are shown in Figures 11(a) and 11(b). The bandwidth correction factor for JM method is calculated with (23) and (52). It should be noted that _λ0,L +_λ0,H ≠1 for real bimodal spectra in Figures 10(a) and 10(b), so, _λ0,L and _λ0,H should be normalized as (14), (15), and (16), while ρ for FC method and MFC method is obtained through (52) and (54).

Figure 11: The bandwidth correction factor of different bimodal methods for (a) group 1 and (b) group 2.

(a) [figure omitted; refer to PDF]

(b) [figure omitted; refer to PDF]

Under the real stress spectra, the zero-order spectral moments corresponding to low frequency and high frequency are very small. JM method provides acceptable damage estimation compared with RFC. The results computed by (54) for FC and MFC method may be not correct in some cases. The incorrect results are mainly caused by (49). However, the proposed analytical solution (see (52)) can give a satisfactory damage prediction and always provide a conservative prediction. More importantly, (52) is more convenient to apply in predicting vibration-fatigue-life than numerical integration.

7. Conclusion

In this paper, bimodal spectral methods to predict vibration-fatigue-life are investigated. The conclusions are as follows:

(1) An analytical solution of the convolution integral of two Rayleigh distributions is developed. Besides, this solution is a general form which is different from that proposed by Jiao and Moan. The latter is only a particular case.

(2) An analytical solution based on bimodal spectral methods is derived. It is validated that the analytical solution shows a good agreement with numerical integration. More importantly, the analytical solution has a stronger attraction than numerical integration in engineering application.

(3) For JM method, the original approximate solution (see (23)) is reasonable in most cases; the analytical solution (see (52)) can give better prediction.

(4) In ideal bimodal processes, JM method and MFC method may overestimate or underestimate the fatigue damage, while in real bimodal processes, they give a conservative prediction. FC method always provides conservative results in any cases. Therefore, FC method can be recommended as a safe design technique in real engineering structures.

Acknowledgments

This work was supported by the National Natural Science Foundation of China (Grant nos. 51409068 and 51209046). This work was also supported by China Postdoctoral Science Foundation (Grant no. 2015M571396). A part of this work was supported by Guangdong Natural Science Foundation (Grant no. 2014A030310194).

References

[1] V. V. Bolotin, Random Vibrations of Elastic systems , Martinus Nijhoff, Hague, The Netherlands, 1984.

[2] T. T. Soong, M. Grigoriu, Random Vibration of Mechanical and Structural Systems , Prentice Hall International, London, UK, 1973.

[3] F. Pelayo, A. Skafte, M. L. Aenlle, R. Brincker, "Modal analysis based stress estimation for structural elements subjected to operational dynamic loadings,", Experimental Mechanics , vol. 55, no. 9, pp. 1791-1802, 2015.

[4] N. W. M. Bishop, F. Sherrat, Finite Element Based Fatigue Calculations , NAFEMS Ltd, East Kilbride, UK, 2000.

[5] C. Braccesi, F. Cianetti, G. Lori, D. Pioli, "Fatigue behaviour analysis of mechanical components subject to random bimodal stress process: frequency domain approach,", International Journal of Fatigue , vol. 27, no. 4, pp. 335-345, 2005.

[6] D. Benasciutti, F. Sherratt, A. Cristofori, "Recent developments in frequency domain multi-axial fatigue analysis,", International Journal of Fatigue , vol. 91, pp. 397-413, 2016.

[7] M. Mrsnik, J. Slavic, M. Boltezar, "Frequency-domain methods for a vibration-fatigue-life estimation--application to real data,", International Journal of Fatigue , vol. 47, pp. 8-17, 2013.

[8] G. Jiao, T. Moan, "Probabilistic analysis of fatigue due to Gaussian load processes,", Probabilistic Engineering Mechanics , vol. 5, no. 2, pp. 76-83, 1990.

[9] T.-T. Fu, D. Cebon, "Predicting fatigue lives for bi-modal stress spectral densities,", International Journal of Fatigue , vol. 22, no. 1, pp. 11-21, 2000.

[10] D. Benasciutti, R. Tovo, "On fatigue damage assessment in bimodal random processes,", International Journal of Fatigue , vol. 29, no. 2, pp. 232-244, 2007.

[11] L. D. Lutes, M. Corazao, S.-L. J. Hu, J. Zimmerman, "Stochastic fatigue damage accumulation,", Journal of Structural Engineering , vol. 110, no. 11, pp. 2585-2601, 1984.

[12] T. Dirlik, Application of computers in fatigue analysis [Ph.D. thesis] , University of Warwick, 1985.

[13] J.-B. Park, J. Choung, K.-S. Kim, "A new fatigue prediction model for marine structures subject to wide band stress process,", Ocean Engineering , vol. 76, pp. 144-151, 2014.

[14] J.-B. Park, C. Y. Song, "Fatigue damage model comparison with formulated tri-modal spectrum loadings under stationary Gaussian random processes,", Ocean Engineering , vol. 105, pp. 72-82, 2015.

[15] D. Benasciutti, C. Braccesi, F. Cianetti, A. Cristofori, R. Tovo, "Fatigue damage assessment in wide-band uniaxial random loadings by PSD decomposition: outcomes from recent research,", International Journal of Fatigue , vol. 91, pp. 248-250, 2016.

[16] C. Braccesi, F. Cianetti, L. Tomassini, "Random fatigue. A new frequency domain criterion for the damage evaluation of mechanical components,", International Journal of Fatigue , vol. 70, pp. 417-427, 2015.

[17] C. Braccesi, F. Cianetti, L. Tomassini, "Validation of a new method for frequency domain dynamic simulation and damage evaluation of mechanical components modeled with modal approach,", Procedia Engineering , vol. 101, pp. 493-500, 2015.

[18] C. Braccesi, F. Cianetti, L. Tomassini, "An innovative modal approach for frequency domain stress recovery and fatigue damage evaluation,", International Journal of Fatigue , vol. 91, pp. 382-396, 2016.

[19] M. A. Miner, "Cumulative damage in fatigue,", Journal Applied Mechanics , vol. 12, pp. A159-A164, 1945.

[20] American Society for Testing and Materials (ASTM), "1049-2005 Standard Practices for Cycle Counting in Fatigue Analysis,"

[21] J. S. Bendat, A. G. Piersol, Random Data. Analysis and Measurement Procedure , A Wiley Interscience Publication, London, UK, 2000.

[22] J. S. Bendat, "Probability functions for random response: prediction of peaks, fatigue damage and catastrophic failures,", NASA, 1964.

[23] P. H. Wirsching, C. L. Light, "Fatigue under wide band random stresses,", ASCE, Journal of the Structural Division , vol. 106, no. 7, pp. 1593-1607, 1980.

[24] API, "Recommended practice for design and analysis of station keeping systems for floating structures,", API RP 2SK , 2005.

[25] M. Abramowitz, I. A. Stegun, Handbook of Mathematical Functions , United States Department of Commerce, 1965.

[26] D. Benasciutti, R. Tovo, "Spectral methods for lifetime prediction under wide-band stationary random processes,", International Journal of Fatigue , vol. 27, no. 8, pp. 867-877, 2005.

[27] C. Han, Y. Ma, X. Qu, M. Yang, P. Qin, "A practical method for combination of fatigue damage subjected to low-frequency and high-frequency Gaussian random processes,", Applied Ocean Research , vol. 60, pp. 47-60, 2016.

[28] R. T. Hudspeth, L. E. Borgman, "Efficient fft simulation of digital time sequences,", ASCE J Eng Mech Div , vol. 105, no. 2, pp. 223-235, 1979.

[29] V. Bouyssy, S. M. Naboishikov, R. Rackwitz, "Comparison of analytical counting methods for Gaussian processes,", Structural Safety , vol. 12, no. 1, pp. 35-57, 1993.

Appendix

Appendix

A. Derivation of (44)

[figure omitted; refer to PDF] By introducing an integral transformation [figure omitted; refer to PDF] (A.1) can be changed as [figure omitted; refer to PDF] in which the standard Normal distribution function can be expressed as a form of an error function. [figure omitted; refer to PDF] where the error function can be written as [figure omitted; refer to PDF] Equation (A.3) then becomes [figure omitted; refer to PDF] The first item for (A.6) on the right side is [figure omitted; refer to PDF] According to the integral transformation shown in Figure 12, the second item for (A.6) on the right side is [figure omitted; refer to PDF] where the definition of [figure omitted; refer to PDF] is shown in Figure 12 and it can be calculated as [figure omitted; refer to PDF] . Therefore, the analytical solution of (A.1) becomes [figure omitted; refer to PDF]

Figure 12: Integral transformation.

[figure omitted; refer to PDF]

B. MATLAB Program for (38)

See Algorithm 1.

Algorithm 1

function Int_eq = Int_M_eq3 (L01, L02, m)

% MATLAB program for Eq. (38).

% L01 represent the zero-order spectral moment for spectra 1

% L02 represent the zero-order spectral moment for spectra 2

% m is fatigue strength exponent

error(nargchk(2,3,nargin))

if nargin < 3 || isempty(m)

m = 0;

end

if m < 0

Int_eq = [figure omitted; refer to PDF] ;

Return

else

if rem(m,1) > 0

Int_eq = [figure omitted; refer to PDF] ;

return

end

end

I0 = sqrt(2*pi)/4+(pi/2-atan(sqrt(L02/L01)))/sqrt(2*pi);

I1 = 0.5 + sqrt(L01/(L01+L02))/2;

if m == 0

Int_eq = I0;

return

end

if m == 1

Int_eq = I1;

return

end

DoubFac_m = DoubleFactorial(m-1);

C1 = DoubFac_m * I1;

C2 = DoubFac_m * I0;

C3 = sqrt(L01/L02)/2/sqrt(2*pi);

C4 = sqrt(2 * L02/(L01+L02));

if mod(m,2)==1

C5 = zeros((m+1)/2,1);

for k = 2 : 1 : (m+1)/2

DoubFac_k = DoubleFactorial(2*k-2);

C5(k) = C4^⋀ (2*k-1) * gamma((2*k-1)/2)/DoubFac_k;

end

Int_eq = C1 + C3 * DoubFac_m * sum(C5);

else

C5 = zeros(m/2,1);

for k = 1 : 1 : m/2

DoubFac_k = DoubleFactorial(2*k-1);

C5(k) = C4^⋀ (2*k) * gamma(k)/DoubFac_k;

end

Int_eq = C2 + C3 * DoubFac_m * sum(C5);

end

return

%-subroutine for caculating double factorial.

function b = DoubleFactorial(a)

% a is a integral number.

b = 1;

for i = a : -2 : 1

b = b*i;

end

return

%-end-