# Chaotic Dynamics of a Mixed Rayleigh–Liénard Oscillator Driven by Parametric Periodic Damping and External Excitations

Authors: Yélomè Judicaël Fernando Kpomahou; Laurent Amoussou Hinvi; Joseph Adébiyi Adéchinan; Clément Hodévèwan Miwadinou
Journal: Complexity
 (2021)

---

## Abstract

In this paper, chaotic dynamics of a mixed Rayleigh–Liénard oscillator driven by parametric periodic damping and external excitations is investigated analytically and numerically. The equilibrium points and their stability evolutions are analytically analyzed, and the transitions of dynamical behaviors are explored in detail. Furthermore, from the Melnikov method, the analytical criterion for the appearance of the homoclinic chaos is derived. Analytical prediction is tested against numerical simulations based on the basin of attraction of initial conditions. As a result, it is found that forω=ν, the chaotic region decreases and disappears when the amplitude of the parametric periodic damping excitation increases. Moreover, increasing of F1 and F0 provokes an erosion of the basin of attraction and a modification of the geometrical shape of the chaotic attractors. For ω≠ν and η=0.8, the fractality of the basin of attraction increases as the amplitude of the external periodic excitation and constant term increase. Bifurcation structures of our system are performed through the fourth-order Runge–Kutta ode 45 algorithm. It is found that the system displays a remarkable route to chaos. It is also found that the system exhibits monostable and bistable oscillations as well as the phenomenon of coexistence of attractors.


---

## Body

## 1. Introduction

Many dynamical systems in various scientific fields are represented in terms of nonlinear ordinary differential equations which appear in a vast range of applications [1–7]. One of the classes of equations widely used is the class of mixed Liénard-type equations that has the following form:(1)x¨+fxx˙2+gxx˙+hx=0,where a dot represents the derivative with respect to time t. f,g, and h are arbitrary functions of x. It is important to notice that, from equation (1), when fx=0, the classical Liénard-type equation is obtained. This latter class of equations plays an important role in science and engineering [8–11], and its mathematical properties have been intensively studied in the literature from mathematical and physical points of view [12–15].Recently, the nonlinear dynamics of resonant RLC series circuit has been described using a general class of mixed Rayleigh–Liénard oscillators [16]. This new class of oscillators contains the class of oscillators given by equation (1). The class of mixed Rayleigh–Liénard oscillators is widely used in physics, electronics, mathematics, biology, engineering, and many other disciplines. In recent years, the nonlinear dynamics of this class of oscillators has been intensively studied, and interesting results such as period-doubling leading to chaotic motion, strange attractors, reverse period-doubling bifurcation, symmetry breaking, antimonotonicity, existence of horseshoe chaos, and so on have been obtained [17–25]. In most of these studies, the Melnikov perturbation method [18, 24, 26, 27] has been widely used to detect chaotic dynamics and to analyze near-homoclinic motion with deterministic or random perturbation. This method is today considered as a powerful analytical tool to provide an approximate criterion for the occurrence of hetero/homoclinic chaos in a wide class of dynamical systems. On the other hand, in these nonlinear dynamics works, the considered nonlinear oscillators are usually subjected to linear and nonlinear damping with constant coefficients. This damping combination is necessary since experiments showed a strong nonlinear dependence of damping with the vibration amplitude of nonlinear vibrating systems. To this end, several phenomenological nonlinear damping models have been introduced into existing fundamental nonlinear oscillators to understand the large amplitude, quasi-periodic and chaotic vibrations of nonlinear vibrating systems [24, 28–32]. In most cases, the phenomenological nonlinear damping does not reproduce with accuracy experimental results. Thus, the problem of finding nonlinear dissipative force capable to reproduce with accuracy experimental results has been recently studied. Amabili [33] used a fractional standard linear solid model and first-order harmonic balance method, and derived for the first time a nonlinear dissipative force of the type x2x˙ for the systems with only cubic stiffness. Recently again, Amabili [34] derived for the first time, from classical viscoelasticity and using the same method, a nonlinear dissipative force containing linear damping x˙, quadratic damping xx˙, and cubic damping x2x˙ terms. In addition, this author showed that in case of small damping, quadratic damping term disappears and the nonlinear damping is of the type x2x˙. In [35], Amabili showed that a mathematical model with linear damping and nonlinear damping of the type x2x˙ for a continuous system is capable to reproduce with accuracy experimental results for geometrically nonlinear vibrations. These results originally introduced by Amabili [34, 35] have been confirmed by Balasubramanian et al. [36]. Therefore, for a better understanding on the dynamical properties of the real systems, the dissipative force to choose must incorporate at least linear, quadratic, and cubic damping terms.Although the problem of identification of nonlinear damping is known, the effect of parametric damping on the dynamical behavior of nonlinear nonautonomous systems constitutes another interesting research topic. In this perspective, Mokni et al. [37] investigated the effect of cubic nonlinear parametric damping on vibration isolation of a sdof system, and the obtained results have revealed that the periodic nonlinear viscous damping enhances vibration isolation more than nonlinear time-dependent damping. On the other hand, the effect of linear parametric damping has been numerically studied by Rajasekar and Lakshmanan [38] in chaos control of an Bonhoeffer–Van der Pol oscillator. Parametric excitation with periodic damping has been examined for energy-harvesting applications [39, 40]. Kwuimy et al. [41] studied the effect of tilted harmonic excitation and parametric damping on the chaotic dynamics in an asymmetric magnetic pendulum. The obtained results have shown that the presence of parametric damping without a periodic fluctuation can enhance or suppress chaos, while a parametric damping with a periodic fluctuation can increase the region of regular motions significantly.It is then important to underline that the study of nonlinear oscillators under the influence of applied periodic force with the presence of the parametric periodic damping has attracted attention of very few investigators. Moreover, the effects of the parametric periodic damping and periodic excitation on the chaotic behavior of a mixed Rayleigh–Liénard oscillator with symmetric fourth-order polynomial double-well potential have not been discussed yet. Therefore, studying chaotic dynamics of such an oscillator must be performed for a better understanding of the dynamical characteristics of this system. In this perspective, we consider the following mixed Rayleigh-Liénard oscillator driven by parametric periodic damping and periodic external excitations:(2)x¨+α0+α1x+ηcosνtx˙+β0x˙2+β1x˙3−ω02x+γx3=ω02F0+F1cosωt,where ω0, γ, α0,1, β0,1, F0,1, η, ν, and ω are constant parameters. Physically, ω0 represents the natural frequency, γ denotes the cubic stiffness parameter, α0 means the linear damping coefficient, and α1 and β0,1 designate the nonlinear damping coefficients. F0 denotes the constant excitation amplitude. η and F1 represent the amplitude of the parametric damping and periodic external excitations, respectively. ν and ω represent, respectively, the frequency of the parametric and external periodic excitations. From equation (2), one can note that the nonlinear dissipative force takes into account at least linear, quadratic, and cubic damping terms required for a good description of the experimental results.Another motivation for our interest in this system is that it has a wide range of applications in various scientific fields. Whenη=α1=β0,1=F0=0, equation (2) serves as a basic model for self-excited oscillations in various disciplines. For α1=β0,1=F0,1=0, equation (2) arises in energy-harvesting applications [39, 40]. On the other hand, the system given by equation (2) belongs to a general class of mixed Rayleigh–Liénard oscillators recently proposed for describing the nonlinear dynamics of resonant RLC series circuit [16]. In nonlinear electronic circuits, note that the memristor represents the fourth basic circuit element. It is described by flux and charge [42, 43]. As a new device with memory, the memristor is used due to its potential applications, in many fields including nonvolatile memories on the nanoscale, neuromorphic systems and chaotic circuits [44, 45]. The study on memristive nonlinear electronic circuits is one of the active topics of research in the past few years, since memristors with versatile nonlinearities can be used in memristor-based oscillators, information encryption, memory, and so on. Given its potential applications, memristors are conveniently introduced into some existing linear or nonlinear electronic circuits to build various novel memristive chaotic circuits [46]. In this way, numerous mathematical models with different complexities have been proposed in the literature [47–52]. On the other hand, the two-dimensional FitzHugh–Nagumo (FHN) model [53, 54] is simplified from the four-dimensional Hodgkin–Huxley model for geometrical explanation of neuronal excitability and spiking under excitation. Several studies have been performed on the significant and complex dynamical aspects of the FHN model [51, 55–57]. Through all these studies, the results showed that the nonautonomous memristive chaotic circuits exhibit rich dynamical behaviors, such as chaos and hyperchaos, coexisting multiple attractors, hidden attractors, and hidden extreme multistability [47–52, 58–61].Therefore, dynamical behaviors of nonlinear electronic circuits must be investigated from theoretical and experimental points of view. Thus, the main objective of this paper is to investigate analytically and numerically the effects of the parametric periodic damping and external excitations on the transitions from regular to chaotic motions. In order to attain our objective, we investigate the equilibrium stability and its evolutions (Section2). We derive the condition of existence of horseshoes chaos using the Melnikov method, and we analyze the effect of certain parameters on chaotic behavior of the system (Section 3). Afterward, we analyze the bifurcation mechanisms and eventual routes to chaos in our system using the fourth-order Runge-Kutta ODE45 algorithm (Section 4). Finally, we end with a conclusion (Section 5).


## 2. Equilibrium Stability and Its Evolutions

The equilibrium points of system (2) can be obtained by solving x˙=y˙=0 [46]. So, these points are expressed as(3)E=x∗,0,where x∗ verifies the following equation:(4)x∗3−ω02γx∗−ω02γF0+F1cosωt=0.Equation (4) describes the dynamics of equilibrium state E of system (2) with time. The Jacobian matrix at the equilibrium state E and for ω=ν is derived from (2) as(5)JE=01ω02−3γx∗2−α0−α1x∗−ηcosωt.The characteristic polynomial of matrix (5) is(6)Ps=s2+mx∗s+nx∗=0,where(7)mx∗=ηγω02F1x∗3+α1−ηF1x∗+α0−ηF0F1,nx∗=3γx∗2−ω02.From equation (6), two eigenvalues can be expressed as(8)s1=−12mx∗−12lx∗,s2=−12mx∗+12lx∗,where lx∗=m2x∗−4nx∗.According to the sign ofmx∗, nx∗, and lx∗, the equilibrium stabilities can be separately analyzed by the eigenvalues s1,2. Thus, the following cases can be enumerated.Case 1.
mx∗=0 and nx∗≻0. In this case, the eigenvalues are both pure imaginary. Therefore, a Hopf bifurcation occurs when x∗≺ω0/3γ and x∗ satisfies the following equation:(9)x∗3+px∗+q=0,where p=ω02F1α1−η/ηγ and q=ω02F1α0−F0η/ηγ. The roots of equation (9) can be derived as [46, 62](10)x∗1=λ1×−q2+Δ3+λ2×−q2−Δ3,(11)x∗2=−q2+Δ3+−q2−Δ3,(12)x∗3=λ2×−q2+Δ3+λ1×−q2−Δ3,where Δ=q/22+p/33, λ1=−1+j3/2, and λ2=−1−j3/2.
According to Cardan discriminant [46, 62], when Δ≻0, one positive real equilibrium point is obtained from (11), since the equilibrium point cannot be a complex number. Thus, the system possesses in this case, one Hopf bifurcation point (HBP). However, two equilibrium points appear in equation (9) when Δ=0. Therefore, two Hopf bifurcation points appear in this case. Now, the system considered possesses three real equilibrium points from (10)–(12) when Δ≺0. Then, three Hopf bifurcation points appear in this condition.Case 2.
mx∗≻0 and nx∗=0. In this case, one obtains s1=−mx∗ and s2=0. Then, the system presents two fold bifurcation points (FBPs) at E=±ω0/3γ,0.Case 3.
mx∗≻0 and nx∗≺0. There are two real roots with opposite signs, implying that the equilibrium state E is an unstable saddle point.Case 4.
mx∗≺0, nx∗≻0, and lx∗≻0. In this case, there are two positive real roots, indicating that the equilibrium state E is an unstable node point.Case 5.
mx∗≺0 and lx∗≺0. There is a pair of complex conjugate roots with positive real parts, which means that the equilibrium state E is an unstable node focus.Case 6.
mx∗≻0, nx∗≻0, and lx∗≺0. One obtains, in this case, a pair of complex conjugate roots with negative real parts, indicating that the equilibrium state E is a stable node focus.Case 7.
mx∗≻0, nx∗≻0, and lx∗≻0. Equation (6) has two negative real roots, which implies that the equilibrium state E is a stable node point.
Forω0=F0=0.25, α0=0.015, γ=1, η=0.8, and ω=ν=0.6, we get numerically the values of the equilibrium state E as well as the Hopf and fold bifurcation points with the time evolution, as shown in Figure 1.
The above analytical and numerical results show that with the evolution of time, the equilibrium point is variable, and its stability has complex transitions between stable points and unstable points via three Hopf bifurcations and two fold bifurcations and then initiates the next period, leading to the generation of complex dynamical behaviors.Figure 1
Equilibrium points vary with the evolution of the time under different values ofF1 and α1: (a) F1=2.2 and α1=−0.25; (b) F1=0.5 and α1=0.025.
(a)(b)


## 3. Melnikov Analysis

### 3.1. Fixed Points and Melnikov Criterion for Homoclinic Chaos

In this section, we first derive the fixed points and the phase portrait corresponding to the system of equation (2) when it is unperturbed. Thus, if we let α0,1=β0,1=F0,1=η=0, the system of equation (2) becomes an unperturbed system and can be rewritten as follows:(13)x˙=y,y˙=ω02x−γx3.The system of equations (13) corresponds to an integrable Hamiltonian system with a potential function(14)Vx=−12ω02x2+14γx4,and the associated Hamiltonian function is given by(15)Hx,y=12y2−12ω02x2+14γx4.From equation (13), we can determine the fixed points and analyze their stabilities. Thus, if γ≺0, the system of equations (13) admits one fixed point 0,0 which is a saddle. It has only a single-well potential. If γ≻0, the system of equation (12) has three fixed points 0,0, x1=−ω0/γ,0, and x2=ω0/γ,0: There are one saddle and two centers. In this case, it has a double-well symmetric potential. We here are interested in double-well potential. Therefore, we fix γ=1 throughout this paper. The potential function and the phase portrait of the system of equations (13) are plotted in Figures 2(a) and 2(b), respectively.Figure 2
(a) Potential function given by equation (14) and (b) corresponding homoclinic orbits for three different values of ω0.
(a)(b)Second, we now suppose that the unperturbed system discussed above is perturbed by a combination of dissipative and external excitations forces. In this section, we investigate analytically the condition for the onset of horseshoes chaos by applying the Melnikov method to equation (2). To this end, scaling α0,1=εα0,1, β0,1=εβ0,1, F0,1=εF0,1, and η=εη, the mixed Rayleigh–Liénard oscillator equation (2) can be rewritten as a first-order differential system in the form(16)x˙=y,y˙=ω02x−γx3+ε−α0−α1x−ηcosνty+ε−β0y2−β1y3+ω02F0+F1cosωt,where ε is a small parameter. It is shown that when ε=0, the system of equations (16) is a Hamiltonian system with a pair of homoclinic trajectories connecting the saddle to itself given by(17)xh±yh±=±ω02γsechω0τ±ω022γsechω0τtanhω0τ,with τ=t−t0. t0 is the cross-section time of the Poincaré map and can be interpreted as the initial time of the forcing time. According to equation (17), the signs refer to the right and left half planes. Both solutions determine the separatrix orbit, since it separates two types of orbits in phase space. ''+'' refers to the homoclinic trajectory with x≻0, and ''-'' refers to the homoclinic trajectory with x≺0. Turning our attention to the perturbed system, that is, ε≠0, we can compute the Melnikov integrals in order to find if the two-dimensional stable and unstable manifolds of the periodic orbits intersect. In this case, the Melnikov function Mt0 is a scalar, and we have(18)Mt0=−α0∫−∞+∞yh2dτ−α1∫−∞+∞xhyh2dτ−β0∫−∞+∞yh3dτ−β1∫−∞+∞yh4dτ+ω02F0∫−∞+∞yhdτ−η∫−∞+∞yh2cosντ+νt0dτ+ω02F1∫−∞+∞yhcosωτ+ωt0dτ,where xh±,yh± = xht,yht.Substituting equation (17) into (18) and evaluating the integral by using the standard integral tables [63], we obtain the Melnikov function:(19)M±t0=A±∓BF1sinωt0−ηCcosνt0,with(20)A±=−4ω02α03γ−16β1ω0735γ2±2ω03F02γ−πω04α182γ3/2−4β0ω05152γ3/2,B=ω02πω2/γcoshπω/2ω01+coshπω/ω0,C=πν2ω02−ν23γsinhπω/ω0.IfMt0 has simple zeros, then a homoclinic bifurcation occurs. Thus, in the case where the amplitude of the parametric periodic damping is equal to zero η=0, the transverse crossings of stable and unstable manifolds exist if Melnikov function (19) has simple zeros, that is to say,(21)A±BF1≺1.Therefore, the necessary condition for which the invariant manifolds intersect is given by(22)F1≻A±B.Ifη≠0, a necessary condition for the onset of Melnikov chaos in our system can be expressed from [64–67] as follows:(23)η≺1∓BF1A±R,with R=A±/CIn order to have a visual information on the range where the horseshoe chaos is possible, we have plotted in Figures3(a)–3(d) the dependence of the amplitude η of the parametric periodic damping excitation on the frequency ω=ν for different values of F1, β0, β1, and F0. From Figure 3(a), we notice that as the amplitude F1 of the external periodic excitation increases, the region where the horseshoes chaos appears decreases. The same observation can be made when F0 decreases (see Figure 3(d)). Furthermore, when the nonlinear damping coefficients β0,1 increase, horseshoes chaos domain decreases (see Figures 3(b) and 3(c)).Figure 3
Melnikov threshold curvesΓ+ for horseshoe chaos in the η,ω plane. The basic values are ω0=F0=0.25, α0=0.015α1=0.025, γ=1, F1=0.5, β0=0.01, and β1=0.005.
(a)(b)(c)(d)


### 3.2. Fractal Basin Boundaries

The aim of this section is to verify the validity of the analytical predictions obtained in the previous section and to study the effect ofη, F0, and F1 on basin of attraction of the mixed Rayleigh–Liénard oscillator under study. For this purpose, we numerically solve the equation of motion (2) by using the fourth-order Runge–Kutta ODE45 algorithm with time step 10-3 for generating basins of attraction in a plane element −1.5×1.5,−2×2. For ω0=F0=0.25, α0=0.015, α1=0.025, γ=1, F1=0.5, β0=0.01, β1=0.005, and ω=ν=0.6, the Melnikov threshold for the right plane is given analytically by 1.1397. Therefore, below this critical value, the system under consideration may display a chaotic behavior. However, regular behavior may appear above this critical value. Now, for η=2 chosen in the regular behavior domain predicted by the Melnikov criterion, the basin of attraction shows regular solutions (see Figure 4(d)). When the amplitude η of the parametric periodic damping excitation is chosen in the chaotic behavior region, a fractal structure of the basin of attraction of the initial conditions appears and becomes more and more visible as η decreases (see Figures 4(a)–4(c)). We can conclude that the analytical and numerical predictions are in good agreement.Figure 4
Basins of attraction for motion aroundx=x2 in mixed Rayleigh–Liénard system with double-well potential: (a) η=0.4, (b) η=0.6, (c) η=1, and (d) η=2.
(a)(b)(c)(d)The effects ofF1 and F0 on basin of attraction are investigated through Figures 5 and 6 for η=0.8. From Figure 5, we notice that the basin of attraction becomes less eroding as the amplitude of the external periodic excitation decreases. However, the basin of attraction is destroyed and the fractal behavior becomes more and more visible when the constant term F0 increases (see Figure 6). It is important to point out that the geometric shape of chaotic attractors has completely changed as F1 and F0 increase. When ω≠ν and η=0.8, the fractal behavior increases with increasing F0 and F1 (see Figure 7). Moreover, the geometric shape of chaotic attractors also changes when F1 and F0 increase.Figure 5
Effect ofF1 on basin of attraction with the parameters of Figure 3 and η=0.8: (a) F1=0, (b) F1=1, (c) F1=5, and (d) F1=10.
(a)(b)(c)(d)Figure 6
Effect ofF0 on basin of attraction with the parameters of Figure 3 and η=0.8: (a) F0=0, (b) F0=1, (c) F0=3, and (d) F0=5.
(a)(b)(c)(d)Figure 7
Effects ofF1 and F0 on basin of attraction with ω=1, ν=5−1/2, and η=0.8: (a) F1=0.5, F0=0.25; (b) F0=0.25, F1=1.1; (c) F0=0.25, F1=5; (d) F1=0.5, F0=2; and (e) F1=0.5, F0=5.
(a)(b)(c)(d)(e)


## 3.1. Fixed Points and Melnikov Criterion for Homoclinic Chaos

In this section, we first derive the fixed points and the phase portrait corresponding to the system of equation (2) when it is unperturbed. Thus, if we let α0,1=β0,1=F0,1=η=0, the system of equation (2) becomes an unperturbed system and can be rewritten as follows:(13)x˙=y,y˙=ω02x−γx3.The system of equations (13) corresponds to an integrable Hamiltonian system with a potential function(14)Vx=−12ω02x2+14γx4,and the associated Hamiltonian function is given by(15)Hx,y=12y2−12ω02x2+14γx4.From equation (13), we can determine the fixed points and analyze their stabilities. Thus, if γ≺0, the system of equations (13) admits one fixed point 0,0 which is a saddle. It has only a single-well potential. If γ≻0, the system of equation (12) has three fixed points 0,0, x1=−ω0/γ,0, and x2=ω0/γ,0: There are one saddle and two centers. In this case, it has a double-well symmetric potential. We here are interested in double-well potential. Therefore, we fix γ=1 throughout this paper. The potential function and the phase portrait of the system of equations (13) are plotted in Figures 2(a) and 2(b), respectively.Figure 2
(a) Potential function given by equation (14) and (b) corresponding homoclinic orbits for three different values of ω0.
(a)(b)Second, we now suppose that the unperturbed system discussed above is perturbed by a combination of dissipative and external excitations forces. In this section, we investigate analytically the condition for the onset of horseshoes chaos by applying the Melnikov method to equation (2). To this end, scaling α0,1=εα0,1, β0,1=εβ0,1, F0,1=εF0,1, and η=εη, the mixed Rayleigh–Liénard oscillator equation (2) can be rewritten as a first-order differential system in the form(16)x˙=y,y˙=ω02x−γx3+ε−α0−α1x−ηcosνty+ε−β0y2−β1y3+ω02F0+F1cosωt,where ε is a small parameter. It is shown that when ε=0, the system of equations (16) is a Hamiltonian system with a pair of homoclinic trajectories connecting the saddle to itself given by(17)xh±yh±=±ω02γsechω0τ±ω022γsechω0τtanhω0τ,with τ=t−t0. t0 is the cross-section time of the Poincaré map and can be interpreted as the initial time of the forcing time. According to equation (17), the signs refer to the right and left half planes. Both solutions determine the separatrix orbit, since it separates two types of orbits in phase space. ''+'' refers to the homoclinic trajectory with x≻0, and ''-'' refers to the homoclinic trajectory with x≺0. Turning our attention to the perturbed system, that is, ε≠0, we can compute the Melnikov integrals in order to find if the two-dimensional stable and unstable manifolds of the periodic orbits intersect. In this case, the Melnikov function Mt0 is a scalar, and we have(18)Mt0=−α0∫−∞+∞yh2dτ−α1∫−∞+∞xhyh2dτ−β0∫−∞+∞yh3dτ−β1∫−∞+∞yh4dτ+ω02F0∫−∞+∞yhdτ−η∫−∞+∞yh2cosντ+νt0dτ+ω02F1∫−∞+∞yhcosωτ+ωt0dτ,where xh±,yh± = xht,yht.Substituting equation (17) into (18) and evaluating the integral by using the standard integral tables [63], we obtain the Melnikov function:(19)M±t0=A±∓BF1sinωt0−ηCcosνt0,with(20)A±=−4ω02α03γ−16β1ω0735γ2±2ω03F02γ−πω04α182γ3/2−4β0ω05152γ3/2,B=ω02πω2/γcoshπω/2ω01+coshπω/ω0,C=πν2ω02−ν23γsinhπω/ω0.IfMt0 has simple zeros, then a homoclinic bifurcation occurs. Thus, in the case where the amplitude of the parametric periodic damping is equal to zero η=0, the transverse crossings of stable and unstable manifolds exist if Melnikov function (19) has simple zeros, that is to say,(21)A±BF1≺1.Therefore, the necessary condition for which the invariant manifolds intersect is given by(22)F1≻A±B.Ifη≠0, a necessary condition for the onset of Melnikov chaos in our system can be expressed from [64–67] as follows:(23)η≺1∓BF1A±R,with R=A±/CIn order to have a visual information on the range where the horseshoe chaos is possible, we have plotted in Figures3(a)–3(d) the dependence of the amplitude η of the parametric periodic damping excitation on the frequency ω=ν for different values of F1, β0, β1, and F0. From Figure 3(a), we notice that as the amplitude F1 of the external periodic excitation increases, the region where the horseshoes chaos appears decreases. The same observation can be made when F0 decreases (see Figure 3(d)). Furthermore, when the nonlinear damping coefficients β0,1 increase, horseshoes chaos domain decreases (see Figures 3(b) and 3(c)).Figure 3
Melnikov threshold curvesΓ+ for horseshoe chaos in the η,ω plane. The basic values are ω0=F0=0.25, α0=0.015α1=0.025, γ=1, F1=0.5, β0=0.01, and β1=0.005.
(a)(b)(c)(d)


## 3.2. Fractal Basin Boundaries

The aim of this section is to verify the validity of the analytical predictions obtained in the previous section and to study the effect ofη, F0, and F1 on basin of attraction of the mixed Rayleigh–Liénard oscillator under study. For this purpose, we numerically solve the equation of motion (2) by using the fourth-order Runge–Kutta ODE45 algorithm with time step 10-3 for generating basins of attraction in a plane element −1.5×1.5,−2×2. For ω0=F0=0.25, α0=0.015, α1=0.025, γ=1, F1=0.5, β0=0.01, β1=0.005, and ω=ν=0.6, the Melnikov threshold for the right plane is given analytically by 1.1397. Therefore, below this critical value, the system under consideration may display a chaotic behavior. However, regular behavior may appear above this critical value. Now, for η=2 chosen in the regular behavior domain predicted by the Melnikov criterion, the basin of attraction shows regular solutions (see Figure 4(d)). When the amplitude η of the parametric periodic damping excitation is chosen in the chaotic behavior region, a fractal structure of the basin of attraction of the initial conditions appears and becomes more and more visible as η decreases (see Figures 4(a)–4(c)). We can conclude that the analytical and numerical predictions are in good agreement.Figure 4
Basins of attraction for motion aroundx=x2 in mixed Rayleigh–Liénard system with double-well potential: (a) η=0.4, (b) η=0.6, (c) η=1, and (d) η=2.
(a)(b)(c)(d)The effects ofF1 and F0 on basin of attraction are investigated through Figures 5 and 6 for η=0.8. From Figure 5, we notice that the basin of attraction becomes less eroding as the amplitude of the external periodic excitation decreases. However, the basin of attraction is destroyed and the fractal behavior becomes more and more visible when the constant term F0 increases (see Figure 6). It is important to point out that the geometric shape of chaotic attractors has completely changed as F1 and F0 increase. When ω≠ν and η=0.8, the fractal behavior increases with increasing F0 and F1 (see Figure 7). Moreover, the geometric shape of chaotic attractors also changes when F1 and F0 increase.Figure 5
Effect ofF1 on basin of attraction with the parameters of Figure 3 and η=0.8: (a) F1=0, (b) F1=1, (c) F1=5, and (d) F1=10.
(a)(b)(c)(d)Figure 6
Effect ofF0 on basin of attraction with the parameters of Figure 3 and η=0.8: (a) F0=0, (b) F0=1, (c) F0=3, and (d) F0=5.
(a)(b)(c)(d)Figure 7
Effects ofF1 and F0 on basin of attraction with ω=1, ν=5−1/2, and η=0.8: (a) F1=0.5, F0=0.25; (b) F0=0.25, F1=1.1; (c) F0=0.25, F1=5; (d) F1=0.5, F0=2; and (e) F1=0.5, F0=5.
(a)(b)(c)(d)(e)


## 4. Bifurcation Mechanisms and Transition to Chaos

The aim of this section is to investigate some bifurcation mechanisms in our system when certain system parameters evolve. For this purpose, we numerically solve the equation of motion (2) under consideration using the fourth-order Runge–Kutta ODE45 algorithm with time step 10-2. The initial conditions used to perform the numerical simulations are X0=0.5,0.5 and Y0=−0.5,−0.5. The bifurcation diagram and its Lyapunov exponent as function of the amplitude η are plotted in Figure 8. Through this figure, various bifurcations such as period doubling, quasiperiodic, periodic windows, intermittency, and chaos are exhibited by system (2). By comparing the two sets of data, one can notice from Figure 8 that system (2) displays the monostability and bistability phenomena as well as the phenomenon of coexistence of attractors. These phenomena are illustrated in Figure 9 for several different values of η.Figure 8
Bifurcation diagrams and Lyapunov exponents versus the amplitudeη of the parametric periodic damping excitation with the parameters of Figure 1(a) for ω=ν=0.2 and F1=0.1. The specific initial conditions used are X0=0.5,0.5 (blue color) and Y0=−0.5,−0.5 (red color). Bifurcation diagrams are obtained by scanning the parameter η upwards (blue color) and downwards (red color).Figure 9
Various phase portraits of system (2) illustrating the monostability and bistability oscillations and the coexistence of attractors for different values of η with the parameters of Figure 8.When the amplitudeF1 is used as control parameter, the dynamics of system (2) becomes rich. One can see that the system presents periodic windows, reverse periodic windows, periodic bubble, intermittency, and chaotic oscillations (see Figure 10). In addition, we also notice that system (2) exhibits the coexistence of attractors as well as the monostability and bistability phenomena. In order to have an idea about these attractors, various phase portraits of system (2) are plotted in Figure 11 for different values of F1. We have also investigated through Figure 12 the effect of the constant term F0 on the bifurcation diagram when F1=0.1 and ω=ν=0.2. We clearly see that system (2) exhibits periodic windows, reverse periodic windows, intermittency, and chaos. We also notice that the chaotic behavior is more dominant than periodic behavior when F0 evolves. On the other hand, the system under consideration exhibits monostability and bistability phenomena as well as the phenomenon of coexistence of attractors.Figure 10
Bifurcation diagrams and Lyapunov exponents versus the amplitudeF1 with the parameters of Figure 1(a) for ω=0.2 and specific initial conditions: X0=0.5,0.5 (blue color) and Y0=−0.5,−0.5 (red color). Bifurcation diagrams are obtained by scanning the parameter F1 upwards (blue color) and downwards (red color).Figure 11
Various phase portraits of system (2) showing the monostable and bistable oscillations and the coexistence of attractors for different values of F1 with the parameters of Figure 10.Figure 12
Bifurcation diagrams and Lyapunov exponents versus the constant termF0 with the parameters of Figure 1(a) for ω=ν=0.2 and F1=0.1. The specific initial conditions used are X0=0.5,0.5 (blue color) and Y0=−0.5,−0.5 (red color). Bifurcation diagrams are obtained by scanning the parameter F0 upwards (blue color) and downwards (red color).Whenω=0.2, ν=5−1/2, and F1=0.1, the bifurcation diagram and its corresponding Lyapunov exponents as F0 evolves are shown in Figure 13. We observe that system (2) transits from a chaotic motion to a period-3 motion followed by a chaotic motion. Afterward, the chaotic behavior disappears and a period-2 orbit takes place in the system. We also notice that the periodic behavior is more dominant than chaotic behavior. By comparing the two sets of data, one can notice from Figure 13 that system (2) shows the monostability and bistability phenomena as well as the phenomenon of coexistence of attractors.Figure 13
Bifurcation diagrams and Lyapunov exponents as a function of the constant termF0 with the parameters of Figure 1(a) for ω=0.2, ν=5−1/2, and F1=0.1 and specific initial conditions X0=0.5,0.5 (blue color) and Y0=−0.5,−0.5 (red color). Bifurcation diagrams are obtained by scanning the parameter F0 upwards (blue color) and downwards (red color).


## 5. Conclusions

This paper deals with chaotic dynamics in mixed Rayleigh–Liénard oscillator driven by parametric periodic damping and external excitations. The equilibrium point and its evolutions with the time are theoretically investigated. It is found that the equilibrium stability has complex dynamical transitions between stable and unstable points via Hopf/fold bifurcations. Moreover, The Melnikov perturbation method is used to derive the analytical criterion for the appearance of chaos in the system under consideration. A convenient demonstration of the use and accuracy of the method is obtained from the basin of attraction. The effects of the parametric periodic damping and external excitations are investigated. It is found that forω=ν, the horseshoe chaos decreases and disappears as the amplitude of the parametric periodic damping increases. For η=0.8, the basin of attraction is destroyed and the fractal behavior becomes more and more visible as the constant term F0 increases. However, the fractality of the basin of attraction decreases when the amplitude F1 of the external periodic excitation increases. Moreover, the geometric shape of chaotic attractors is modified when F0 and F1 increase. For ω≠ν and η=0.8, the fractality of the basin of attraction becomes visible when the amplitude of the external periodic excitation and the constant term increase. Bifurcation diagrams, Lyapunov exponents, and phase portraits are plotted using the fourth-order Runge–Kutta ode 45 algorithm, and period-doubling, periodic windows, reverse periodic windows, periodic bubble, intermittency, and chaotic behavior are obtained. It is important to underline that the system under consideration can vibrate from a chaotic motion to a period-3 motion followed by a chaotic motion. In addition, we also notice that the system presents monostable and bistable oscillations as well as the coexistence of attractors.



---

*Source: 6631094-2021-03-23.xml*