Stochastic Bifurcation of a Strongly Non-Linear Vibro-Impact System with Coulomb Friction under Real Noise

Abstract

This manuscript investigated the response of a strongly non-linear vibro-impact (VI) system with Coulomb friction. The impact model is used with classical impact. The excitation is modelled by real noise. First, the VI system is converted into a simplified system without any barrier by non-smooth transformation (symmetric transformation). The stochastic averaging method is adopted to obtain the theoretical stationary probability function of the VI system. Next, the Duffing Van der Pol VI system with Coulomb friction is used to verify the validity of the proposed theoretical method compared with numerical simulations. Moreover, the influence of bandwidth, noise intensity, and friction amplitude are further analyzed in detail on the probability density function (PDF) of distribution of the VI system. The P-bifurcation is studied by a qualitative change of friction amplitude and restitution coefficient on the stationary probability distribution, which indicated that these parameters can arouse the emergence of stochastic P-bifurcation.

1. Introduction

1.1. Background

In the field of science and engineering, vibro-impact (VI) systems can model some relevant problems such as ship roll motion against icebergs, rotor-stator rubbing in rotating machinery, or impact interaction in pipe-baffle interfaces due to seismic excitations and so on [1]. Moreover, the VI systems displayed some interesting phenomena, such as grazing bifurcation [2], torus bifurcation [3], and chatter and sticking [4]. Some modelling techniques are developed to model impact, such as the Herz contact law, for classical impact with an instantaneous velocity jump. Because excessive impacts may cause damage to some structural components in practice, in order to avoid harmful impact or some random vibro-impact, research on the response and bifurcation of vibro-impact has received much attention.

1.2. Literature Survey

Namachachivaya [5,6] studied the dynamic behaviors of VI systems under random perturbation. Huang [7] investigated the stationary response for stochastically excited multi-degree-of-freedom VI systems. Xu [8] studied the stochastic response of VI systems with inelastic contact. Rong [9] applied multiple scales to study the amplitude-frequency response of VI systems. Yang [10] investigated the random response of Rayleigh VI systems excited by Poisson white noise. Zhu [11,12] studied the random response of VI systems through the exponential-polynomial closure method.

As another representative non-smooth factor, friction, which exists in the mechanical system, has been extensively researched in the last two decades. For example, Green [13] analyzed the effect of friction on the response of an electromagnetic energy harvester with Duffing-type non-linearities. Sun [14] applied generalized cell mapping to investigate the random response of the Coulomb friction system. Pankaj [15] studied the stochastic bifurcation of a Duffing oscillator with Coulomb friction under Poisson white noise. Sun [16] researched the reliability of a non-linearly damped oscillator with Coulomb friction under Gaussian white noise excitation. If the system response is required to be below a critical level for safety consideration, a control device may be needed. Rigatos [17,18] studied the sensorless control and Bryson [19] considered the optimal control; Pappalardo [20,21] investigated the algorithm controlling the mechanical system.

1.3. Formulation of the Problem

In addition, in some mechanical systems, such as rotor-stator interactions in turbines, interacting gear teeth, riveted or bolted supports in structural systems [22,23], vibro-impact, and friction dynamics arise simultaneously when some rigid subparts of mechanical systems move separately in partial-motion restraining mechanisms. The friction and impact between gear wheels and railways usually induces the “chatter” phenomenon, which tends to generate undesired noise and excessive vibrations which could influence performance and even cause damage to the system.

Since random factors in most real-world circumstances may influence the dynamic behavior of VI systems, some researchers studied VI systems under random excitation. For mathematical simplicity, the random driving force is generally assumed as zero mean Gaussian white noise. However, white noise has an idealization of the correlation time scale, while practical physical models need to consider finite correlation times, so the white noise does not describe random disturbances in nature. Wide-band noise as a finite correlation of time stochastic excitation has received increasing attention. Zhu [24] studied the dynamic behaviors of strongly non-linear systems under broad noise by averaging method; Lin and Zhu [25,26] researched the response and stability for wide-band noise for an excited viscoelastic system.

1.4. Contribution of this Study

As mentioned above, the random vibration of VI systems has attracted much attention in the past. However, the investigation hasn’t considered friction in the VI system. Moreover, the existence of friction in the VI system is unavoidable, and the random excitation modeled by real noise is more reasonable. Bifurcation analysis under different system parameters is essential to avoid the appearance of undesirable dynamic behaviors. Thus, it is significant to investigate the stochastic response and bifurcation of the VI system with Coulomb friction subjected to real noise excitation.

1.5. Organization

In this paper, the stochastic response and bifurcation of a strongly non-linear VI system with Coulomb friction under real noise excitation is investigated. The model of classical impact with instantaneous velocity jump is used in this paper. In the first step, the VI system is converted into the simplified system without any barrier by non-smooth transformation. In the second step, the stationary probability density function (PDF) is obtained through the stochastic averaging method. The effects of bandwidth, noise intensity of wideband noise, and friction parameters on stationary probability density are considered. Lastly, stochastic P-bifurcation is studied by the variation of friction parameters and restitution coefficient on the stationary probability distribution.

2. System Description and Non-Smooth Transformation

2.1. System Description

Consider the strongly non-linear VI system with Coulomb Friction (Figure 1) subjected to real noise excitation

$$\{\begin{array}{}\mathrm{(1)}& \ddot{x}+\epsilon c(x,\dot{x})\dot{x}+g(x)+\mu \mathrm{sgn}(\dot{x})={\epsilon }^{\raisebox{1ex}{$1$}\!\left/ \!\raisebox{-1ex}{$2$}\right.}{f}_i(x,\dot{x}){\xi }_i(t)x>b,\hfill \mathrm{(2)}& {\dot{x}}_{+}=-r{\dot{x}}_{-}x=b.\hfill \end{array}$$

where $\ddot{x},\dot{x}$ and $x$ denote the acceleration, velocity, and displacement. Here $c(x,\dot{x})$ is the damping coefficient, $g(x)$ is the strongly non-linear stiffness and is an odd function. $b$ denotes distance of the unilateral barrier from the static equilibrium position of the system and we assigned zero to $b$ in this paper. ${\dot{x}}_{+}$ and ${\dot{x}}_{-}$ denote impact velocities before and after impact. $\epsilon$ is a small scale parameter, and $r(0<r\le 1)$ is the restitution coefficient. ${\xi }_i(t)$ represents real noise with correlation function $R_{ij}(\tau )$. Coulomb friction $\mu \mathrm{sgn}(\dot{x})$ has constant amplitude $\mu$ and $\mathrm{sgn}(\cdot )$ denotes sign function

$$\mathrm{sgn}(\dot{x})=\{\begin{array}{l}1\dot{x}>0\\ 0\dot{x}=0\\ -1\dot{x}<0\end{array}.$$

Although Coulomb friction is intuitively simple, its discontinuity added difficulty in trying to analytically evaluate the system response.

2.2. Non-Smooth Transformation

The procedure of non-smooth transformations mentioned in Ref. [27] are introduced as follows

$$x=\left|y\right|=x_1,\dot{x}=\dot{y}\mathrm{sgn}(y)=x_2$$

(3)

The transformation of Equation (3) maps the domain $x>0$ of original phase plane $\left(x,\dot{x}\right)$ onto the whole phase plane $\left(y,\dot{y}\right)$; the transformed equations of new variables can be written as:

$$\{\begin{array}{l}\ddot{y}+\epsilon c(y,\dot{y})\dot{y}+g(y)+\mu \mathrm{sgn}(\dot{x})={\epsilon }^{\raisebox{1ex}{$1$}\!\left/ \!\raisebox{-1ex}{$2$}\right.}{f}_i(y,\dot{y}){\xi }_i(t)\mathrm{sgn}(y)y(t_{•})\ne 0,\\ {\dot{y}}_{+}=r{\dot{y}}_{-}y(t_{•})=0\end{array}$$

(4)

Through the new equation, one can find that the jump of transformed velocity $\dot{y}$ changes proportional to $(1-r)$ instead of $(1+r)$ for the original velocity $\dot{x}$. By virtue of $y(t_{\ast })=0$ and $\delta (t-t_{\ast })=\left|\dot{y}(t)\right|\delta (y(t))$ at impact, Equation (2) can be considered as an additional impulse through Dirac delta function and the impulsive term can be written as

$$({\dot{y}}_{-}-{\dot{y}}_{+})\delta (t-t_{\ast })=(1-r)\dot{y}\left|\dot{y}\right|\delta (y)$$

(5)

The term $(1-r)\dot{y}\left|\dot{y}\right|\delta (y)$ tends to depict the impact losses of the system. Thus, Equation (1) and (2) can be combined into one equation by adding the impulsive damping term $(1-r)\dot{y}\left|\dot{y}\right|\delta (y)$ through the transformation (5)

$$\ddot{y}+\epsilon (f(y,\dot{y}))\dot{y}+\mu \mathrm{sgn}(\dot{y})+g(y)={\epsilon }^{\raisebox{1ex}{$1$}\!\left/ \!\raisebox{-1ex}{$2$}\right.}{f}_i(y,\dot{y}){\xi }_i(t)\mathrm{sgn}(y),$$

(6)

$$f(y,\dot{y})=c(y,\dot{y})+\frac{1}{\epsilon }(1-r)\left|\dot{y}\right|\delta (y),U(y)={\displaystyle \underset{0}{\overset{y}{\int }}g(u)}du.$$

3. Stochastic Averaging Procedure

We can assume that the solution of Equation (6) has the following form [25]:

$$\{\begin{array}{l}y(t)=A(t)\cos \theta (t),\\ \dot{y}(t)=-A(t)\nu (A,\theta )\sin \theta (t),\\ \theta (t)=\varphi (t)+\phi (t)\end{array}$$

(7)

where $A,\theta ,\varphi ,\phi$ are all random processes, $\nu (A,\theta )$ is the instantaneous frequency of the oscillator and has the following form:

$$\nu (A,\theta )=\frac{d\varphi }{dt}=\sqrt{\frac{2(U(A)-U(A\cos \theta ))}{A^2{\sin }^2\theta }=}{b}_0(a)+{\displaystyle \sum _{r=1}^{\infty }{b}_r(A)}\cos r\theta ,$$

(8)

Substituting Equation (7) into Equation (6), we can observe the transformation equations of amplitude $A(t)$ and phase angle $\phi (t)$:

$$\begin{array}{c}\frac{dA}{dt}=\epsilon F_1+{\epsilon }^{1/2}{G}_{1k}(A,\phi ){\xi }_k(t),\frac{d\phi }{dt}=\epsilon F_2+{\epsilon }^{1/2}{G}_{2k}(A,\phi ){\xi }_k(t),\\ F_1=\frac{-A^2{\nu }^2{\sin }^2\theta }{g(A)}f(A\cos \theta ,-A\nu \sin \theta )-\frac{A\nu \sin \theta }{g(A)}\frac{\mu }{\epsilon }\mathrm{sgn}(-A\nu \sin \theta ),\\ F_2=\frac{-A{\nu }^2\sin \theta \cos \theta }{g(A)}f(A\cos \theta ,-A\nu \sin \theta )-\frac{\nu \cos \theta }{g(A)}\frac{\mu }{\epsilon }\mathrm{sgn}(-A\nu \sin \theta ),\\ G_{1k}=-\frac{A\nu \sin \theta }{g(A)}{f}_k(A\cos \theta ,-A\nu \sin \theta ),G_{2k}=-\frac{\nu \cos \theta }{g(A)}{f}_k(A\cos \theta ,-A\nu \sin \theta ).\end{array}$$

According to the Stratonovich-Khasminskii limit theorem [28,29], the process $A(t)$ can be approximated as a Markov process and is governed by the following averaged $\mathrm{It}\widehat{\mathrm{o}}$ stochastic differential equation:

$$dA=m(A)dt+\sigma (A)dB(t)$$

(9)

where

$$m(A)=\epsilon {\langle F_1+{\int }_{-\infty }^0\left[{\frac{\partial G_{1k}}{\partial A}|}_t{G_{1j}|}_{t+\tau }+{\frac{\partial G_{1k}}{\partial \theta }|}_t{G_{2j}|}_{t+\tau }\right]R_{kj}(\tau )d\tau \rangle }_t,$$

(10)

$${\sigma }^2(A)=\epsilon {\langle {\int }_{-\infty }^{+\infty }{G}_{1k}{|}_t{G}_{1j}{|}_{t+\tau }{R}_{kj}(\tau )d\tau \rangle }_t.$$

(11)

$B(t)$ denotes a unit wiener process, ${\langle \rangle }_t$ represents the time averaging.

To obtain the explicit expressions of the averaged $m(A)$ an ${\sigma }^2(A)$ in Equation (10) and (11), one can expand $\epsilon F_1$ and ${\epsilon }^{1/2}{G}_{ik}$ into Fourier series with respect to $\theta$:

$$\{\begin{array}{l}\epsilon F_1(A,\Gamma )=F_{10}+{\displaystyle \sum _{i=1}^{\infty }}(F_{1i}^{(s)}\cos i\theta +F_{1i}^{(t)}\sin i\theta )\\ {\epsilon }^{1/2}{G}_{jk}(A,\Gamma )=G_{jk0}+{\displaystyle \sum _{i=1}^{\infty }}(G_{jki}^{(s)}\cos i\theta +F_{jki}^{(t)}\sin i\theta )\\ j=1,2;k=1,2,\cdots ,m\end{array},$$

(12)

Using Equation (12) and time averaging, we obtained the drift and diffusion coefficients given by:

$$\begin{array}{c}m(A)=F_{10}(A)+\frac{\pi }{2}{\displaystyle \sum _{i=1}^{\infty }\{\left[\frac{d{G}_{1ki}^{(s)}}{dA}{G}_{1li}^{(s)}+\frac{d{G}_{1ki}^{(t)}}{dA}{G}_{1li}^{(t)}+i(G_{1ki}^{(t)}{G}_{2li}^{(s)}-G_{1ki}^{(s)}{G}_{2li}^{(t)})\right]}{S}_{kl}(i\omega (A))\\ +\left[\frac{d{G}_{1ki}^{(s)}}{dA}{G}_{1li}^{(t)}-\frac{d{G}_{1ki}^{(t)}}{dA}{G}_{1li}^{(s)}+i(G_{1ki}^{(t)}{G}_{2li}^{(t)}+G_{1ki}^{(s)}{G}_{2li}^{(s)})\right]I_{kl}(i\omega (A))\}+\pi \frac{d{G}_{1k0}}{dA}{G}_{1l0}{S}_{kl}(0),\\ {\sigma }^2(A)=2\pi G_{1k0}{G}_{1l0}{S}_{kl}(0)+\pi {\displaystyle \sum _{i=1}^{\infty }[\left(G_{{}_{1ki}}^{(s)}{G}_{{}_{1li}}^{(s)}+G_{{}_{1ki}}^{(t)}{G}_{{}_{1li}}^{(t)}\right)}{S}_{kl}(i\omega (A))+\\ \left(G_{{}_{1ki}}^{(s)}{G}_{{}_{1li}}^{(t)}-G_{{}_{1ki}}^{(t)}{G}_{{}_{1li}}^{(s)}\right)I_{kl}(i\omega (A))]\end{array}$$

where

$$S_{kl}(\omega )=\frac{1}{\pi }{\int }_{-\infty }^0{R}_{kl}(\tau )\cos \omega \tau d\tau I_{kl}(\omega )=\frac{1}{\pi }{\int }_{-\infty }^0{R}_{kl}(\tau )\sin \omega \tau d\tau .$$

4. Response Probability Density Functions

The PDF for amplitude is governed by the following Fokker-Planck-Kolmogorov (FPK) equation:

$$\frac{\partial p}{\partial t}=-\frac{\partial }{\partial A}\left[m(A)p\right]+\frac{1}{2}\frac{{\partial }^2}{\partial A}\left[{\sigma }^2(A)p\right],$$

(13)

where the boundary condition is:

$$\begin{array}{c}p=\mathrm{finite}\mathrm{as}A=0,\\ p\to 0,\partial p/\partial A\to 0\mathrm{as}A\to \infty .\end{array}$$

According to Zhu [24], $N$ is a normalization constant and the stationary solution of FPK equation is

$$p(A)=\frac{N}{{\sigma }^2(A)}\exp \left[{\displaystyle {\int }_0^A\frac{2m(u)}{{\sigma }^2(u)}du}\right],$$

(14)

The stationary probability of $E$ can be derived from the Equation (14) as follows:

$$p(E)=p(A)\left|\frac{dA}{dE}\right|={\frac{p(A)}{g(A)}|}_{A=U^{-1}(E)},$$

(15)

where $A=U^{-1}(E)$ is the inverse function of $U(A)=E$, the stationary probability density of velocity and displacement can be further obtained from $p(E)$ as follows:

$$p(y,\dot{y})={\frac{p(E)}{T(E)}|}_{E=(1/2){\dot{x}}^2+U(x)},$$

(16)

where

$$T(E)={\frac{2\pi }{b_0(A)}|}_{A=U^{-1}(E)}.$$

Through the transformation (3), the joint stationary PDF of the original displacement $x_1$ and velocity $x_2$ can be expressed as

$$p(x_1,x_2)=2{\tilde{p}}_{y_1,y_2}(x_1,x_2)$$

The stationary PDFs of displacement $x_1$ and velocity $x_2$ are derived, respectively:

$$p_{x_1}(x_1)={\displaystyle {\int }_{-\infty }^{+\infty }p(x_1,s)}ds,p_{x_2}(x_2)={\displaystyle {\int }_{-\infty }^{+\infty }p(x_2,s)}ds.$$

5. Example

The Duffing-Van der Pol VI system with Coulomb Friction are often met in science and engineering as a representative model of a strongly non-linear system under real noise excitation, which has the following form:

$$\{\begin{array}{l}\ddot{x}+(c_2{x}^2-c_1)\dot{x}+{\omega }_0^{2}x+\alpha x^3+\mu \mathrm{sgn}(\dot{x})=f\xi (t)x>0,\\ {\dot{x}}_{+}=r{\dot{x}}_{-}x=0\end{array}$$

(17)

where $c_1,c_2,\alpha ,f$ and ${\omega }_0$ are positive parameters, $\xi (t)$ is the real noise and the spectral density has the following low-pass type

$$S(\omega )=\frac{D}{\pi }\frac{1}{{\omega }^2+{\alpha }_1^2},$$

where ${\alpha }_1$ and $D$ are the bandwidth and noise intensity of the real noise. With the above spectral density the random excitation $\xi (t)$ can be generated from a first filter

$$\dot{\xi }(t)+\alpha \xi (t)=W(t)$$

where $W(t)$ is a white noise with noise intensity $2D$. According to Equation (6), the equation of system (17) can be transformed as the following stochastic system

$$\ddot{y}+(c_2{y}^2-c_1+(1-r)\left|\dot{y}\right|\delta (y))\dot{y}+\mu \mathrm{sgn}(\dot{y})+{\omega }_0^{2}y+\alpha y^3=f\xi (t)\mathrm{sgn}(y),$$

The instantaneous frequency $\nu (A,\theta )$ in Equation (8) can be approximated as the following finite sum.

$$\nu (A,\theta )=b_0(A)+b_2(A)\cos 2\theta +b_4(A)\cos 4\theta +b_6(A)\cos 6\theta ,$$

(18)

$$\begin{array}{c}{b}_0(A)=\sqrt{({\omega }_0^2+\frac{3}{4}\alpha A^2)}(1-\frac{{\eta }^2}{16}),b_2(A)=\sqrt{({\omega }_0^2+\frac{3}{4}\alpha A^2)}(\frac{\eta }{2}+\frac{3{\eta }^3}{64}),\\ b_4(A)=\sqrt{({\omega }_0^2+\frac{3}{4}\alpha A^2)}(-\frac{{\eta }^2}{16}),b_6(A)=\sqrt{({\omega }_0^2+\frac{3}{4}\alpha A^2)}(\frac{3{\eta }^3}{64}).\end{array}$$

where

$$\eta =(\alpha a^2/4)/({\omega }_1^2+3\alpha a^2/4).$$

Equation (18) indicates the average frequency is $\omega (a)=b_0(a)$. Thus, applying the stochastic averaging method, we obtained drift and diffusion coefficients of Equation (9) defined as

$$\begin{array}{l}m(A)=\frac{A^2}{8g(A)}\left[\left[\epsilon c_1(4{\omega }_0^2+5\alpha A^2/2)-\epsilon c_2{A}^2({\omega }_0^2+3\alpha A^2/4)]\right]\right]\\ -\frac{A^2}{8g(A)}\left[8\pi (1-r)({\omega }_0^2+\alpha A^2/2)(b_0(A)-b_2(A)+b_4(A)-b_6(A)\right]\\ -\frac{A\mu }{\pi g(A)}\left[2b_0(A)-\frac{2}{3}{b}_2(A)-\frac{2}{15}{b}_4(A)-\frac{2}{35}{b}_6(A)\right]+\frac{\pi A}{8g(A)}\times \\ \{(4b_0{d}_0-2b_2{d}_0-2b_0{d}_2+b_2{d}_2)S(\omega (A))+(b_2{d}_2-b_4{d}_2-b_2{d}_4+b_4{d}_4)S(3\omega (A))\\ +(b_4{d}_4-b_6{d}_4-b_4{d}_6+b_6{d}_6)S(5\omega (A))+b_6(A)S_1(7\omega (A))d_6\}+\frac{\pi A}{8g^2(A)}\times \\ [(4b_0^2-b_2^2)S(\omega (A))+3(b_2^2-b_4^2)S(3\omega (A))+5(b_4^2-b_6^2)S(\omega 5(A))+7b_6^{2}S(7\omega (A))],\\ {\sigma }^2(A)=\frac{\pi A^2}{4g^2(A)}\{{(2b_0-b_2)}^{2}S(\omega (A))+{(b_2-b_4)}^{2}S(3\omega (A))\\ +{(b_4-b_6)}^{2}S(5\omega (A))+b_6^{2}S(7\omega (A)).\end{array}$$

where

$$d_i=\frac{d}{dA}\left[b_i(A)\frac{A}{g(A)}\right],(i=0,2,4,6).$$

5.1. The Effect of System Parameters on Response

In order to verify the accuracy of the proposed analytical method, the analytical stationary probability density was compared with the numerical results based on the original system (17). When system parameters were chosen as $c_2=0.05,c_1=0.025,\alpha =0.6,r=0.96,$ ${\omega }_0=1.0$, the effect of the magnitude parameter ${\alpha }_1$, friction amplitude $\mu$, and noise intensities $D$ on the responses of the system (17) were analyzed, respectively.

When $\mu =0.001,{\alpha }_1=2.0$ is fixed, all stationary PDFs of amplitude $A,$ displacement $x_1$, and velocity $x_2$ for different noise intensities $D$ are shown respectively in sub-captions (a–c) of Figure 2. These figures indicated that the influence of noise intensity on the response of system (17) is existent. Concretely, the increase of noise intensity $D$ leads to the reduction of all stationary PDFs of amplitude, displacement, and velocity.

When $\mu =0.001,D=0.01$ is fixed, all stationary PDFs of amplitude $A,$, displacement $x_1$, and velocity $x_2$ for different bandwidths ${\alpha }_1$ are shown respectively in sub-captions (a–c) of Figure 3. A larger bandwidth ${\alpha }_1$ corresponds to a higher peak of all stationary PDFs in Figure 3. When $D=0.01,{\alpha }_1=2.0$ is fixed, the stationary PDFs of amplitude $A,$ displacement $x_1$, and velocity $x_2$ for different friction parameters are shown respectively in sub-captions (a–c) of Figure 4. Concretely, the increase of friction amplitude $\mu$ leads to higher peaks of stationary PDFs.

By virtue of these figures the accuracy of the analytical method is corroborated by comparison with numerical results. As indicated, the bandwidth, friction amplitude, and noise intensities play a very constructive role in the response of the original system (17).

5.2. Stochastic Bifurcations

In this section, when the parameters $c_2=0.05,$ $c_1=0.025,\alpha =0.6,D=0.01,$ ${\omega }_0=1.0,{\alpha }_1=2.0$ are fixed, friction amplitude $\mu$ and restitution coefficient $r$ are adopted to analyze stochastic P-bifurcation.

Figure 5 shows the influence of the friction amplitude $\mu$ on the stochastic system when the restitution coefficient $r=0.975$ is fixed. The analytical results of joint PDFs on the displacement $x_1$ and velocity $x_2$ under different friction parameters $\mu$ are exhibited in sub-captions (a–c) of Figure 5. As shown in Figure 5a, there is a singular peak in the equilibrium position for the case of $\mu =0.008$. While $\mu$ is decreased to $\mu =0.001$, one observes the shape of the crater in Figure 5c. This phenomenon demonstrates that stochastic P-bifurcation occurred according to the concept of stochastic bifurcation [30]. Figure 5d exhibited the section graphs of analytical results of joint PDFs on the surface $x_1=0$ for different values $\mu$. We can conclude that joint PDFs exhibit a transformation from one peak to three peaks when the friction amplitude decreases slowly. Stochastic P-bifurcation occurs as the friction amplitude $\mu$ is close to $\mu \approx 0.038$.

As indicated, Figure 6a,b exhibited analytical results of joint PDFs of the displacement $x_1$ and velocity $x_2$ under different restitution coefficient $r$. The transformation of the joint PDFs from a singular peak to the shape of crater in the equilibrium position indicated stochastic P-bifurcation occurred as the restitution coefficient $r$ varied from $r=0.96$ to $r=0.975$.

Figure 6c exhibited the section graphs of analytical results of joint PDFs on the surface $x_1=0$ for different restitution coefficients $r$. In Figure 6c, the shape of joint PDFs changed from one peak to two peaks when the restitution coefficient $r$ increases gradually, which mean the appearance of stochastic P-bifurcation occurred. P-bifurcation occurs as the restitution coefficient $r$ is close to $r\approx 0.968$.

Based on the above analysis, it is obvious that the change of friction amplitude $\mu$ and restitution coefficient $r$ can arouse the occurrence of stochastic P-bifurcation.

6. Conclusions

This paper is devoted to researching the stationary response and bifurcation of a non-linear VI system with Coulomb friction excited by real noise. By virtue of the non-smooth transformation, the VI system is converted into a simplified system without barriers. Furthermore, the analytical solutions are derived by the stochastic averaging method. Effects of parameters such as excitation intensity, bandwidth, and friction amplitude are further analyzed in detail on the PDF distribution of the VI system. The analytical results obtained by proposed methods are verified by comparing with the simulated result of original systems. The stochastic P-bifurcation is explored by analyzing the stationary response. The joint PDFs indicate that friction amplitude $\mu$ and restitution coefficient $r$ can arouse the emergence of stochastic P-bifurcation.