Djaffar Boussaa
email: boussaa@lma.cnrs-mrs.fr

Robert Bouc
email: bouc@lma.cnrs-mrs.fr

Elastic Perfectly Plastic Oscillator Under Random Loads: Linearization and Response Power Spectral Density

Keywords: nonlinear random vibration, hysteretic systems, linearization with random parameters

The randomly-excited elastic-perfectly-plastic oscillator-hysteretic bilinear oscillator with zero secondary stiffness-has been extensively researched. The vast majority of the work on that system has investigated the time-domaine statistics of the response. No studies have focused on the power spectral densities. This study specifically examines the system's velocity power spectral density-the system's displacement is nonstationary-under wide-band random excitations by means of a statistical linearization/stochastic averaging technique that is developed in one existing and two new procedures, one with constant, the other with amplitude-dependent parameters. The three procedures are evaluated against Monte Carlo simulation and classical Gaussian linearization in terms of the velocity average power, the peak frequency of the power spectral density, its peak value, its bandwidth, and its overall shape around the main frequency. The best predictions were yielded by the new procedure with amplitude-dependent parameters, which combines an extended amplitudephase transformation, a linearization with random parameters into a Maxwell system, a correlationfunction-based criterion, the conditional power spectral density concept, and a power conservation correction step. The new procedures can be adapted to apply to other hysteretic systems.

Introduction

This study deals with the response of the symmetric elastic perfectly plastic (EPP) oscillator subjected to wide-band random excitations. More precisely, the system under consideration is the hysteretic bilinear oscillator with the same elastic stiffness and the same yield limit in traction and in compression and zero secondary stiffness.

The vast majority of the work on this extremely nonlinear system has focused on the timedomain statistics of various response quantities . In comparison, the response power spectral densities (PSDs) have drawn little attention despite their practical importance, as in structural fatigue life assessment [START_REF] Boussaa | Fatigue-seismic ratcheting interactions in pressurized elbows[END_REF][START_REF] Nies | Spectral method in multiaxial random fatigue[END_REF] and floor spectrum modeling [START_REF] Politopoulos | Some aspects of floor spectra of 1dof nonlinear primary structures[END_REF]. In [START_REF] Bouc | Drifting response of hysteretic oscillators to stochastic excitation[END_REF], where the focus is on the nonstationarity of the displacement response, the authors have obtained an approximation of the value of the velocity PSD at zero frequency, which is the asymptotic rate of the variance of the displacement response. The rest of the spectrum has remained unstudied. The aim of this study is to obtain an analytical estimate to the velocity PSD, which can be used to estimate its peak frequency, its peak value, its bandwidth, its overall shape around the main frequency, and the velocity average power.

The methodological approach taken combines statistical linearization with stochastic averaging [START_REF] Roberts | Random vibration and statistical linearization[END_REF][START_REF] Socha | Linearization Methods for Stochastic Dynamic Systems[END_REF]. As far as the estimation of response PSDs are concerned, linearization procedures fall into two classes. In the first class, the nonlinear system is replaced by a linear system, the parameters of which are constant and determined by minimizing an error function in some sense. The PSD of the nonlinear system is then approximated by that of the equivalent linear system. The consensus is that this class of procedures gives correct estimates of the second-order response moments but may fail to predict the response PSD of strongly nonlinear systems [START_REF] Bouc | The power spectral density of response for a strongly non-linear random oscillator[END_REF][START_REF] Soize | Stochastic linearization method with random parameters for sdof nonlinear dynamical systems: prediction and identification procedures[END_REF]. In the second class, as initiated in the pioneering works of Crandall [START_REF] Crandall | Random forcing of nonlinear systems[END_REF] and Miles [START_REF] Miles | An approximate solution for the spectral response of duffing's oscillator with random input[END_REF], the nonlinear system is replaced by a family of linear systems with random parameters related to the system's response. The response PSD of the nonlinear system is then approximated by the average, in some sense, of the response PSDs of the family of equivalent linear systems. This latter class is regarded as an improvement on the first class as it has successfully predicted the response PSDs of various systems where linearization procedures of the first class failed [START_REF] Politopoulos | Some aspects of floor spectra of 1dof nonlinear primary structures[END_REF][START_REF] Bouc | The power spectral density of response for a strongly non-linear random oscillator[END_REF][START_REF] Soize | Stochastic linearization method with random parameters for sdof nonlinear dynamical systems: prediction and identification procedures[END_REF][START_REF] Fogli | Spectral response of a stochastic oscillator under impacts[END_REF][START_REF] Bellizzi | Analysis of multi-degree of freedom strongly non-linear mechanical systems with random input: Part II: Equivalent linear system with random matrices and power spectral density matrix[END_REF][START_REF] Rüdinger | Spectral density of oscillator with bilinear stiffness and white noise excitation[END_REF]. Two further points should be mentioned about the second class. First, to yield good results, it has often required a post-processing correction step that ensures that some quantities, reliable estimates of which exist, are conserved during approximations. Second, the research to date has tended to focus on nonlinear oscillators which, unlike the EPP system, have (i) a stiffness that increases with displacement (i.e., hardening oscillators) and (ii) a stiffness and damping that are separable. To the authors' knowledge, only one linearization procedure of either class has been proposed to linearize the EPP system. The first one is a (constant-parameter) Gaussian linearization [START_REF] Roberts | Random vibration and statistical linearization[END_REF][START_REF] Asano | An alternative approach to the random response of bilinear hysteretic systems[END_REF][START_REF] Spanos | Third-order statistical linearization-based approach to derive equivalent linear properties of bilinear hysteretic systems for seismic response spectrum analysis[END_REF]]. The second one, which is based on a mean-square criterion, is the random-parameter linearization procedure described in [START_REF] Bouc | Drifting response of hysteretic oscillators to stochastic excitation[END_REF]. Their ability to predict the velocity PSD, however, has not been assessed so far.

This study develops two new linearization procedures, one being a constant-parameter non-Gaussian linearization with a mean-square criterion and the other being an amplitude-dependentparameter linearization with a correlation-function-based criterion. The equivalent linear system is sought for as a Maxwell system. This choice has advantages. It preserves the dimension of the ordinary differential equation describing the displacement: The displacement response is described by a differential equation of order 3 for both the EPP and the Maxwell systems. It also preserves an important feature of the displacement response: The displacement response of both the EPP and the Maxwell systems exhibit drifting.

The remaining part of the paper is organized as follows. Section 2 provides the problem statement and sets up the machinery needed to develop the PSD approximation approach. Section 3 illustrates the use of the machinery by addressing the linear Maxwell system. Section 4 discusses the linearization procedures based on a mean-square criterion: the new constant-parameter linearization and (a variant of) the linearization described in [START_REF] Bouc | Drifting response of hysteretic oscillators to stochastic excitation[END_REF]. Section 5 develops the new amplitude-parameter linearization built using a criterion based on correlation functions. Section 6 provides the velocity PSD for each linearization. Section 7 compares the predictions of the three procedures between them, with those of the classical Gaussian statistical linearization, and with those of Monte Carlo simulations in terms of several PSD features. In Section 8, an average power conservation correction is proposed as a post-processing step to improve on the predictions of the amplitude-dependent linearization procedures.

Background

This section gives the equations governing the motion of the EPP oscillator. These equations are then reformulated in terms of new variables. This leads to a new system that lends itself to various approximations, such as stochastic averaging. The number of the new variables exceeding that of the old variables, closing relations are introduced. Moreover, because this study concentrates on the region around the velocity PSD's main frequency, it suffices to consider a white noise input. The qualitative and quantitative conclusions reached in this study are expected to hold for more general wide-band, stationary, physical inputs.

Below, the autocorrelation function R η and the PSD S η for a wide-sense stationary random process η are taken to be related by the following Fourier transform pair:

S η (ω) = ∞ -∞ R η (τ ) exp(-iωτ ) dτ, R η (τ ) = 1 2π ∞ -∞ S η (ω) exp(iωτ ) dω. (1)

Problem statement

The non-dimensionalized stochastic differential equation governing the motion of the EPP system has the form [START_REF] Bouc | Drifting response of hysteretic oscillators to stochastic excitation[END_REF] ẍt

+ z t = √ σ ẇt , (2)
żt = ẋt (1 -H (|z t | -1) H (z t ẋt )) , (3)
where the overdot denotes differentiation with respect to time t; x is the displacement; z is the restoring force (Figure 1); √ σ is the excitation intensity; H is the Heaviside function (H(x) = 0 if x < 0, 1 otherwise); w is the Wiener process; and ẇ, formal time-derivative of w, is the unit Gaussian white noise (S ẇ = 1). The parameter is introduced as an order of smallness. No damping other than the damping coming from the friction element is assumed. However, the methodology developed below can be readily adapted to apply in the presence of a linear viscous damping element in parallel with the spring/friction-element pair. The aim of this study is to estimate the PSD of the velocity response ẋ.

Extended amplitude-phase transformation

The first change of variables is the following amplitude-phase transformation:

x t = t 0 u τ dτ + a t cos Φ t , (4) ẋt
= u t -a t Ω(a t ) sin Φ t , (5)
Φ t = t 0 Ω(a τ ) dτ + ϕ t . (6)
Eq. ( 4) is the decomposition of the displacement, x, into a drift term t 0 u τ dτ and an oscillatory term a t cos Φ t . The displacement is nonstationary, and so is assumed the drift term, whereas the oscillatory term is regarded as stationary [START_REF] Bouc | Drifting response of hysteretic oscillators to stochastic excitation[END_REF]. The decomposition of the displacement translates into decomposition (5) of the velocity, ẋ, into the drift rate u and the oscillatory velocity -aΩ sin Φ. The velocity is stationary and so are assumed the drift rate and the oscillatory velocity. For brevity, the amplitude of the oscillatory term, a, is referred to as the amplitude. The yet unknown function Ω(a) is an amplitude-dependent frequency to be made precise later, and ϕ is a phase. Differentiating (4) with respect to t and equating the result to (5) gives ȧ cos Φa φ sin Φ = 0.

(

) 7
The second change of variables is a decomposition of z, namely,

z t = z(a t , Φ t ) + v t , (8)
where v is a fluctuation, and z is the periodic response to the quasi-static cyclic straining a cos Φ with a and ϕ regarded as fixed. That is, z is the unique 2π-periodic function with respect to Φ that solves the following differential equation:

Ω(a) ∂ z ∂Φ = -aΩ(a) sin Φ (1 -H (|z| -1) H (-zaΩ(a) sin Φ)) , (9)
which, owing to the positivity of a and Ω, simplifies to

∂ z ∂Φ = -a sin Φ (1 -H (|z| -1) H (-z sin Φ)) . ( 10
)
One can check that the solution to this equation is as follows. If 0 ≤ a ≤ 1, then

z = a cos Φ. ( 11
) Otherwise, if a > 1, then z(a, Φ) =          a(cos Φ -1) + 1, if 0 < Φ ≤ Φ 0 ; -1, if Φ 0 < Φ ≤ π; a(cos Φ + 1) -1, if π < Φ ≤ π + Φ 0 ; 1, if π + Φ 0 < Φ ≤ 2π. (12)
Here, Φ 0 = arccos (1 -2/a) .

The first cosine and sine coefficients of z, respectively denoted by h c and h s , play an important role below. They are given by

h c (a) = a 2π (2Φ 0 (a) -sin (2Φ 0 (a))), if a > 1; a, if a ≤ 1; (13)
and

h s (a) =    4 π 1 a -1 , if a > 1; 0, if a ≤ 1. (14)

Closing relations

The motion of the system was initially described by two processes, x and z, and two equations, (2) and (3). After the changes of variables ( 4)-( 6) and the decomposition (8), the system becomes described by the five processes a, ϕ, Φ, u, and v. The differential equations satisfied by these processes are now derived.

Substituting the velocity ẋ in [START_REF] Grossmayer | Elastic-plastic oscillators under random excitation[END_REF] and decomposition (8) of z into (2) leads to

-ȧ (Ω + aΩ ) sin Φ -a φΩ cos Φ + r + ( u + v) + h c (a) -aΩ 2 cos Φ = √ σ ẇ, (15)
where the prime ( ) denotes differentiation with respect to a and r is the sum of the first sine harmonic and the higher harmonics of the Fourier series of z, i.e.,

r = z -h c cos Φ. ( 16
)
Inspection of (15) suggests the following:

Ω 2 (a) = h c (a) a , (17) u
+ v = 0. ( 18
)
Substituting ( 17) and ( 18) into (15) and using ( 6) and [START_REF] Spencer | Reliability of randomly excited hysteretic structures[END_REF] shows that a, ϕ, and Φ are governed by the following (Stratonovich) stochastic differential system:

ȧ = r ΩD sin Φ - √ σ ΩD (sin Φ) ẇ, (19) φ
= r aΩD cos Φ - √ σ aΩD (cos Φ) ẇ, (20)
Φ = Ω + r aΩD cos Φ - √ σ aΩD (cos Φ) ẇ, (21)
where

D = 1 + (aΩ /Ω) (sin Φ) 2 .
The system ( 19)-( 21) is in a form suitable for application of the principle of stochastic averaging. Now substituting decomposition (8) of z into the constitutive equation ( 3) yields the differential equation in v:

v = (y + u) 1 -H |z + v| -1 H (z + v) (y + u) - ∂ z ∂φ (Ω + φ) - ∂ z ∂a ȧ. ( 22
)
where ẋ = y + u, and y is the oscillatory velocity, namely,

y = -aΩ sin Φ, (23)
so that Eq. ( 18), ( 19)- [START_REF] Bensoussan | An analytical approach for the growth rate of the variance of the deformation related to an elasto-plastic oscillator excited by a white noise[END_REF], and ( 22) form the stochastic differential system governing the motion of the EPP system in terms of the a, ϕ, Φ, u, v processes.

Averaged equations

Because of the whiteness of the input, applying the averaging principle here amounts to averaging with respect to Φ the drift vector and diffusion tensor of the Itô equation corresponding to the Stratonovich system ( 19)- [START_REF] Bensoussan | An analytical approach for the growth rate of the variance of the deformation related to an elasto-plastic oscillator excited by a white noise[END_REF]. The joint process (a t , ϕ t ) (respectively (a t , Φ t )) is replaced by a Markov diffusion joint process over R + × [0, 2π], which, by an abuse of notation, is also denoted (a t , ϕ t ) (respectively (a t , Φ t )). The application of the stochastic averaging principle to the case in which the fast variable Φ has an integral dependence on a, as in ( 6), has been discussed in [START_REF] Bouc | The power spectral density of response for a strongly non-linear random oscillator[END_REF][START_REF] Bellizzi | Analysis of multi-degree of freedom strongly non-linear mechanical systems with random input. Part I: Non-linear modes and stochastic averaging[END_REF] (See also [START_REF] Fogli | The dynamics of a stochastic oscillator with impacts[END_REF] for a more formal mathematical treatment of this point).

Under the first harmonic assumption (i.e., z ≈ h c cos Φ + h s sin Φ or r ≈ h s sin Φ), the applica-tion of stochastic averaging to the system ( 19)-( 21) yields

ȧ = h s √ γ + γ Ω + σ 2 6γ -2γ 2 -3aγ 16aΩ 2 γ 5/2 + √ σ √ 2Ωγ 3/4 ẇ1 , (24) φ
= √ σ √ 2aΩγ 1/4 ẇ2 , (25)
Φ = Ω + √ σ √ 2aΩγ 1/4 ẇ2 , (26)
where ẇ1 and ẇ2 are independent unit white noises, and γ = 1 + aΩ /Ω. It can be shown that γ ≥ 1/4. The (stationary) joint probability density functions (pdfs) of (a, ϕ) ∈ R + × S 1 and (a, Φ) ∈ R + × S 1 , where S 1 is the unit circle, are then deduced from the Fokker-Planck equation as

p aϕ (a, ϕ) = 1 2π p a (a) = p aΦ (a, Φ), (27)
with

p a (a) = CaΩ 3/2 γ 3/4 H(a) exp 4 σ 2 a 0 h s (α)Ω(α) γ 1 + √ γ dα , ( 28
)
where C is a normalization constant.

It can be shown (as in the appendix of [START_REF] Roberts | The yielding behaviour of a randomly excited elasto-plastic structure[END_REF]) that the expectation of h s (a) with respect to p a (a), which is related to the average energy loss by cycle due to hysteresis, is of the order σ 2 .

The Maxwell system

This section has been included for several reasons: it presents the Maxwell system, which is the linear system used in all the linearizations to follow; it illustrates how the decompositions and quantities introduced in the previous section can be used to analyze that system; and, more importantly, it delineates a decomposition that highlights the respective contribution of the drift rate, u, and oscillatory velocity, y, to the PSD of the velocity, ẋ = u+y. One advantage of opting for the Maxwell system over the more commonly used Kelvin model is that the displacement response of the former shares one essential feature-nonstationarity-with the EPP system whereas that of the latter does not.

Velocity PSD of the Maxwell system

The response of the Maxwell system under white noise is governed by the following stochastic differential equations:

d ẋ dt + z = √ σ ẇ, ( 29
)
dz dt + Az = B ẋ, (30)
where A and B are two given constants such that A > 0, B > 0, and B -A 2 > 0.

The pair ( ẋ, z) is (asymptotically) zero-mean and stationary. The PSD of ẋ is given by

S ẋ(ω) = σ 2 ω 2 + A 2 (ω 2 -B) 2 + A 2 ω 2 . ( 31
)
It follows that S ẋ(0) > 0 , and therefore that the displacement response of the system exhibits drift [START_REF] Crandall | Accumulated slip of a friction-controlled mass excited by earthquake motions[END_REF]. It also follows that

E ẋ2 = 1 2π ∞ -∞ S ẋ(ω) dω = σ 2 2 1 A 1 + A 2 B . (32)

Decomposition of the velocity PSD

It is now shown that the decompositions and quantities introduced in the previous section make it possible to recover [START_REF] Fogli | Spectral response of a stochastic oscillator under impacts[END_REF] and [START_REF] Bellizzi | Analysis of multi-degree of freedom strongly non-linear mechanical systems with random input: Part II: Equivalent linear system with random matrices and power spectral density matrix[END_REF] and to obtain a natural decomposition S ẋ = S y + S u .

For this, the first step is to obtain the dynamical equations for a, ϕ, Φ, y, z, u, and v in the case of the Maxwell system. The periodic response z satisfies the following analog of Eq. ( 9):

Ω ∂ z ∂Φ + Az = By. (33)
Now, combining the identities Φ = Ω + φ, ẋ = y + u, and dz/dt = (∂ z/∂Φ) Φ + (∂ z/∂a) ȧ with [START_REF] Bensoussan | An ultra weak finite element method as an alternative to a Monte Carlo method for an elasto-plastic problem with noise[END_REF], and ( 19)-( 21) with D = 1, and ( 29)- [START_REF] Miles | An approximate solution for the spectral response of duffing's oscillator with random input[END_REF] gives, after some algebra, the equations governing the oscillatory components y and z as

dy dt + z = √ σ ẇ, (34)
dz dt + Az = By -rA + A √ σ ẇ, (35)
and those governing the perturbations u and v as

du dt + v = 0, (36)
dv dt + Av = Bu + rA -A √ σ ẇ. ( 37
)
The solutions of both systems can be expressed as the superposition of the response to noise terms and the response to rA terms. Because rA appears with opposite signs in the two systems, its contribution to ẋ = u + y vanishes. The same can be said about its contribution to z = v + z. In the rest of the section, (y, z) and (u, v) refer to the response to the noise terms alone. The next step is computing the PSDs. Let S u and S y denote the PSDs of u and y, respectively, and S uy their cross-spectrum. From ( 34)- [START_REF] Fogli | The dynamics of a stochastic oscillator with impacts[END_REF] with the term rA omitted, it may be concluded that

S u (ω) = σ 2 A 2 (ω 2 -B) 2 + A 2 ω 2 , ( 38
)
S y (ω) = σ 2 ω 2 (ω 2 -B) 2 + A 2 ω 2 , (39)
(S uy (ω)) = 0, [START_REF] Spanos | On the determination of the power spectrum of randomly excited oscillators via stochastic averaging: An alternative perspective[END_REF] where denotes the real part. Hence

S ẋ(ω) = S u (ω) + S y (ω), (41)
which gives back Eq. (31).

Energy considerations

From ( 38) and [START_REF] Bouc | Modèle mathématique d'hystérésis [mathematical model for hysteresis[END_REF], it follows that the stationary variance of the drifting rate, E u 2 , and the variance of the oscillatory component, E y 2 , are given respectively by

E u 2 = A σ 2 2B
, and

E y 2 = σ 2 2A . ( 42
)
Adding Eu 2 to Ey 2 gives back [START_REF] Bellizzi | Analysis of multi-degree of freedom strongly non-linear mechanical systems with random input: Part II: Equivalent linear system with random matrices and power spectral density matrix[END_REF]. The value of E y 2 in (42), obtained by integrating S y (ω) over R, coincides with the value obtained by averaging y 2 = (-aΩ sin Φ)

2 over the joint pdf ( 27)- [START_REF] Soize | Stochastic linearization method with random parameters for sdof nonlinear dynamical systems: prediction and identification procedures[END_REF] given by stochastic averaging, namely,

σ 2 2A = 1 2π 2π 0 ∞ 0 (-aΩ sin Φ) 2 p a (a) da dΦ. ( 43
)
To show this, note that p a in (28) reduces in the case of the Maxwell system to the Rayleigh distribution

p a (a) = a q 2 exp - a 2 2q 2 with q 2 = σ 2 2AΩ 2 = σ 2 2A (B -A 2 ) . (44)
The decomposition of S ẋ into S y and S u used here is applied to the velocity PSD of the EPP system within the context of linearization in Section 8. The equality between the value of E y 2 estimated by stochastic averaging and that obtained by integrating S y is not expected to hold in the EPP case, since then the estimate to E y 2 based on S y depends on the adopted linearization whereas the estimate based on the joint pdf of (a, Φ) is not. Imposing the equality of the two values, however, seems desirable. It strengthens the consistency of the approach as it enforces a conservation of the average power of the oscillatory velocity y. It also offers a natural means to introduce and determine a correction factor, as is described in Section 8.

Parameters of the equivalent Maxwell system via mean-square criteria

This section proposes two pairs of parameters (A, B) for the equivalent Maxwell system ( 29)-( 30) obtained via mean-square criteria. The first pair consists of constant parameters; the second, of amplitude-dependent parameters.

Constant parameters

The constant pair is new and is based on the following criterion: min

A,B ∞ 0 1 2π 2π 0 G (y, z) 2 dΦ p a (a) da, (45)
where

G (y, z) = y 1 -H (|z| -1) H (zy) + Az -By (46)
with y and z being the oscillatory components given by ( 23) and ( 11)-( 12), respectively. The criterion amounts to minimizing the expectation of the squared gap G using the joint pdf of the pair (a, Φ), given by ( 27)- [START_REF] Soize | Stochastic linearization method with random parameters for sdof nonlinear dynamical systems: prediction and identification procedures[END_REF]. From [START_REF] Paparizos | Some observations on the random response of an elasto-plastic system[END_REF], it follows that the problem is equivalent to min

A,B ∞ 0 1 2π 2π 0 Ω ∂ z ∂Φ + Az -By 2 dΦ p a (a) da. ( 47
)
Solving the problem for A and B gives

A = - 1 D ∞ 0 a 2 Ω(a) 4 p a (a) da ∞ 0 aΩ(a)h s (a)p a (a) da , (48)
B = 2 D ∞ 0 a 2 Ω(a) 4 p a (a) da ∞ 0 z 2 p a (a) da , (49)
where

z 2 = (1/2π) 2π 0 z (a, Φ) 2 dΦ and D = 2 ∞ 0 a 2 Ω 2 p a (a) da ∞ 0 z 2 p a (a) da - ∞ 0 aΩ(a)h s (a)p a (a) da 2 . (50)

Amplitude-dependent parameters

The amplitude-dependent pair, perhaps the first such pair that comes to mind, is described in [START_REF] Bouc | Drifting response of hysteretic oscillators to stochastic excitation[END_REF] and is defined by the following minimization problem:

min A,B 1 2π 2π 0 Ω ∂ z ∂Φ + Az -By 2 dΦ. ( 51
)
It follows that

A(a) = - aΩ 3 h s 2 z 2 -h 2 s and B(a) = Ω 2 2 z 2 2 z 2 -h 2 s . (52)

Linearization parameters via a correlation-function-based criterion

This section obtains a new pair of amplitude-dependent linearization parameters from a correlationfunction-based criterion. Because it is not clear how to directly apply this criterion to the EPP system due to the nonsmoothness of its governing equations, a simpler, less nonsmooth hysteretic system is first introduced as a nonlinear approximating system to the EPP system in a sense that will be made precise below. The simpler system is then linearized according to the correlationfunction-based criterion.

The approximating system

The approximating system used here is the following:

d ẋ dt + z = √ σ ẇ, ( 53
) dz dt + A * | ẋ| z = B * ẋ, (54)
where A * and B * are the loop parameters. The constitutive model ( 54) is the simplest Bouc-Wen (BW) constitutive model [START_REF] Bouc | Modèle mathématique d'hystérésis [mathematical model for hysteresis[END_REF]. It is rate independent, like the EPP system, but is less nonsmooth.

In keeping with the spirit of this study, A * and B * are sought as functions of a. Their determination is carried out along the same lines as the determination of A and B in the amplitude-dependent parameter linearization in Subsection 4.2. A gap between the EPP and BW systems is defined as

G (y, z) = y 1 -H (|z| -1) H (zy) + A * |y| z -B * y, (55)
which can be rewritten as

G (y, z) = Ω ∂ z ∂Φ + A * |y| z -B * y. (56)
Now, a variant of criterion ( 51) is used that requires the vanishing of the first two Fourier coefficients of G, namely,

2π 0 G (y, z) cos Φ sin Φ dΦ = 0 0 . ( 57
)
Under the additional z ≈ h c cos Φ + h s sin Φ, solving (57) for A * and B * yields

A * (a) = - 3π 4aΩ h s aΩ and B * (a) = Ω 2 + 2 h s aΩ 2 . (58)
A mean-square criterion as in (51) and/or the full expression ( 11)-( 12) for z can be used in (57), but this leads to unwieldy expressions for the parameters. Moreover, let ∂ z ∂Φ + A * |a sin Φ| z = -B * a sin Φ (59) be the analogue for the BW system of [START_REF] Borsoi | Approche probabiliste de la ruine d'un oscillateur élasto-plastique. In 2ème Colloque National AFPS 'Génie parasismique et aspects vibratoires dans le génie civil[END_REF]. It can be shown that the first order Galerkin approximation of the periodic solution of (59) coincide with the first order approximation of the periodic solution of [START_REF] Paparizos | Some observations on the random response of an elasto-plastic system[END_REF]. It follows that the energy dissipated over a cycle by the EPP system and by the BW system is the same.

Linearization parameters of the approximating system

The linearization of the BW system is carried out using correlation functions. From (53)-(54), it follows that correlation functions c ẋ ẋ(t, s) and c z ẋ(t, s) satisfy the following differential equations:

∂c ẋ ẋ ∂t (t, s) = -c z ẋ(t, s), (60)
∂c z ẋ ∂t (t, s) = B * c ẋ ẋ(t, s) -A * E [| ẋ(t)| z(t) ẋ(s)] , (61)
where t > s ≥ 0. Because of the term E [| ẋ(t)| z(t) ẋ(s)], the following decomposition is introduced to make progress analytically:

E [| ẋ(t)| z(t) ẋ(s)] ≈ αE [z(t) ẋ(s)] + βE [ ẋ(t) ẋ(s)] , (62)
= αc z ẋ(t, s) + βc ẋ ẋ(t, s),

where α and β are parameters to be determined. Note that if the response ( ẋ, z) of the BW system were Gaussian, then the decomposition would hold true, and using properties of Gaussian variables would make it possible to obtain α and β. That response being non-Gaussian, further approximation is required to determine them. This will be addressed shortly, but first note that decomposition (63) is more pertinent and better justifiable for a term of the form E [| ẋ(t)| z(t) ẋ(s)] than it is for a term of the form E [ ẋ(t)H (|z(t)| -1) H (z(t) ẋ(t)) ẋ(s)], whence the introduction of the BW system. Parameters α and β are obtained following the approach described in the previous section. A gap is introduced,

G (y, z) = y |y| z -αyz -βy 2 , ( 64
)
and associated with the following criterion:

These relations, which clearly show that A and B depend on a, are to be compared with those given by (52), which under the first-harmonic approximation (namely, z

2 ≈ h 2 c + h 2 s /2) reduce to A(a) = - h s aΩ and B(a) = Ω 2 + h s aΩ 2 . ( 68
)
The approximation introduced in decomposition (63), along with criterion (64), turns out to be a linearization approach. The coefficients A and B thus obtained are used as the linearization parameters of the linearized EPP system.

Velocity PSD of the EPP system

Damping ratio

The damping ratio is defined as

ξ = A 2 √ B . (69)
It is a constant in the case of a linearization with constant parameters, and a function of a in the case of a linearization with amplitude-dependent parameters. In the latter case, the average damping ratio is defined as

ξ = ∞ 0 A(a) 2 B(a) p a (a) da, (70)
where p a is given by (28).

General form of the velocity PSD

The non-Gaussian constant-parameter linearization procedure readily gives the velocity PSD, S ẋ, as [START_REF] Fogli | Spectral response of a stochastic oscillator under impacts[END_REF], and the corresponding velocity average power, E ẋ2 , as Eq. ( 32), with A and B given by Eq. ( 48) and (49), respectively.

The amplitude-dependent-parameter linearizations lead to the same expression [START_REF] Fogli | Spectral response of a stochastic oscillator under impacts[END_REF], except that now parameters A and B depend on a, where a is considered as a random variable (assumption motivated by the fact that process a t is slow) independent of ẇ with pdf p a given in [START_REF] Soize | Stochastic linearization method with random parameters for sdof nonlinear dynamical systems: prediction and identification procedures[END_REF]. The quantity S ẋ thus obtained is regarded as the "conditional velocity PSD given a" [START_REF] Bouc | The power spectral density of response for a strongly non-linear random oscillator[END_REF][START_REF] Spanos | On the determination of the power spectrum of randomly excited oscillators via stochastic averaging: An alternative perspective[END_REF]. To reflect this, the corresponding velocity PSD is denoted as

S ẋ(ω | a) = σ 2 ω 2 + A(a) 2 (ω 2 -B(a)) 2 + A(a) 2 ω 2 . ( 71
)
The velocity PSD, S ẋ(ω), is then defined as the average S ẋ(ω | a) over the pdf of a:

S ẋ(ω) = ∞ 0 S ẋ(ω | a) p a (a) da. ( 72
)
An immediate difficulty arises with S ẋ(ω | a) in (71), which is singular at ω = 1 for all 0 ≤ a ≤ 1, because then A(a) = 0 and B(a) = 1.

One approach to removing the singularity consists in replacing A(a) in the denominator of (71) by a parameter Ā(a) defined from the average damping ratio in (70) as follows:

Ā(a) = 2 ξ B(a). ( 73
)
This amounts to modifying (71) to

S ẋ(ω | a) = σ 2 ω 2 + A(a) 2 (ω 2 -B(a)) 2 + Ā(a) 2 ω 2 . ( 74
) It follows that E ẋ2 | a = σ 2 2 1 Ā(a) 1 + A(a) 2 B(a) . (75)
Finally,

E ẋ2 = σ 2 2 ∞ 0 1 Ā(a) 1 + A(a) 2 B(a) p a (a) da. ( 76
)

Linearization procedure vs Monte Carlo simulation

Below, the procedure based on the classical Gaussian linearization with constant parameters is referred to as the Gaussian linearization procedure; the procedure based on the non-Gaussian linearization with constant parameters, as the non-Gaussian linearization procedure; the procedure based on the linearization with amplitude-dependent parameters and mean-square criterion, as the direct linearization procedure; the procedure based on the linearization with amplitude-dependent parameters and a correlation-function-based criterion, as the indirect procedure.

General trends of the main frequency in simulation velocity PSDs

Examination of the simulation PSDs for 0.1 ≤ √ σ ≤ 0.9 reveals the following general trends (Figure 2):

• The peak frequency shifts to the left with increasing input level.

• The peak of the velocity PSD decreases with increasing input level, although the rate of decrease slows markedly for √ σ ≥ 0.5.

• The velocity PSD widens around its peak frequency with increasing input level.

The various procedures are now evaluated on how well they reproduce these trends, qualitatively and quantitatively. Particular attention is paid to the following quantities: average power, peak frequency, peak value, and bandwidth.

Velocity average power

The indirect linearization and the non-Gaussian linearization procedures gave accurate estimates of the velocity average power E ẋ2 over the whole input level range investigated. Although the Gaussian linearization procedure yielded acceptable predictions, it failed to capture the trend of the simulation data. The direct linearization procedure significantly overestimated the velocity average power over the whole input level range investigated (Figure 3A).

Velocity PSD peak frequency

The shift to the left of the peak frequency with increasing input level was reproduced qualitatively by the four procedures. The quantitative agreement varied from procedure to procedure. The Gaussian linearization procedure gave good predictions for 0.1 ≤ √ σ ≤ 0.7. For √ σ > 0.7, that procedure failed to reproduce the trend of the data. The non-Gaussian linearization procedure yielded acceptable results for 0.1 ≤ √ σ ≤ 0.3. For √ σ > 0.3, it significantly overestimated the frequency shift. The direct linearization procedure underestimated the frequency shift for 0.1 < √ σ ≤ 0.7, and grossly overestimated it for √ σ > 0.7. The indirect linearization procedure reproduced the trend of the simulation data over the whole input level range investigated. It outperformed the Gaussian linearization procedure for √ σ ≥ 0.7 (Figure 3B).

Velocity PSD peak

The direct linearization procedure significantly overestimated the peak of the velocity PSD over the whole input level range investigated. The Gaussian linearization procedure also overestimated the peak over a good portion of the range and failed to reproduce the trend of the data. The non-Gaussian linearization procedure gave good predictions for 0.1 ≤ √ σ ≤ 0.5. For √ σ > 0.5, its predictions deteriorated with increasing input level. The indirect linearization procedure is the only procedure whose predictions were in accord with the simulations over the whole input-level range investigated (Figure 3C).

Velocity PSD bandwidth

The bandwidth is defined here as half the ratio of the velocity average power to the peak of the velocity PSD (namely, E ẋ2 /(2 max ω S ẋ(ω))) [41, p. 265].

Although the constant-parameter-linearization procedures gave good predictions of the bandwidth for the lower values of the input level range, they failed to reproduce the trend of the simulation data. The direct linearization procedure underestimated the bandwidth for 0.2 ≤ √ σ ≤ 0.5, and gave good results for 0.6 ≤ √ σ ≤ 0.8. The indirect linearization procedure was the only one to give acceptable predictions over the whole input level range (Figure 3D). Taken together, these results suggest that the indirect linearization procedure was the best performer.

Average power conservation correction

Experience with random vibration of nonlinear systems has shown that the predictions of various approximation methods are often enhanced by implementing a mechanism that ensures that some pertinent quantities, estimates of which exist and are deemed reliable, are conserved in some sense during approximations [START_REF] Bhartia | Associate linear system approach to nonlinear random vibration[END_REF][START_REF] Ziegler | Nonstationary random vibrations of yielding beams[END_REF][START_REF] Irschik | Nonstationary random vibrations of yielding multi-degree-of-freedom systems: Method of effective envelope functions[END_REF]. This has been shown to be the case for linearization procedures with random parameters [START_REF] Bouc | The power spectral density of response for a strongly non-linear random oscillator[END_REF][START_REF] Bellizzi | Analysis of multi-degree of freedom strongly non-linear mechanical systems with random input: Part II: Equivalent linear system with random matrices and power spectral density matrix[END_REF][START_REF] Spanos | Stochastic averaging of preisach hysteretic systems[END_REF]. This section develops an average power conservation correction as a post-processing step to the linearization procedures with amplitude-dependent parameters. An original contribution of this section lies in (i) partitioning the power into the partitioning S ẋ(ω | a) into components S y (ω | a) and S u (ω | a) by analogy with Section 3:

S y (ω | a) = σ 2 ω 2 (ω 2 -B(a)) 2 + Ā(a) 2 ω 2 , (77)
S u (ω | a) = σ 2 A(a) 2 (ω 2 -B(a)) 2 + Ā(a) 2 ω 2 , (78)
and (ii) proposing an average power conservation for each component. The average power conservation correction does not modify the total average power but redistribute it between the two components. It takes advantage of the availability of an estimate provided by stochastic averaging to the average power of the cyclic component.

Average power conservation for the oscillatory velocity PSD

In contrast with what was observed with the Maxwell system (see the remark at the end of Section 3), E y 2 | a differs whether it is computed using (1/(2π))

2π 0 a 2 Ω 2 (sin Φ) 2 dΦ or (1/(2π)) ∞ -∞ S y (ω | a)
dω with S y given in (77). The former value is a natural estimate to E y 2 | a within the framework of stochastic averaging. To ensure that this value is conserved through linearization, S y (ω | a) is scaled by an amplitude-dependent factor G(a) 2 which is determined by enforcing the average power conservation in the sense of stochastic averaging [START_REF] Bouc | The power spectral density of response for a strongly non-linear random oscillator[END_REF]:

1 2π ∞ 0 σ 2 G(a) 2 ω 2 (ω 2 -B(a)) 2 + Ā(a) 2 ω 2 dω = E y 2 | a = 1 2π 2π 0 a 2 Ω 2 (sin Φ) 2 dΦ. ( 79
)
This forces It is desirable to derive a correction mechanism for S u similar to the one implemented for S y in the previous subsection. However, unlike E y 2 | a , for which stochastic averaging provides an estimate (viz., E y 2 | a = a 2 Ω 2 /2), E u 2 | a has no natural a priori estimate. An estimate to E u 2 is, however, provided by the following relation:

G(a) 2 = a 2 Ω(a) 2 Ā(a) σ 2 . ( 80
E u 2 = E ẋ2 -E y 2 , ( 83
)
where E ẋ2 and E y 2 are given by

E ẋ2 = σ 2 2 ∞ 0 1 Ā(a) 1 + A(a) 2 B(a) p a (a) da, (84)
E y 2 = 1 2 ∞ 0 a 2 Ω(a) 2 p a (a) da. ( 85
)
Figure 4 shows the (asymptotic) moments E ẋ2 given in ( 84), E y 2 given in (85), and their difference, E u 2 , as functions of the excitation intensity √ σ in the case of the indirect linearization procedure.

Eq. (83) allows one to scale S u by a factor χ 2 . But this time, the factor χ 2 is independent of a, as the average power conservation is written between a-independent quantities. Writing E u 2 with account taken of the factor χ 2 results in

E u 2 = χ 2 ∞ 0 1 2π ∞ -∞ S u (ω | a) dω p a (a) da, ( 86
) = χ 2 σ 2 ∞ 0 A(a) 2 2 Ā(a)B(a)
p a (a) da.

(87)

The average power conservation amounts to adjusting χ 2 so that

E ẋ2 -E y 2 = χ 2 σ 2 ∞ 0 A(a) 2 2 Ā(a)B(a) p a (a) da. ( 88
)
Hence

S u (ω) = χ 2 σ 2 ∞ 0 A(a) 2 (ω 2 -B(a)) 2 + Ā(a) 2 ω 2 p a (a) da, (89)
where

χ 2 = E ẋ2 -E y 2 σ 2 ∞ 0 A(a) 2 2 Ā(a)B(a)
p a (a) da .

(90)

Consequences on the velocity PSD

Combining (82) and (89) yields the following approximate velocity PSD:

S ẋ(ω) = ∞ 0 Ā(a)a 2 Ω(a) 2 ω 2 + χ 2 σ 2 A(a) 2 (ω 2 -B(a)) 2 + Ā(a) 2 ω 2 p a (a) da. ( 91
)
By design, the average power conservation correction does not alter the velocity average power predicted by either the direct or indirect linearization procedures. There is therefore no point in discussing the application of this correction to the direct linearization procedure, as it will not repair a critical defect of the procedure, namely, its poor estimates to the average power (see Subsection 7.2.). As for the indirect linearization procedure, applying the average power conservation correction (i) enhanced the (already good) predictions for the peak frequency of the PSD, its peak value, and its bandwidth, as shown in figure 5, and (ii) yielded PSDs that are in good agreement with the simulation PSDs around the main frequency, as shown in Figure 6.

Conclusions

Two original linearization procedures were developed for the EPP system. Their predictions, along with those of the classical Gaussian linearization and of the linearization described in [START_REF] Bouc | Drifting response of hysteretic oscillators to stochastic excitation[END_REF], were evaluated against Monte Carlo simulation in terms of the velocity average power, the peak frequency, the peak value, the bandwidth, and the overall shape around the main frequency. The following conclusions can be drawn:

1. The classical Gaussian linearization procedure gave good predictions of the peak frequency for small-strength inputs and performed acceptably well in estimating the velocity average power. But it performed very poorly in estimating the peak of the velocity PSD. 2. The non-Gaussian linearization procedure performed well in predicting the velocity average power over the whole input level range, the peak of the PSD and the bandwidth for smallstrength inputs. However, the procedure performed poorly in estimating the peak frequency. 3. The direct linearization procedure, despite being the most natural, poorly predicted the velocity average power, the peak frequency and the peak of the velocity PSD. The application of the average power conservation correction did not improve the predictions of the procedure. 4. The indirect linearization procedure, which combined the stochastic averaging principle, an approximating nonlinear hysteretic system with amplitude-dependent loop parameters, a linearization with amplitude-dependent parameters using a correlation-function-based criterion, the conditional PSD notion, and an average power conservation correction, was the best performer. It performed well in estimating the velocity average power and the PSD peak over the whole input level range; it performed acceptably in estimating the bandwidth over the whole input level range and the peak frequency over almost the whole input level range. The procedure qualitatively and quantitatively reproduced the simulation velocity PSDs in the main frequency band, over an input level range of almost one order of magnitude.

5. A pertinent estimation to the equivalent damping of the EPP system was obtained (Eq. ( 70)). 6. All the quantities occurring in the three linearization procedures are based on the two first Fourier coefficients of the quasi-static periodic response of the system to a constant-amplitude periodic input. To apply any of the three procedures to a given randomly excited symmetric nonlinear hysteretic system, it suffices to have the two first Fourier coefficients of the corresponding periodic response and proceed as described here for the EPP system.

List of Figures

y

  |y| zαyzβy 2 2 dΦ. (65) Solving (65) under the first-order approximation, z ≈ h c cos Φ + h s sin Φ, gives α(a) = 32 15π aΩ and β(a) = -32 45π h s . (66) Inserting decomposition (63) into (61) and eliminating c z ẋ from (60)-(61) gives the velocity correlation equation of a Maxwell system with parameters A = A * α and B = B * -A * β. Combining these with (58) and (66) yields

2 (ω 2 - 0 Ā(a)a 2 Ω(a) 2 ω 2 (ω 2 -) 8 . 2 .

 22022282 ) Now, combining (77) and (80) gives S y (ω | a) = Ā(a)a 2 Ω(a) 2 ω B(a)) 2 + Ā(a) 2 ω 2 . (81) Finally, averaging over the pdf of a yields S y (ω) = ∞ B(a)) 2 + Ā(a) 2 ω 2 p a (a) da. (82Average power conservation for the drift rate PSD

1 1Figure 1 :ω

 11 Figure 1: Restoring force diagram.

Figure 2 :Figure 3 :

 23 Figure 2: Two-sided simulation velocity PSDs in the circular frequency range [0, 2] for √ σ = 0.1 to 0.9 with increment 0.2.

Figure 4 :Figure 5 :Figure 6 :

 456 Figure 4: The moments E ẋ2 given in (84), E y 2 given in (85), and their difference, E u 2 as functions of the excitation intensity √ σ for the indirect linearization.

1

  Restoring force diagram. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 30 2 Two-sided simulation velocity PSDs in the circular frequency range [0, 2] for of the average power E ẋ2 (A), of the peak frequency of the velocity PSD S ẋ (B), of the peak value of the velocity PSD S ẋ (C), and of the bandwidth (D) as determined from Monte Carlo simulations ( ) and as predicted by the amplitude-dependent parameter indirect linearization after application of the average power correction ( ). . . . . . . . . . . . . . . . . . . . . . . . . . . . . 34 6 Two-sided velocity PSD S ẋ after application of the average power correction

				√ σ =
	3	0.1 to 0.9 with increment 0.2. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 31 Dependence on √ σ of the average power E ẋ2 (A), of the peak frequency of the
		velocity PSD S ẋ (B), of the peak value of the velocity PSD S ẋ (C), and of the band-
		width (D) as determined from Monte Carlo simulations ( ) and as predicted by the
		four procedures: (1) constant parameter Gaussian stochastic linearization (Classi-
		cal statistical linearization) (	), (2) constant parameter non Gaussian stochas-
		tic linearization (	), (3) amplitude-dependent parameter direct linearization
		(	), (4) amplitude-dependent parameter indirect linearization (	). . . . . 32
	4 5	The moments E ẋ2 given in (84), E y 2 given in (85), and their difference, E u 2 , as functions of the excitation intensity √ σ for the indirect linearization. . . . . . . 33 Dependence on √ σ