# Random Response and Crossing Rate of Fractional Order Nonlinear System with Impact

Authors: Yang Zheng; Jiayin Wu; Ming Xu
Journal: Shock and Vibration
 (2021)

---

## Abstract

The random response and mean crossing rate of the fractional order nonlinear system with impact are investigated through the equivalent nonlinearization technique. The random additive excitation is Gaussian white noise, while the impact is described by a phenomenological model, which is developed from the actual impact process experiments. Based on the equivalent nonlinearization technique, one class of random nonlinear system with exact probability density function (PDF) solution of response is selected. The criterion of the appropriate equivalent nonlinear system is the similarity with the original system on the damping, stiffness, and inertia. The more similar, the higher the precision. The optimal unknown parameters of the equivalent random nonlinear system in the damping and stiffness terms are determined by the rule of smallest mean-square difference. In the view of equivalent nonlinearization technique, the response of the original system is the same as that of the equivalent system with the optimal unknown parameters in analytical solution manner. Then, the mean crossing rate is derived from stationary PDF. The consistence between the results from proposed technique and Monte Carlo simulation reveals the accuracy of the proposed analytical procedure.


---

## Body

## 1. Introduction

The impact can be found frequently in the structural and mechanical engineering, where there exists clearance between two bodies [1, 2]. Due to the discontinuity, the dynamical behavior of the vibration system comprising impact is different from the usual smooth system and becomes extremely complex [3–5]. The impact process is usually modeled as a velocity jump at the time of impact occurring, which is expressed by x˙+=−rx˙− with 0≤r≤1. The rebound velocity x˙+ equals to r by impact velocity x˙− and has an opposite direction of movement. The mapping of the velocity before and after impact is not continuous both in the size and orientation. Many research studies are devoted to reveal the dynamical behaviors of vibro-impact system [6, 7]. To deal with the velocity jump at collision moment, several effective transformations of state variables are always helpful [7–9]. Based on these transformations, the vibro-impact system with viscoelastic damping and nonzero offset barrier is studied by the stochastic averaging [10]. Multivalue response of a nonlinear vibro-impact system under narrow-band excitation is investigated [11]. The triple-valued response under a certain case is observed, which may have two or four steady-state solutions. Other research studies with respect to the vibro-impact system can be found in [12–15].Based on the impact experiments of elastic-plastic structures, a phenomenological impact model is developed. Due to time duration of the elastic-plastic deformation process, the impact model of the velocity jump at collision moment is not appropriate in this case [16, 17]. By using this phenomenological impact model, the stochastic averaging method has been adopted to obtain the stationary probability density function (PDF) of the vibro-impact system subjected to Gaussian white noises [18]. Due to the limitation on slight nonlinearity and weak excitations density, the equivalent nonlinearization technique has a remarkable advantage, i.e., it can be applicable to systems which are strongly nonlinearity. This technique has been developed and used to complex nonlinear system in Hamiltonian framework [19–23]. Considering the significant advantage, the equivalent nonlinearization technique has been successfully adopted to evaluate the reliability and the reliability-based design of inelastic structure, and a simple equivalent nonlinear system that retains the dynamic characteristics of the first two modes and the global yielding behavior has been developed to replace the original system [24]. It can be expected that the equivalent nonlinearization technique is applicable to deal with the random system-incorporated inelastic impact described by the phenomenological impact model.Besides, the mean crossing rate, which always relates to predicting the extreme response statistics and the system reliability, is another important quantity which should be considered. The evaluation of the mean crossing rate has been of great interest among researchers. Naess and Karlsen [25] calculated the level crossing rate of second-order stochastic Volterra systems. Naess et al. [26] applied Monte Carlo simulation to predict the extreme response statistic of floating offshore structures subjected to random seas. Beck and Melchers [27] investigated time variant reliability of uncertain structure by using the ensemble crossing rate. As illustrated in [28], the mean crossing rate for a stationary stochastic process can be determined by the joint probability density of the stochastic process and its derivative process. It means that the equivalent nonlinearization technique is applicable to derive not only the system response but also the mean crossing rate.The random response and mean crossing rate of the nonlinear system with fractional order stiffness and impact are studied by using the equivalent nonlinearization technique. The paper is organized as follows. In Section2, the fractional order nonlinear system with impact is established and the impact model is described. In Section 3, an equivalent nonlinearization technique is adopted and the approximate solution of system response is derived. In Section 4, the mean crossing rate is introduced and the analytical formula is expressed by the response PDF of system. In Section 5, examples are given to shown the proposed analytical technique and different system parameters are discussed. The conclusions are drawn in Section 6.


## 2. Fractional Order Nonlinear System with Impact

The present paper is concerned with a vibro-impact system with right-side barrier under additive excitations described by Gaussian white noise, as shown in Figure1. The motion equation governing the mechanical behavior of the system is(1)X¨+2ζX˙+k1X+k3XXα−1+gX,X˙=Wt,where X is the system displacement and overhead dot indicates differentiation with respect to time t, ζ is the viscous damping coefficient, and Wt is Gaussian white noise in sense of Stratonovitch [29] with zero mean and correlation function Rτ=2Dδτ.Figure 1
Schematic of a vibrating system with a barrier on the right side.The system stiffness here is fractional order nonlinear. The physical valueα is larger than 1 and not an integer for many materials due to the nonlinear stress-strain relationship [30]. This nonlinear stiffness has been studied by Cvetivanin and Zukuvic [31, 32]. In Figure 1, we consider the elastic force obeys the fractional order law. k1 and k3 in system (1) are the linear and nonlinear coefficients, respectively. Specially, parameter k1 can be positive for a beam, while negative for a beam with sufficient axial load [33].Impact forcegx,x˙, which depends on both the velocity and displacement, is governed by the phenomenological impact model when collision occurs [17]. This impact model is efficient to describe the collision for the elastic-plastic materials. The distance of clearance between the mass and barrier is δr.The impact forcegx,x˙ will be given for two different cases (the dissipative collision and the conservative collision), respectively. The criterion to distinguish the dissipative and conservative collision is the relation of the maximal right-side displacement a1 of the impact mass and the elastic limit x0. The dissipative collision for a1>x0, while conservative collision for a1<x0


### 2.1. For Dissipative Collision

The impact forcegx,x˙ can be explicitly derived for the loading phase g+ and unloading phase g− as follows:(2a)g+=k2x−δrp,x∈δr,a1,0,x∈−a2,δr,forx˙>0loading phase,(2b)g−=s2|x−δr−Δr|q,x∈δr+Δr,a1,0,x∈−a2,δr+Δr,forx˙<0unloading phase,with(3a)Δr=spa1−x0,(3b)s2=Fma1−δr−Δrq,(3c)Fm=k2a1−δrp,where sp,p,andq are material parameters that are determined through experiments, Δr is the plastic deformation induced by the dissipative impact, k2 and s2 represent the stiffness coefficients of the loading phase and unloading phase, respectively, and a2 represents the maximal left-side displacement of the impact mass. The relation of impact force versus relative displacement for dissipative collision is depicted in Figure 2.Figure 2
The constitutive relationship of the phenomenological impact model.


### 2.2. For Conservative Collision

The impact forcegx,x˙ becomes simple and only depends on the system displacement. The impact forces gx,x˙ can be written as(4)g=k2x−δrp,x∈δr,a1,0,x∈−a2,δr.


## 2.1. For Dissipative Collision

The impact forcegx,x˙ can be explicitly derived for the loading phase g+ and unloading phase g− as follows:(2a)g+=k2x−δrp,x∈δr,a1,0,x∈−a2,δr,forx˙>0loading phase,(2b)g−=s2|x−δr−Δr|q,x∈δr+Δr,a1,0,x∈−a2,δr+Δr,forx˙<0unloading phase,with(3a)Δr=spa1−x0,(3b)s2=Fma1−δr−Δrq,(3c)Fm=k2a1−δrp,where sp,p,andq are material parameters that are determined through experiments, Δr is the plastic deformation induced by the dissipative impact, k2 and s2 represent the stiffness coefficients of the loading phase and unloading phase, respectively, and a2 represents the maximal left-side displacement of the impact mass. The relation of impact force versus relative displacement for dissipative collision is depicted in Figure 2.Figure 2
The constitutive relationship of the phenomenological impact model.


## 2.2. For Conservative Collision

The impact forcegx,x˙ becomes simple and only depends on the system displacement. The impact forces gx,x˙ can be written as(4)g=k2x−δrp,x∈δr,a1,0,x∈−a2,δr.


## 3. Equivalent Nonlinearization Technique

As mentioned above, the equivalent nonlinearization technique has some advantages compared with the stochastic averaging technique, especially for the effectiveness to the strong nonlinearity and its concision. To utilize the equivalent nonlinearization technique, the first and critical step is the selection of the equivalent nonlinear system family. The original system solution is approximately expressed by that of the equivalent system family. The more similar, the higher the accuracy. There are no rigorous rules for selecting the equivalent system family. General speaking, an efficient rule is to make inertia, stiffness, and damping of the original and equivalent system close. The above rules, which confine the selection of the equivalent nonlinear system, guarantee the accuracy of the results from the equivalent nonlinearization technique.The phenomenological impact model is depicted for two different cases. In the conservative collision case, the loading path coincides with the unloading path, while in the dissipative collision case, the loading path is deviated from the unloading path beyond the elastic limit and the nonzero area encircled by the two paths indicates the dissipative mechanism of the inelastic impact. In other words, in case of the dissipative collision case, the impact is not only to store the potential energy due to its elasticity but also to dissipate energy as a damping. So, in the selection of the equivalent nonlinear system family, the conservative and dissipative components of inelastic impact should be considered.It is reasonable to reflect the conservative component of inelastic impact based on the loading path:(5)g¯x=k2x−δrp,x≥δr,0,x<δr.In the case of the conservative collision, equation (5) is the exact expression of the impact process, but it is not suitable for the dissipation collision. In order to reflect the inconsistency of loading and unloading paths, one intuitive and convenient selection of equivalent stiffness to reflect the conservative component of inelastic impact is 1+Seg¯x, in which Se is a correction coefficient.Obviously, the linear damping is the simplest reflection of the dissipative component of the inelastic impact, but is not appropriate since the dissipative component will not play a role in conservative collision and depend on the system state in dissipative collision. To avoid the irrationality of linear damping, quasi-linear damping with damping coefficient depending on the system states is a good alternative and can be selected asζeHX˙, where H=Hx,x˙=x˙2/2+Gx and the potential energy Gx=∫0xk1x+k3xxα−1+1+Seg¯xdx. The selection of ζeHX˙ has two remarkable advantages: (i) due to the similar properties between the impact process and hysteretic behavior [34], ζeHX˙ can reflect the dissipative component well and (ii) the solvability of the equivalent nonlinear system including ζeHX˙ is guaranteed.Based on the above analysis, the original system is replaced by the nonlinear system family,(6)X¨+2ζ+ζeHX˙+k1X+k3XXα−1+1+Seg¯X=Wt,which has the stationary PDF with normalization constant N as(7)psx,x˙=Nexp−1D∫0e2ζ+ζeudu|e=x˙2/2+Gx.In equivalent nonlinear system (6), the constant Se and the function ζeH are both unknown and should be determined through some criterions. Unfortunately, to derive the optimal function ζeH is almost impossible. A practical way is to expand the function ζeH in a power series:(8)ζeH=∑i=1∞biH−H∗i,H≥H∗,0,H<H∗,where H∗ denotes the system energy critical value of elastic impact and inelastic impact, i.e., H∗=1/2k1x02+1/α+1k3x02x0α−1+k2/1+px0−δr1+p, and bi are undetermined constant parameters. Now, the equivalent system of the original system comes down to determine unknown parameters Se and bi,i=1,2,….The criterion of minimizing the mean-square valueEe2 is selected to derive the unknown parameters Se and bi, where e=ζeHX˙+1+Seg¯X−gX,X˙ and Eg is the expectation operator. Thus, the optimal parameters Seandbi can be solved from the following equations:(9)∂∂SeEe2=0,∂∂biEe2=0,i=1,2,…Because of the infinite number ofbi, equation (9) is unsolvable and should be truncated. If only b1 is kept and bii≥2 are ignored, parameters Se and bi can be solved from equation (9) and can be written as(10)b1=EH−H∗X˙gX,X˙ℏH−H∗EH−H∗2X˙2ℏH−H∗,Se=EgX,X˙g¯X−Eg¯2XEg¯2X,where hg represents Heaviside step function.In equations (9) and (10), the joint PDF depends on the unknown parameters Se and bi, so the undetermined parameters Se and bi should be obtained through iterative technique.


## 4. Stochastic Response and Mean Crossing Rate Evaluation

As discussed in the above section, by substituting the parametersSeandbi into the PDF in equation (7), the joint PDF of the original system can be approximated to that of the equivalent system. Then, the PDF of original system states X and X˙, the marginal PDF, can be written directly.Except the stochastic responses analysis, another quantity of particular importance is the mean crossing rate, which associates with the extreme response prediction, system safety, and reliability. The mean crossing rate, denoted asνact, can be described as the mean number of a stochastic process Yt crossing a critical value ac per unit time. The mean crossing rate νact which usually depends on the time except Yt is a stationary stochastic process. Particularly, in the special case of Yt which is a mean-square differentiable stationary stochastic process with continuous time and state, the mean crossing rate of Yt can be defined as [28, 35](11)νac=∫−∞∞y˙pac,y˙ⅆy˙,where υac is independent of the time and ac represents the critical excursion value of system displacement, pac,y˙=py,y˙|y=ac, in which py,y˙ is the stationary joint probability density of stochastic variances Yt and Y˙t. Similar to the mean crossing rate, the mean up-crossing and down-crossing rates, which represent the mean number of a stochastic process Yt up-crossing and down-crossing a critical value ac of system displacement per unit time, respectively, can be defined as(12a)νac+=∫0∞y˙pac,y˙ⅆy˙,(12b)νac−=∫−∞0−y˙pac,y˙dy˙.Sincepy,y˙ is an even function of stochastic variance Y˙, the mean up-crossing and down-crossing rate satisfies(13)νac+=νac−=12νac.Considering the relation in equation (13), only the mean up-crossing rate is discussed in this paper.


## 5. Numerical Results and Discussion

Some numerical calculations are carried out to validate the proposed analytical technique. System parameters are selected ask1=1, k2=5,k3=1,α=1.5,δr=0.1,x0=0.125,2D=0.04,ζ=0.01,sp=0.2,p=1.39, and q=1.8 [17], unless otherwise mentioned. Figures 3(a) and 3(b) relate the displacement PDF px and the velocity PDF px˙ under different fractional orders α. It can be clearly seen that the analytical results (solid lines) agree with the Monte Carlo simulation (MCS) results of the original system (1) (circles). As the fractional order α increases from α=1.5 to α=3.5 and α=5.5, the displacement PDF is to be more flat and the velocity PDF is slightly changed. It means that the displacement is sensitive to the fractional order α, while the velocity is not. Figures 4(a) and 4(b) plot the displacement and velocity PDFs of the nonlinear system under different Gaussian white noise intensities 2D. The bigger intensity 2D means large displacement and velocity response, which is consistent with our intuition. The influences of the nonlinear stiffness coefficients k3 are shown in Figures 5(a) and 5(b). The high nonlinear stiffness k3 restricts the system displacement response. However, such restriction has little effects on the velocity response, where the curves of the velocity PDF almost coincide. Then, the bistable system (k1=−1) is studied. The displacement and velocity PDFs of the nonlinear system under different fractional order α is given in Figures 6(a) and 6(b). The displacement PDF has two peaks, while velocity PDF only has one peak. The system with small fractional order α (i.e., α=1.5) has a flat PDF. Figures 7(a) and 7(b) are the curves of the displacement and velocity PDFs. The responses increases with the increased excitation intensity 2D. The influences of the nonlinear stiffness coefficients k3 on the displacement and velocity PDFs are shown in Figures 8(a) and 8(b). Strong nonlinear stiffness coefficient k3 reduces the displacement response.Figure 3
The PDF of the nonlinear system under a different fractional orderα. (a) The displacement PDF; (b) the velocity PDF (k1=1, solid lines: the analytical results, and circles: MCS results).
(a)(b)Figure 4
The PDF of the nonlinear system under different Gaussian white noise intensities2D. (a) The displacement PDF; (b) the velocity PDF (k1=1, solid lines: the analytical results, and circles: MCS results).
(a)(b)Figure 5
The PDF of the nonlinear system under different nonlinear stiffness coefficientsk3. (a) The displacement PDF; (b) the velocity PDF (k1=1, solid lines: the analytical results, and circles: MCS results).
(a)(b)Figure 6
The PDF of the nonlinear system under a different fractional orderα. (a) The displacement PDF; (b) the velocity PDF (k1=−1, solid lines: the analytical results, and circles: MCS results).
(a)(b)Figure 7
The PDF of the nonlinear system under different Gaussian white noise intensities2D. (a) The displacement PDF; (b) the velocity PDF (k1=−1, solid lines: the analytical results, and circles: MCS results).
(a)(b)Figure 8
The PDF of the nonlinear system under different nonlinear stiffness coefficientsk3. (a) The displacement PDF; (b) the velocity PDF (k1=−1, solid lines: the analytical results, and circles: MCS results).
(a)(b)Besides the system responses, the mean up-crossing rate plays an important role in the system reliability. Figure9 gives relation of the mean up-crossing rate υac+ and the critical value of excursion ac and fractional order α. The results from the proposed technique agrees well with those from the Monte Carlo simulation(MCS). To give a clear view, Figure 10 plots the curve of the mean up-crossing rate υac+ with critical value of excursion ac under α=1.5. The mean up-crossing rage monotonously decrease with the excursion ac. The mean up-crossing rate υac+ of the bistable system is plotted in Figures 11 and 12. The mean up-crossing rate υac+ increases first and then decreases with the increased critical value of excursion ac. The comparison between the proposed technique and the Monte Carlo simulation (MCS) shows that the accuracy is not high, but acceptable. These can be clearly observed in Figure 12, where the fractional order α equals to 5.5.Figure 9
The mean up-crossing rateυac+ vs. critical value of excursion ac and fractional order α (k1=1). (a) The analytical results. (b) MCS results.
(a)(b)Figure 10
The mean up-crossing rateυac+ vs. critical value of excursion ac under α=1.5 (k1=1, solid lines: the analytical results, and circles: MCS results).Figure 11
The mean up-crossing rateυac+ vs. critical value of excursion ac and fractional order α (k1=−1). (a) The analytical results. (b) MCS results.
(a)(b)Figure 12
The mean up-crossing rateυac+ vs. critical value of excursion ac under α=5.5 (k1=−1, solid lines: the analytical results, and circles: MCS results).


## 6. Conclusions

The random response and mean crossing rate of the nonlinear system with fractional order stiffness and impact has been investigated through equivalent nonlinearization technique. The impact is described by an empirical model developed from impact experiment of elastic and plastic materials. By using the equivalent nonlinearization technique, the original vibro-impact system is equivalently replaced by a nonlinear system. Through minimizing the mean-square value of the system difference, one optimal system is chosen from the equivalent nonlinear system family. Then, the joint PDF of system displacement and velocity are analytically obtained through the equivalent nonlinear system. The agreement between the analytical results and MCS validates effectiveness of the proposed technique. The proposed technique is also adopted to derive the mean crossing rate of the vibro-impact system, and the acceptable precision of the results illustrates the effectiveness of the proposed technique to evaluate the mean crossing rate. It is necessary to emphasize that, in comparison with the stochastic averaging technique, the present technique can extend the applicable range of system parameters to which the stochastic averaging technique is invalid.



---

*Source: 6666843-2021-05-27.xml*