Statistics of extreme values of stochastic processes generated by rotating machines and associated probabilistic reliability model

The random stress processes recorded on the mechanical structures of rotating machines have non-Gaussian structure. These processes are logically composed of a deterministic periodic process, centred on the harmonics of the rotor rotation frequency, on which is routinely superposed a Gaussian zero-mean random process. This article first presents the mathematical principles of extreme value models, adapted to their specific nature. The proposed probabilistic models are then compared with each other in order to examine their degree of similarity and conservatism. On the basis of Gumbel’s theory, adopted as an asymptotic approach to extreme values, a statistical model of extreme values of sine plus noise composite random processes is proposed and discussed.


Introduction
The random stress processes recorded on the mechanical structures of rotating machines are non-Gaussian in structure [1], and consequently necessitate the use of probabilistic dimensioning approaches, adapted to their specific nature.The non-Gaussian nature of the vibration environments produced by a large number of carriers in the civil and military sectors on embedded equipment is also due to high-level rotating processes (engines in the case of motor vehicles [2][3][4] and/or tracks in the case of tracked vehicles [5][6][7] and/or rotors and propellers in the case of helicopters [8] or propeller-driven aircraft [9]).
These technical considerations mean that structures at extreme stresses can no longer be dimensioned on the basis of conventional probabilistic theories, developed historically in the context of stationary, ergodic and Gaussian processes, and widely incorporated into FE (Finite Element) design tools [10].On the basis of previous work on the dimensioning of mechanical structures subjected to the vibratory excitations produced by tracked vehicles [11], an analytical and numerical approach to the local maximum values induced by this type of sine plus noise composite process is proposed.This local statistical approach is then used to examine the statistics of the overall maxima of such processes, which by nature have a Corresponding author: b.colin@nexter-group.fr a great impact on the dimensioning of mechanical structures when the fundamental modes of such structures are locked on one of the rotational harmonics of the rotating machine under study.
First, the mathematical principles of the extreme value models, adapted to the non-Gaussian probabilistic nature of this type of sine plus noise composite process, are set out.The proposed probabilistic models are then compared with each other in order to examine their degree of similarity and conservatism.On the basis of Gumbel's theory, adopted as an asymptotic approach to extreme values, a statistical model of extreme values of sine plus noise composite random processes is proposed and discussed.

Nomenclature
In addition, for analytical simplification reasons, the narrow-band random process b (t, ξ) is supposed to be Gaussian, with mean value equal to zero and variance equal to σ 2 , whose law of density of probability is given by the following expression, namely: For successful dimensioning of structures subjected to this type of narrow-band sine plus noise composite random stress process x (t, ξ), the probability density function of the local maxima of such a stress process should first be established.To accomplish this, we show that it is then possible to use the analytical expression of the envelope Z x (t, ξ) process and to use the analytical expression of the phase φ x (t, ξ) process of this narrow-band process x (t, ξ), centred on f 0 ≈ f 1 Therefore, we obtain: (3) From these technical considerations it becomes clear that the analytical model of the envelope of the narrow-band random process considered is accessible using the theory of analytical signals and thus of complex envelopes [14], through the Hilbert transform.This principle was applied by the author in 1992 to establish the expression for the probability density of the local maxima of narrow-band sine plus noise composite random processes [15].

Rice's analytical model
The analytical model of this local maxima probability density is defined as [15]: where σ 2 is the variance of the narrow-band random stress process b (t, ξ), whose the mean frequency is f 0 and the mean value is zero.and S is the amplitude of the sinusoidal random process s (t, ξ), whose the frequency is f 1 and the mean value is zero.
The dimensionless coefficient a 2 0 represents the severity of the sine plus noise composite process x (t, ξ) = s (t, ξ)+b (t, ξ), as it is the square of the signal-to-noise ratio rms values (named also SNR:Signal-to-Noise Random, in anglo-saxon technical literature).
This local maxima probability density P Zx (z), shown in Figure 1 below, is thus the product of two simple analytical functions: And where I 0 (a) is the modified Bessel function of the first kind, of dimensionless argument a = zS/σ 2 , which is a symmetrical function as shown in Figure 2 below.The probability density P Ray (z; σ)is the Rayleigh distribution relating to the case of narrow-band pure random processes and h(a 2 0 , a) is the multiplicative function of the Rayleigh distribution relating to the probability density of the extrema of the sine plus noise narrow-band composite random processes.This probability density P Zx (z), defined by expression (4), is also referred to in the literature as the Rice function, after its originator.This probability function was defined during his work on radio communication [16], and is a generalisation of the Rayleigh function used to describe the behaviour of a radio signal propagating along several paths before being received by a radar antenna.

Dimensionless expression of the Rice function
In order to better define the behaviour of the Rice function, it is useful to present the function in its dimensionless form, by normalising the random amplitude  Z x (ξ) of the envelope Z x (t, ξ) with respect to the rms value σ of the gaussian noise b (t, ξ) alone.
The dimensionless random variable is then defined by expression (7) below: Rice's dimensionless expression can then be written in the following form: for η 0 (8) At this stage it is important to note that the entire Rice distribution shifts to the right (towards high normalised maxima values) as the dimensionless coefficient a 2 0 of the process increases (case of a sinusoidal random process of increasing rms value for a Gaussian random process of constant variance).
It can be seen that the zero-order modified Bessel function of the first kind is an even function, of value one for a zero argument value.The latter increases continuously on a positive support, for which the analytical expression is: exp(x cos θ)dθ for x ∈ (9)

Asymptotic expression of the Rice function
In order to successfully study the dimensioning of structures at the extreme stresses induced by these narrow-band sine plus noise composite random processes, the analytical expressions of the asymptotic functions of the local maxima must be defined.It is effectively very useful to have analytical knowledge of the mathematical form of the distribution tail of the random stress process local maxima, in order to define a posteriori the form of the probability density function (p.d.f) of the highest maximum generated by the sine plus noise composite process considered over its entire load time T .Consequently, in order to obtain the normalised asymptotic expression of the Rice function, it is possible to use the asymptotic expression [17] of the zero-order modified Bessel function of the first kind, whose the analytical expression is defined by: Inserting expression (10) in expression (8) of the normalised probability density of the Rice function, the normalised asymptotic expression of the Rice function, P asymp η (η), can then be established: It can be seen that this dimensionless asymptotic function is an exponential decay function and depends only on the dimensionless coefficient a 2 0 , defined in Equation (5) above.

Local maxima distribution function
Following the determination of the probability density function (p.d.f) of the local maxima of the composite stress random process, it is now essential to be able to obtain the expression for the local maxima distribution function, in order to access the probability of exceeding a given threshold z 0 = βσ (see Fig. 3) by the envelope Z x (t, ξ) of the composite process considered x (t, ξ).Based on Figure 4 below and on the definition of the cumulative function repartition P Zx (z 0 ) of the random variable Z x (ξ) defined by expression (12) below, it can be shown that this cumulative function repartition (c.f.r) is written using Q 1 (α, β) defined by expression (13) below.
With for analytic expression of the first-order Marcum function: At this stage, it is important to note that, for a given risk of exceeding the threshold, the normalised threshold level β is higher the greater the value of the dimensionless coefficient a 2 0 of the sine plus noise random process.

Series expressions of the first-order Marcum function
The first-order Marcum function Q 1 (α, β) is by definition between 0 and 1 and decreases continuously, as shown in Figure 4 above.It is expressed in the form of an unbounded integral (13), which can be put in the form of convergent infinite series [15,18] as below: The general term of the infinite series ( 14) of the firstorder Marcum function includes two characteristic functions well known to mathematicians: factorial function k! which is expressed on the basis of the gamma function Γ (k + 1), as follows: complementary gamma function γ 1 (a, x) which is expressed on the basis of the gamma function Γ (a) and the incomplete gamma function γ (a, x), as follows: Also, it is possible to define a convergent infinite serie, as follows: after [18] (17 The general term of the infinite series (17) of the firstorder Marcum function includes the nth-order modified Bessel function of the first kind, defined as follows: To manage the calculation error on the numerical values of the first-order Marcum function [20], these convergent series formulations ( 14) and ( 17) are used in programming, in the scientific tools developed by the author for fatigue design of structures mounted on rotating machines (for example engine equipment and equipment on tracked vehicles [11]).

Asymptotic expression of the first-order Marcum function
With the aim of being able to use a fast calculation method employing elementary mathematical functions to evaluate the extreme stresses generated by a sine plus noise composite process, a simple analytical form of the asymptotic expression of Q 1 (α, β) is now needed.Using the asymptotic expression of the zero-order modified Bessel function of the first kind (10), inserted into the integral expression (13), the asymptotic expression Q asymp 1 (α, β) of the first-order Marcum function can be defined as follow: Then, on the basis of the theoretical work of Nuttall [21], it is possible to write: where Φ(x) is the distribution function of the standard normal distribution.Consequently, it can be seen that the first asymptotic expression Q asymp 1 (α, β) of the firstorder Marcum function is proportional to the distribution function of the standard normal distribution, which is an elementary mathematical function, expressed by: ) Knowing that error function erf(x) is an integral function related to the complementary error function erfc(x) expressed as follows: And, From these technical considerations, the second asymptotic function of the first-order Marcum function can also be written using the complementary error function (Nuttall-Colin):

Design of structures such that the local stress maximum occurs on average once over the load time T of the process
By definition, the local maximum stress ZL max of the composite random process x (t, ξ) likely to occur on average once over the given load time T of the random process x (t, ξ), whose the mean frequency is equal to f 0 ≈ f 1 , is given by the expression: ZL max = σ β a 2 0 , n with α = 2a 2 0 and n = f 1 T fixed (25) The value β a 2 0 , n , associated with the local maximum stress of the composite process for any number n of load cycles, then corresponds to the root of the characteristic function g (α, β, n) defined by: For the pure gaussian random process (a 2 0 = 0), it is possible to write the following expression, in using the expressions ( 13) and (26): The following trivial relationship is also obtained: The expressions can thus be used to plot the graph of the normalised local maxima of a sine plus noise composite random process, as shown in Figure 5.
The value of β thus represents the normalised local maximum stress, which corresponds in fact to the threshold of the local maxima of the composite stress process considered, and the risk of exceeding this threshold is equal to the reciprocal of the mean number of load cycles n = f 1 T .The function β a 2 0 , n for different dimensionless values of a 2 0 as a function of n ∈ [1.3E + 04 cycles] is shown in Figure 5.
Having previously established the asymptotic function of the first-order Marcum function in Section 3.2, It shows that the values obtained are very similar and that the errors generated are less than 1%, giving credibility to the asymptotic calculation of βasymp , which requires a much shorter calculation time than that needed to estimate β.This asymptotic method will be used for the tailoring tests in composite environments, for example, as in such cases the number n of load cycles is always high and βasymp is calculated for each degree of freedom considered (ranging from 1000 to 2000 d.o.f) by vibration environment situation to be characterised in terms of extreme stresses.
Companies involved in the dimensioning of mechanical structures of rotating machines generally apply conventional processing techniques, using the power spectral density (PSD) concept and assuming normal distributions of the stochastic processes studied, even though they are not appropriate given the statistical nature of the composite processes.To identify the design errors that can be made in the mechanical strength of structures when using a stochastic approach that is not appropriate for the statistical nature of the physical stress processes occurring in rotating machines, in what follows a comparison is proposed between the stochastic approach developed previously ( ZL max ) and the approach using Gaussian processes ( ẐL max ), applying constant rms values in these two probabilistic approaches.When all calculations are completed, Figure 6 and the expressions (29) and (30) below are obtained, with regard to this model of comparison (of design at extreme values) between the Sine plus Random Composite approach (C S+R) and the iso-energy gaussian random approach, (GR IE) in expression (29): And logically, inserting (27) in (29), the following is obtained: 1 (for gaussian random process) for any value of n (30) At this stage, it is important to note that a Gaussian random approach such as PSD is always more conservative than the sine plus random composite approach in terms of local maximum stress, however severe the composite random stress process (narrow-band process) applied.This justifies the use of probabilistic structure dimensioning models, appropriate for the statistics of the excitation processes, with a view to minimum necessary dimensioning, enabling the design to be optimised.

Design of structures in terms of the overall stress maximum with defined risk δ of exceeding the threshold (controlled-risk dimensioning)
By definition, ZG max (δ) at imposed risk δ, associated with the composite process x (t, ξ), corresponds to the At this stage, it can be seen that it is now possible to determine the design value ZG max (δ) to which the designer has to design the structure, for a risk δ that the designer defines according to the level of reliability required by the customer, as follow [23]: For the case of 'standard' structures and for the case of 'secured' structures, the dimensionless coefficients [− ln [− ln (1 − δ)]] appearing in expression (33) of the overall maximum at imposed risk δ are given in Table 2.
The comparative approach between the secured design approach with regard to the local maximum ZL max and  the overall maximum approach ZG max (δ) with a defined risk δ is then described by the basic expression (34) below (in terms of Gumbel's asymptotic approximation), and by Figures 7 and 8, associated with a risk δ of 1% and 0.1%, respectively, of exceeding the threshold.At this stage, it is possible to compare composite Sine plus Random model at fixed risk δ (δ-C S+R) with the Composite Sine plus Random Model (C S+R), which corresponds to overall maxima approach divided by local maxima approach for a Composite Sine plus Random process.
Note: This ratio depends on a 2 0 , n and δ (fixed at 1% or 0.1%) Figures 7 and 8 above show clearly that the overall maxima approach with defined risk δ of exceeding the threshold leads logically to securing of the design of the structures subjected to random stress processes compared with a dimensioning approach using the local maximum stresses, used routinely by designers in dynamic dimensioning.
This work shows that the overall maxima approach with defined risk δ of exceeding the threshold leads logically to securing of the design of structures subjected to random stress processes compared with a dimensioning approach using the local maximum stresses, used routinely by designers in dynamic dimensioning.
In the case of sine plus random composite processes, it can be seen that the higher the value of the dimensionless coefficient a 2 0 , the lower the level of conservatism of the overall maxima approach for a given risk δ and number of cycles n.
Lastly, it is also shown that, for a stress random process generated by a rotating machine, the lower the stipulated risk δ, the higher the severity of the overall approach.And this conclusion constitutes an logical and very relevant result in a security context of design.

Empirical design model of structures subjected to sine plus noise composite random process
In practice, industrial firms concerned by dynamic dimensioning of structures subjected to sine plus noise composite processes are led to extract the sinusoidal nature from the composite process by adapting Kalman filtering techniques [24].Extraction of the periodic character of the process in the time domain by Kalman filtering is a way of characterising the random nature of the composite process correctly, by subtracting the extracted sinusoidal process from the composite process considered.On the basis of this extraction, designers are led to estimate the extreme stresses generated by each of these two individual processes treated separately, using deterministic time-based methods for the sinusoidal process and spectral methods for the random process [25].
Given the lack of availability of computer codes capable of handling the complex nature of such sine plus noise composite processes, the conservative dimensioning attitude then consists in summing the extreme stresses of these two individual processes, characterised as described above, in order to dimension their structures for the extreme stresses.With a view to describing the conservative character of this empirical dimensioning method, based on the summing of the individual effects, it appeared of interest to demonstrate the usefulness of having dimensioning models taking account of the non-Gaussian nature of the composite processes.
Consequently, in what follows the relationship is established between the local maximum stress ( ZL max ) of the sine plus noise composite process (C S+R) and the summing of the individual stresses (S+R) associated with the sinusoidal (S) and narrow-band Gaussian random (σ 2 ln(f 1 )) processes according to (27).The empirical model of dimensioning to the individual extreme stresses (S+R) normalised with respect to the sine plus noise composite model (C S+R) leads to the following analytical The results given in Figure 9 above show that use of the empirical model of summing the individual effects leads to over-estimation of the local extreme stress of the sine plus composite process, however severe the composite process in terms of load time T and dimensionless ratio a 2 0 .Consequently the empirical approach of summing the effects is observed to be a conservative approach, with a limited severity increase of 20%, which enables its continuing use in the design of tracked running gear components, for example, at Nexter Company (Land defence systems French integrator).

Conclusion
Modelling of random load processes in terms of imposed displacement and/or imposed force is necessary for optimum dimensioning of all the equipment of industrial products subjected to high load levels.In terms of dimensioning to extreme stresses, many of the models used are based on the theory of Gaussian random processes, applying conventional spectral approaches now incorporated into all finite element program codes (FE codes).
However, in the case of complex processes, such as the ones described in this paper, encountered in the rotating machinery sector, the sine plus noise composite nature requires more sophisticated calculation approaches that are not widespread in the specialist literature.In fact, although the probabilistic nature of local extreme stresses can be supported positively by the probabilistic model of S.O Rice, the models developed by the author show the need to enrich it through an overall asymptotic approach, applying the extreme value theory (EVT).These approaches of dimensioning to overall maximum stresses are all the more important in that they provide designers with an effective dimensioning model, able to manage equipment dimensioning criticality by specifying a risk of exceeding limits, as explained in Section 5.

Fig. 6 .
Fig. 6.Composite model normalised with respect to the Gaussian noise model, at iso-energy.

Table 2 .
Coefficient [− ln [− ln (1 − δ)]] for δ = 1% and δ = 0.1%.Risk Delta − ln(− ln(1 − Delta)) of the composite process generated over the load time T .This stress is associated with a risk δ defined by the designer, according to the level of criticality of the product.In general, the designers define a risk of exceeding the threshold δ of 1% in the case of 'standard' structures and 0.1% in the case of 'secured' structures.As the probability density function of the composite random process local maxima is described asymptotically by an exponential decrease, it can be shown that the distribution function of the highest maxima Y n (ξ) = MAX k=1 to n [z (t k , ξ)] over the duration T of the composite random process tends to a Gumbel function [22]: F Yn (y) − −−− → high n exp [− exp (−α n (y − y n ))] where y n = σ βasymp with α n = nP asymp Zx α n is the value of the maximum nth-order extreme intensity of the local maximum Z x (ξ), and where y n is the maximum nth-order characteristic value of Z x (ξ).The following expressions are obtained for α n and y n : y n = ZL max = σ β where β

Fig. 9 .
Fig. 9. Composite model (C S + R) normalised with respect to the empirical model of summing of the individual maximum stresses (S + R).

Table 1 .
(24)arison between β a 2 0 , n and βasymp a 2 0 , n obtained by the Nuttall-Colin asymptotic approach.canalsobeusedto obtain the numerical values of the normalised local maximum stress βasymp , when the number n of cycles is high.Using the Nuttall-Colin asymptotic expression(24), the comparison in Table1can be made between the values β and βasymp for values of a 2 0 ∈ [1, 10] and n ∈ [1E + 02, 1E + 07].