# Dynamic Behaviors of a Two-Degree-of-Freedom Impact Oscillator with Two-Sided Constraints

Authors: Songtao Li; Qunhong Li; Zhongchuan Meng
Journal: Shock and Vibration
 (2021)

---

## Abstract

The dynamic model of a vibroimpact system subjected to harmonic excitation with symmetric elastic constraints is investigated with analytical and numerical methods. The codimension-one bifurcation diagrams with respect to frequency of the excitation are obtained by means of the continuation technique, and the different types of bifurcations are detected, such as grazing bifurcation, saddle-node bifurcation, and period-doubling bifurcation, which predicts the complexity of the system considered. Based on the grazing phenomenon obtained, the zero-time-discontinuity mapping is extended from the single constraint system presented in the literature to the two-sided elastic constraint system discussed in this paper. The Poincare mapping of double grazing periodic motion is derived, and this compound mapping is applied to obtain the existence conditions of codimension-two grazing bifurcation point of the system. According to the deduced theoretical result, the grazing curve and the codimension-two grazing bifurcation points are validated by numerical simulation. Finally, various types of periodic-impact motions near the codimension-two grazing bifurcation point are illustrated through the unfolding diagram and phase diagrams.


---

## Body

## 1. Introduction

In mechanical engineering, there exist the vibroimpact phenomena widely, and systems interacting via impact have been extensively studied in recent years. Over the past years, scholars have mainly dedicated themselves to study the bifurcation phenomenon in smooth systems. Recently, a lot of work has gone into investigating nonsmooth bifurcations [1–8] of dynamical systems. The focus of investigations has gradually begun to change from a unilateral constraint system [9–11] to a multiconstraint system [12–17]. The impact oscillators can be divided into rigid or elastic impact oscillators according to the hardness of constraint. Shaw and Holmes [18] studied the motions of a system with a single piece rigid stop by using a one-dimensional mapping. Jiang and Wiercigroch [19] developed the concept of discontinuity geometry of rigid impact oscillators into the elastic impact oscillators, and the geometry analysis methods are applied to study the mechanisms of grazing bifurcations of system with unilateral soft constraint. Ing et al. [20] applied the methods of theoretical analyses and simulation to study different bifurcation scenarios for an impact oscillator with one-sided amplitude constraint. Du et al. [21] studied the intermittent chaos control method of a symmetrical collision system with the two-degree-of-freedom elastic double-impact system. The controlling idea of the Hopf bifurcation was applied to the system to provide a new control method for the chaos control of such system. Gritli et al. [22, 23] considered a state-feedback controller and applied the Linear Matrix Inequality (LMI) approach to control the motion of the system subject to norm-bounded parametric uncertainties. Shen et al. [24, 25] proposed a discrete feedback controller to suppress grazing-induced motion. Gritli and Belghith [26–29] considered the dynamic behavior of the model under the OGY-based state-feedback control, and the bifurcation phenomena were carried out via the bifurcation diagrams. Chávez et al. [30] studied a single-parameter and two-parameter bifurcations of the system that consists of two identical Duffing oscillators interacting via soft impact with the aid of COCO.As a special bifurcation of nonsmooth systems, the grazing bifurcation is a qualitative transition between impact and nonimpact motion [31–46]. Nordmark [31, 32] studied the dynamic behavior near the grazing trajectory and a nonlinear mapping is obtained. The normal form mapping for grazing bifurcation in an n-dimensional impact system is obtained in [35]. The existence condition and stability of the grazing periodic orbit are derived in a two degree-of-freedom vibroimpact system with one-sided constraint in [36]. Yin et al. [37–42] analyzed and measured the degenerate grazing bifurcations of the impact systems. The dynamic behavior of the grazing bifurcation is investigated by using the Poincare composite mapping technology, and the conditions of the degenerate grazing bifurcation are found. And the relationship between the observed bifurcations and degenerate grazing point is presented by a two-parameter continuation. Kowalczyk et al. [43] made a simple classification of codimension-two bifurcations of nonsmooth systems. Dankowicz and Zhao [44] deduced specific formulae for the local map on the vicinity of the codimension-two points and presented the bifurcation behavior by numerical simulations. The codimension-two grazing bifurcations in single-degree-of-freedom impact oscillators are studied and the dynamic response near the bifurcation points is presented in [45]. Xu et al. [46] investigated the codimension-two grazing bifurcation in n-degree-of-freedom impact oscillator with bilateral constraints by using a classical approach of discontinuity mappings.Based on previous studies, we discuss the bifurcation behavior of the system with a two-sided elastic constraint by using the path-following method. And then, abundant and complex bifurcation behaviors are exhibited, such as grazing bifurcation, saddle-node bifurcation, fold bifurcation, period-doubling bifurcation. Much attention has been paid to analyze codimension-two bifurcation of a system with rigid impact, little work has been studied on the analysis of such bifurcation in the soft impact system. And then, on the basis of discontinuity mappings of the rigid vibroimpact system, the zero-time-discontinuity mappingsZDM1, ZDM2 of the elastic constraint vibroimpact system are deduced. And the compound mapping P is applied to obtain the existence conditions and specific mathematical expressions of codimension-two grazing bifurcation point.This paper is discussed as follows: in section2 we introduce a dynamic model of a two-degree-of-freedom vibroimpact system with two-sided soft constraints. The codimension-1 bifurcation analysis of the system is presented in section 3. In section 4, we derive a Poincare mapping near the grazing bifurcation by using the zero-discontinuity-mapping method; the numerical simulations according to the deduced theoretical result are presented. In Section 5, the conclusions are given.


## 2. Physical Model and Equations of Motion

### 2.1. Physical Model

A physical model for the two-degree-of-freedom impact oscillator with massesM1 and M2 is described in Figure 1. The masses M1 and M2 are connected via the linear springs and linear viscous dampers, and the mass M1 is limited by the symmetrical elastic constraint corresponding to two discontinuity surfaces D1 and D2. The excitations on the masses M1 and M2 are harmonic forces with amplitudes P1,P2, respectively. The excitation frequency and the phase are the same for the masses, where the frequency is taken as the controlling parameter in the following: when the displacement X1 of the oscillator is Dor−D, it will collide with the right (or left) elastic constraint. Damping in this model is considered as proportional damping of the Rayleigh type, which implies C1/K1=C2/K2.Figure 1
The model of the two-degree-of-freedom vibroimpact system with symmetrical constraints.


### 2.2. Differential Equation of Motion

The equation of motion of the physical model in Figure1 can be expressed as follows:(1)M1X¨1+K1X1−X2+C1X˙1−X˙2+Fc=P1sinΩT+τ,(2)M2X¨2+K2X2+C2X˙2+K1X2−X1+C1X˙2−X˙1=P2sinΩT+τ,Fc=K0X1+D,X1≤−D,0,−D<X1<D,K0X1−D,X1≥D.Introduce the following dimensionless variables and time:(3)μm=M2M1,μk1=K2K1,μc=C2C1,μk0=K0K1,f2=P2P0,t=TK1M1,ω=ΩK1M1,ζ=C12K1M1,d=DK1P0,xi=XiK1P0,i=1,2.Systems (1) and (2) can be expressed in a normalized form as(4)Mx¨+2Zx˙+Kx+fc=F,in which(5)M=100μm,K=1−1−11+μk1,Z=ζK=ζ−ζ−ζζ×1+μc,F=1−f2f2,fc=μk0x1+d,0T,x1≤−d,0,0T,−d<x1<d,μk0x1−d,0T,x1≥d,where the dot ‘⋅’ denotes the differentiation with respect to the nondimensional time t. Let x=x1,x2T; it satisfies that x˙=v1,v2T and M and K denote the matrix of mass and the matrix of stiffness in nondimension form, respectively. The state-space discontinuity surfaces D1 and D2 can be expressed as D1=x|hLCx=x1+d=0, D2=x|hRCx=x1−d=0, where hLCx denote the distance to the left constraint D1 and hRCx denote the distance to the right constraint D2.


### 2.3. Grazing Periodic Motion

The general solution for nonimpact motions of equation (4) is described in the following form:(6)xit=∑j=12φije−ηjtajcosωdjt+Ajsinωt+τ+Bjcosωt+τ,xi.t=∑j=12φije−ηjtbjωdj−ajηjcosωdjt−bjηj+ajωdjsinωdjt+Ajωcosωt+τ−Bjωsinωt+τ,where φij are the elements of the canonical modal matrix Ψ, ηj=ξωj2, ωdj=ωj2−ηj2, f˜1,f˜2T=ΨT1−ppT, aj and bj are the constants of integration that are determined by the initial conditions and system parameter, Aj and Bj are the constants of amplitude, and their expressions are as follows:(7)Aj=f˜jωj2−ω2ωj2−ω22+2ηjω2,Bj=f˜j−2ηjωωj2−ω22+2ηjω2.The considered system consists of bilateral symmetrical elastic stops; the initial conditions and periodic motion conditions of the orbit with double grazing bifurcation satisfy the following formula:(8)x10=x10=d,x20=x20,v10=v10=0,v20=v20,x1πω=−d,x12pπω=x10=d,x22pπω=x20,v12pπω=v10=0,v22pπω=v20.By substituting the above conditions into the general solution, we can see that(9)d2=d12+d22≠0,τ=arctand1d2,where d1=A1φ11+A2φ12,d2=B1φ11+B2φ12.With appropriate system parameters, the motion of the system is periodic. The periodic motion of the system is represented byp,q,n, the number of periods is represented by n, p denotes the number of collisions between the mass M1 and the left side constraint, and q denotes the number of collisions between the mass M1 and the right side constraint.


## 2.1. Physical Model

A physical model for the two-degree-of-freedom impact oscillator with massesM1 and M2 is described in Figure 1. The masses M1 and M2 are connected via the linear springs and linear viscous dampers, and the mass M1 is limited by the symmetrical elastic constraint corresponding to two discontinuity surfaces D1 and D2. The excitations on the masses M1 and M2 are harmonic forces with amplitudes P1,P2, respectively. The excitation frequency and the phase are the same for the masses, where the frequency is taken as the controlling parameter in the following: when the displacement X1 of the oscillator is Dor−D, it will collide with the right (or left) elastic constraint. Damping in this model is considered as proportional damping of the Rayleigh type, which implies C1/K1=C2/K2.Figure 1
The model of the two-degree-of-freedom vibroimpact system with symmetrical constraints.


## 2.2. Differential Equation of Motion

The equation of motion of the physical model in Figure1 can be expressed as follows:(1)M1X¨1+K1X1−X2+C1X˙1−X˙2+Fc=P1sinΩT+τ,(2)M2X¨2+K2X2+C2X˙2+K1X2−X1+C1X˙2−X˙1=P2sinΩT+τ,Fc=K0X1+D,X1≤−D,0,−D<X1<D,K0X1−D,X1≥D.Introduce the following dimensionless variables and time:(3)μm=M2M1,μk1=K2K1,μc=C2C1,μk0=K0K1,f2=P2P0,t=TK1M1,ω=ΩK1M1,ζ=C12K1M1,d=DK1P0,xi=XiK1P0,i=1,2.Systems (1) and (2) can be expressed in a normalized form as(4)Mx¨+2Zx˙+Kx+fc=F,in which(5)M=100μm,K=1−1−11+μk1,Z=ζK=ζ−ζ−ζζ×1+μc,F=1−f2f2,fc=μk0x1+d,0T,x1≤−d,0,0T,−d<x1<d,μk0x1−d,0T,x1≥d,where the dot ‘⋅’ denotes the differentiation with respect to the nondimensional time t. Let x=x1,x2T; it satisfies that x˙=v1,v2T and M and K denote the matrix of mass and the matrix of stiffness in nondimension form, respectively. The state-space discontinuity surfaces D1 and D2 can be expressed as D1=x|hLCx=x1+d=0, D2=x|hRCx=x1−d=0, where hLCx denote the distance to the left constraint D1 and hRCx denote the distance to the right constraint D2.


## 2.3. Grazing Periodic Motion

The general solution for nonimpact motions of equation (4) is described in the following form:(6)xit=∑j=12φije−ηjtajcosωdjt+Ajsinωt+τ+Bjcosωt+τ,xi.t=∑j=12φije−ηjtbjωdj−ajηjcosωdjt−bjηj+ajωdjsinωdjt+Ajωcosωt+τ−Bjωsinωt+τ,where φij are the elements of the canonical modal matrix Ψ, ηj=ξωj2, ωdj=ωj2−ηj2, f˜1,f˜2T=ΨT1−ppT, aj and bj are the constants of integration that are determined by the initial conditions and system parameter, Aj and Bj are the constants of amplitude, and their expressions are as follows:(7)Aj=f˜jωj2−ω2ωj2−ω22+2ηjω2,Bj=f˜j−2ηjωωj2−ω22+2ηjω2.The considered system consists of bilateral symmetrical elastic stops; the initial conditions and periodic motion conditions of the orbit with double grazing bifurcation satisfy the following formula:(8)x10=x10=d,x20=x20,v10=v10=0,v20=v20,x1πω=−d,x12pπω=x10=d,x22pπω=x20,v12pπω=v10=0,v22pπω=v20.By substituting the above conditions into the general solution, we can see that(9)d2=d12+d22≠0,τ=arctand1d2,where d1=A1φ11+A2φ12,d2=B1φ11+B2φ12.With appropriate system parameters, the motion of the system is periodic. The periodic motion of the system is represented byp,q,n, the number of periods is represented by n, p denotes the number of collisions between the mass M1 and the left side constraint, and q denotes the number of collisions between the mass M1 and the right side constraint.


## 3. Codimension-One Bifurcation Analysis

In this section, we will use the COCO-toolbox to study the dynamic behavior, giving the path-following analysis with respect to the frequency of the excitation of the system parameters. The state space of the system is divided into three regions according to the motion state of the oscillator. And then event functions are defined to describe the critical conditions of each mode. We takeu=x1,x˙1,x2,x˙2TεR4, as phase variable, and the values of the other nondimensional parameters are fixed as λ=μm,d,μk1,μc,ξ,μk2,f2,τ=0.25,0.1,0.1,0.25,0.2,100,0.5,0.9276 except for the frequency of the excitation.The segments needed for numerical simulation in COCO are shown in Table1.Table 1
The modes, vector fields, and event functions used for numerical simulation in COCO.
IntervalModeVector fieldEvent functionx1≤−dLeft contact (LC)u˙=fLCt,u,λhLC=x1+d−d<x1<dFree flight (FF)u˙=fFFt,u,λhLC=x1+d,hRC=x1−dx1≥dRight contact (RC)u˙=fRCt,u,λhRC=x1−dIn order to better show the dynamic behaviors of the oscillators under the case of elastic constraint, we take the frequency of excitation as a controlling parameter. Because the constraints on both sides of the system are symmetrical, we take the left contact time as an example. And we discuss the bifurcations of the system as the contact time between the massM1 and the left-hand side spring varies with the frequency of external excitation.In the following bifurcation diagram, the point on the solid line corresponds to a stable orbit, while the solid line on the phase diagram represents a stable solution. The vertical coordinate LTOC in the following bifurcation diagram represents the contact time between the massM1 and the spring on the left.In Figure2, we show the result of the numerical continuation of the period-1 solution with respect to the frequency ω. We start continuation from the larger value of ω; at this moment, there is no interaction between the oscillator M1 and the stop. Thus, the nonimpacting solution corresponds to a horizontal line with the LTOC being zero in Figure 2. From here, we can find a stable (0, 0, 1) solution. In the direction of decreasing ω, we find a grazing bifurcation point where the orbit makes tangential contact with both the left and right constraint surfaces (denoted as DGR1) at ω≈2.768 as shown in Figure 3(a), after which the impacting motion occurs. Very close to the grazing point, a saddle-node bifurcation is detected. Here, a Floquet multiplier of the periodic solution crosses the unit circle from the inside, and therefore the stability is lost. The orbit goes to the direction of increasing ω until the second saddle-node bifurcation SN2 is encountered at ω≈6.021. And the solution obtains stability; after this point, the orbit turns to the direction of decreasing ω. If we trace this stable branch, the third saddle-node bifurcation SN3 occurs at ω≈2.639, and the orbit loses stability again. If ω is decreased further, another double grazing bifurcation point DGR2 shown in Figure 3(b) is found. As shown in Figure 2(c), in the direction of increase of ω, the saddle-node bifurcation point SN4 and period-doubling bifurcation point PD1 are detected in turn. And the stability of the periodic solution changes. Along the unstable (2, 1, 1) orbit, we detect the period-doubling bifurcation point PD2 and the saddle-node bifurcation point SN5; finally, the stability is lost. This unstable impacting solution is traced further as shown in Figure 2(d). A saddle-node bifurcation point SN6 occurs, and the orbit regains stability. When ω is increased along the stable (2, 1, 1) orbit, the period-doubling bifurcation PD3 is found, and the stable orbit becomes unstable (2, 1, 1) orbit again. By tracing this unstable branch, we can detect the left grazing point where the solution is tangent to the left constraint (denoted as LGR) in (1, 1, 1) orbit at ω≈2.132. If ω is increased further, a period-doubling bifurcation PD4 in (1, 1, 1) orbit is encountered at ω≈2.435, and the period-doubling bifurcation makes the (1, 1, 1) orbit stable. If we trace this period orbit in the direction of increasing the frequency ω, the saddle-node point SN3 is detected. After this point, the branch of (1, 1, 1) orbit turns to the direction of decreasing ω until a period-doubling point PD5 appears, and then the solution loses the stability. As shown in Figure 2(e), as ω is decreased further along the unstable (1, 1, 1) branch, a grazing bifurcation where the orbit makes tangential contact with right constraint (denoted as RGR) is detected, after which the impact between the oscillator M1 and the right spring occurs, and then the unstable (1, 2, 1) branch emerges. In the direction of decreasing the frequency ω, a period-doubling bifurcation PD6 in (1, 2, 1) orbit makes the (1, 2, 1) branch stable again. When we decrease ω further, a saddle-node bifurcation SN7 in (1, 2, 1) orbit is detected, and the solution loses the stability. The branch has now switched to the direction of increasing the parameter ω. Along this branch, we find another saddle-node bifurcation SN8 in (1, 2, 1) orbit at ω≈2.373. If ω is decreased further, a period-doubling bifurcation PD7 is detected. Finally, the solution loses stability.Figure 2
The result of the numerical continuation of the period-1 solution with respect to the frequencyω. The black solid and black dashed orbits denote stable and unstable solutions, respectively. Bifurcation points are marked by black dots.Figure 3
(a) The phase diagram at the double grazing bifurcation point DGR1 of the period-1 orbit detected atω≈2.76811438. (b) The phase diagram at the double grazing bifurcation point DGR2 of the period-1 orbit detected at ω≈1.16136663.By tracing the period-1 solutions via the continuation platform COCO, we detect eight period-doubling points which give rise to four period-2 branches. The path-following analysis of the period-2 branch with respect to the frequencyω is shown in Figure 4.Figure 4
The result of the numerical continuation of the period-2 solution with respect to the frequencyω. The solid and dashed orbits denote stable and unstable solutions, respectively. Bifurcation points are marked by black dots. (a) Two orbits of period-2 originated from the period-doubling bifurcation points PD1 and PD3 of the period-1 orbit, respectively. (b) Two orbits of period-2 originated from the period-doubling bifurcation points PD5 and PD7 of the period-1 orbit, respectively.In Figure4(a), starting from the period-doubling bifurcation point PD1 in (2, 1, 1) orbit, and an unstable (4, 2, 2) orbit can be found through period-doubling bifurcation. In the vicinity of the period-doubling bifurcation point PD1, we detect a grazing point LGR1 after which the (4, 2, 2) orbit becomes (3, 2, 2) orbit. As ω is increased, a period-doubling bifurcation PD and a saddle-node bifurcation SN in (3, 2, 2) orbit are detected respectively, and the stability of periodic solution changes. If the frequency is reduced to ω≈2.064 a grazing bifurcation LGR2 in (3, 2, 2) orbit appears, where the solution makes tangential contact with the left stop, after which the branch of unstable (3, 2, 2) orbit becomes the branch of unstable (4, 2, 2) orbit. Later, as ω is increased, we find another saddle-node bifurcation SN and period-doubling bifurcation PD in (4, 2, 2) orbit; finally, the solution becomes stable. Along this stable branch, finally, it returns to the period-doubling bifurcation point PD2 of the period-1 orbit.We start continuation at the period-doubling bifurcation point PD3 of the period-1 branch, and a period-doubling bifurcation PD is detected whenω is increased. Below this point, we find a (4, 2, 2) branch, which terminates at the grazing bifurcation point LGR3. Along this (3, 2, 2) orbit, a period-doubling bifurcation point PD and a saddle-node bifurcation point SN are detected, respectively, and the stability of periodic solution changes. As ω is increased, we find another grazing bifurcation point LGR4 after which the (3, 2, 2) orbit becomes (2, 2, 2) orbit. Finally, the period-2 branch returns to the period-doubling point PD4 of period-1 orbit.The results of the numerical continuation of other period-2 impacting solutions with respect to the frequencyω are presented in Figure 4(b). Starting from the period-doubling bifurcation point PD5 of the period-1 orbit, we can detect a stable branch until the period-doubling bifurcation is reached at ω≈2.386, and then the orbit loses stability. When we trace this branch further in the direction of decreasing the frequency ω, a grazing bifurcation RGR1 is found, where the (2, 2, 2) orbit becomes (2, 3, 2) orbit. As ω is decreased, a saddle-node bifurcation occurs. Below this value, we find a small window of stability, which ends at the period-doubling bifurcation point PD, and then the branch loses stability. Along this unstable branch, another grazing bifurcation point RGR2 is detected where the (2, 3, 2) orbit becomes the (2, 4, 2) orbit. Finally, the orbit returns to the period-doubling bifurcation point PD6 of the period-1 orbit. A similar behavior can be observed around another period-2 orbit, which is shown in the upper part of Figure 4(b).


## 4. Codimension-Two Grazing Bifurcation

### 4.1. The Zero-Time-Discontinuity Mapping

Now suppose thatx∗1 and x∗2 are symmetrical grazing points of a state-space trajectory with the constraints D1 and D2, such that(10)hLCx∗1=0,hRCx∗2=0,v1∗1=hP1x∗1=∂hLCx∗1∂t=0,v1∗2=hP2x∗2=∂hRCx∗2∂t=0,a1∗1=∂2hLCx∗1∂t2<0,a2∗2=∂2hRCx∗2∂t2<0.We define two Poincare surfacesΠ1 and Π2: Π1=xεR5|θ=θ∗, Π2=xεR5|θ=θ∗+π, where the state vector x=x1,v1,x2,v2,θT, θ=ωt+τ, the modulus of θ being 2π, θ∗=τ∗, such that x∗1 and x∗2 correspond to a transversal intersection of a state-space trajectory with Poincare surfaces Π1 and Π2, respectively, and Φx∗1,t+T/2=Φx∗2,t, Φx∗2,t+T/2=Φx∗1,t. Thus, we can introduce two local Poincare maps P1,smooth:Π1⟶Π2 and P2,smooth:Π2⟶Π1, such that P1,smoothx∗1=x∗2 and P2,smoothx∗2=x∗1.According to the discontinuity mapping of the rigid constraint system introduced by Nordmark [4], the discontinuity mapping is extended to the elastic unilateral constraint system. As shown in Figure 5, there is a grazing trajectory with the discontinuous boundary. Consider the point x1 near the grazing point x∗; therefore, it can be seen that the trajectory x1,x2,x3,x4 passes through the discontinuous boundary. From the reverse direction of its trajectory, the time required to move from the point x3 to the point x0 is the same as that from the point x1 to the point x3. So, we assume that there is an instantaneous jump from x1 to x0; then, the zero-time-discontinuity mapping can be defined as the mapping from x1 to x0. For the specific process of discontinuity-mapping expression of the elastic constraint system, see [4].Figure 5
The sketch map of the zero-time-discontinuity mapping near the grazing trajectory.Based on the previous analysis, we extend the discontinuity mapping of elastic constraint systems from unilateral constraint to bilateral constraints. The two corresponding zero-time-discontinuity mappingsZDM1 and ZDM2 near two grazing points x∗1 and x∗2 can be written, respectively, as follows:(11)ZDM1=x,hminLCx≥0,x∗1+83−2a1∗1hminLC3/2Q,hminLCx<0,ZDM2=x,hmaxRCx<0,x∗2+83−2a1∗2hmaxRC3/2Q,hmaxRCx≥0,where Q=0,−μk0,0,0,0T,hmaxRCx represents the maximum displacement across the right collision surface, hminLCx represents the minimum displacement across the left collision surface.The two zero-time-discontinuity mappingsZDM1 and ZDM2 are defined in the vicinity of the grazing points x∗1 and x∗2, respectively, so the Poincare mapping P near the grazing trajectory has the following form:(12)P=P2,smooth∘ZDM2∘P1,smooth∘ZDM1.


### 4.2. The Poincare Mapping

The mappingP1,smooth can be expanded near the grazing point x∗1 as follows:(13)P1,smoothx=P1,smoothx∗1+P1,smooth,xx∗1x−x∗1+h.o.t=x∗2+P1,smooth,xx∗1x−x∗1+h.o.t.Similarly, we expand the event functionhLC near x∗1 and then hLC have the following form:(14)hLCx=hLCx∗1+hxLCx∗1x−x∗1+h.o.t.As the pointx∗1 is on the impact surface D1, it needs to satisfy hLCx∗1=0, such that (14) is expressed as(15)hLCx=hxLCx∗1x−x∗1+h.o.t,And, consequently,(16)P1=P1,smoothx°ZDM1x=x∗2+P1,smooth,xx∗1x−x∗1+h.o.t,hminLCx≥0,x∗2+P1,smooth,xx∗183−2a1∗1hminLC3/2Q+h.o.t,hminLCx<0.Similar to the above analysis, we have(17)P2=P2,smoothx°ZDM2x=x∗1+P2,smooth,xx∗2x−x∗2+h.o.t,hmaxRCx≤0,x∗1+P2,smooth,xx∗283−2a1∗2hmaxRC3/2Q+h.o.t,hmaxRCx>0.


### 4.3. Stability and Codimension-Two Grazing Bifurcation

#### 4.3.1. Conditions of Codimension-Two Grazing Bifurcation

If the starting point near the grazing point, which starts from either the impact side or nonimpacting side, converges to the grazing point after the iteration of mappings (16) and (17), the grazing periodic orbit is stable. For an impact point x in the neighborhood of the grazing point x∗1, it satisfies hxLCx∗1x−x∗1<0. After the iteration of the mapping P1, it reaches the nonimpacting side, which satisfies(18)hxRCx∗2P1x−x∗2=hxRCx∗2P1,smooth,xx∗183−2a1∗1hminLC3/2Q<0;that is,(19)hxRCx∗2P1,smooth,xx∗1Q<0.And after the iteration of the mappingP2, the oscillator collides with the discontinuity surface D1 again, which satisfies(20)hxLCx∗1P2P1x−x∗1=hxLCx∗1P2,smooth,xx∗2P1,smooth,xx∗183−2a1∗1hminLC3/2Q<0;that is,(21)hxLCx∗1P2,smooth,xx∗2P1,smooth,xx∗1Q<0.This shows that if the oscillator collides with the impact surfaceD1 again after iteration, the impact will be persisted, which results in a large stretch of the orbit in a certain direction. Therefore, the grazing periodic orbit is unstable.Based on the above analysis, if(22)hxRCx∗2P1,smooth,xx∗1P2,smooth,xx∗2P1,smooth,xx∗1j−1Q<0,hxLCx∗1P2,smooth,xx∗2P1,smooth,xx∗1iQ<0,2≤j≤i,it means that the oscillator changes from the impact side to the nonimpacting side after several iterations but finally collides with the impact surfaceD1 again and then the impact will be perpetuated. Therefore, the grazing periodic orbit loses stability.Similarly, if(23)hxRCx∗2P1,smooth,xx∗1Q>0,hxRCx∗2P1,smooth,xx∗1P2,smooth,xx∗2Q>0,the oscillator collides with the impact surface D2 again after iteration, and the impact will be perpetuated; then, the grazing periodic orbit is unstable.According to the above analysis, if(24)hxRCx∗2P1,smooth,xx∗1P2,smooth,xx∗2P1,smooth,xx∗1iQ>0,hxRCx∗2P1,smooth,xx∗1P2,smooth,xx∗2j+1Q>0,1≤i≤j,that indicates the oscillator changes from the impact side to the nonimpacting side after several iterations but finally impacts the impact surface D2 again and the impact will be perpetuated. Thus, the grazing periodic orbit loses stability.In the same way, if(25)hxRCx∗2P1,smooth,xx∗1Q>0,hxLCx∗1P2,smooth,xx∗2Q<0,then the oscillator impacts the impact surface D1 and D2 after iteration. If(26)hxRCx∗2P1,smooth,xx∗1P2,smooth,xx∗2P1,smooth,xx∗1iQ>0,hxLCx∗1P2,smooth,xx∗2P1,smooth,xx∗1P2,smooth,xx∗2jQ<0,1≤i≤j,this means that the oscillator collides again with the discontinuity surfaces D1 and D2 after several iterations, and the impact will be perpetuated. Thus, the grazing periodic orbit loses stability.According to the above analysis, four special cases can be obtained as follows:(27)hxLCx∗1P2,smooth,xx∗2P1,smooth,xx∗1iQ=0,hxLCx∗1P2,smooth,xx∗2P1,smooth,xx∗1P2,smooth,xx∗2jQ=0,hxRCx∗2P1,smooth,xx∗1P2,smooth,xx∗2P1,smooth,xx∗1iQ=0,hxRCx∗2P1,smooth,xx∗1P2,smooth,xx∗2j+1Q=0,0≤i,0≤j.Based on the definition in the literature [24], the points corresponding to the above four cases are the codimension-two grazing bifurcation points. Let ξn=0n=1,2,3,4 express the condition that the codimension-two grazing bifurcation points satisfy; then, it can be written in four cases as follows:(H1)ξ1=hxLCx∗1P2,smooth,xx∗2P1,smooth,xx∗1iQ=0(H2)ξ2=hxLCx∗1P2,smooth,xx∗2P1,smooth,xx∗1P2,smooth,xx∗2jQ=0(H3)ξ3=hxRCx∗2P1,smooth,xx∗1P2,smooth,xx∗2P1,smooth,xx∗1iQ=0(H4)ξ4=hxRCx∗2P1,smooth,xx∗1P2,smooth,xx∗2j+1Q=00≤i,0≤jNext, we take the third case (H3) as an example to simplify the calculation formula:(28)P1,smooth,xx∗1P2,smooth,xx∗2P1,smooth,xx∗1n=P11n+1P12n+1P13n+1P14n+1P15n+1P21n+1P22n+1P23n+1P24n+1P25n+1P31n+1P32n+1P33n+1P34n+1P35n+1P41n+1P42n+1P43n+1P44n+1P45n+1P51n+1P52n+1P53n+1P54n+1P55n+1.Thus,(29)ξ3=10000P11n+1P12n+1P13n+1P14n+1P15n+1P21n+1P22n+1P23n+1P24n+1P25n+1P31n+1P32n+1P33n+1P34n+1P35n+1P41n+1P42n+1P43n+1P44n+1P45n+1P51n+1P52n+1P53n+1P54n+1P55n+10−μk0000=−μk0P12n+1.Similarly, we can simplify the calculation formulae of codimension-two bifurcation points for the other three cases (H1), (H2), and (H4). Since the system is symmetrically constrained, we can get thatξ1 and ξ4, ξ2 and ξ3 are the same results, respectively, after simplification. The corresponding results are as follows:(30)ξ1=ξ4=−μk0P12n,ξ2=ξ3=−μk0P12n+1,where P12n and P12n+1 denote the (1, 2) elements of matrices P2,smooth,xx∗2P1,smooth,xx∗1n and P1,smooth,xx∗1P2,smooth,xx∗2P1,smooth,xx∗1n, respectively.


#### 4.3.2. Numerical Simulation

According to the theoretical results deduced above, the existence condition of codimension-two grazing bifurcation points isξn=0(n = 1, 2, 3, 4). If the first case (H1) is taken as an example, then the codimension-two grazing bifurcation points correspond to ξ1=−μk0P12n=0. Taking μk=0.25, μm=0.25, ζ=0.2, n = 1, when ω=0.1376,0.159,0.2091,0.3096,0.6189,…, we have ξ1=0. The grazing curve is shown in Figure 6. The points for which ξ1=0 are indicated by ‘∗’, which are the codimension-two grazing points. It can be seen from Figure 6 that the number of points of codimension-two grazing bifurcation increases with the decrease of frequency ω.Figure 6
The diagram of codimension-two grazing bifurcation points. The black solid line represents the grazing curve, and the points on the line represent the codimension-two grazing bifurcation points.Whenn=1, taking μk=0.25, μm=0.25, ζ=0.2, a codimension-two grazing bifurcation point corresponds to ω=0.3096, d=10.88, which satisfies the existence condition of grazing bifurcation. The grazing periodic trajectory is shown in Figure 7.Figure 7
A grazing periodic trajectory of oscillatorM1.In order to further study the dynamic behavior of the system in the region near a codimension-two grazing bifurcation point, firstly fixedμk=0.25, μm=0.25, ζ=0.2, μk0=100, and then ω=0.3096, d=10.88 corresponds to the codimension-two grazing bifurcation point. Changing the values of the parameters ω, d in the vicinity of ω=0.3096, d=10.88, different types of periodic-impact motions can be obtained by selecting a large number of parameters for numerical simulation. And finally, according to the dynamic behavior of the system, the parameter plane is divided into several regions. The unfolding diagram near the codimension-two grazing bifurcation point is drawn in Figure 8, where the dashed line in the diagram represents the grazing curve.Figure 8
Schematic of unfolding diagram near a codimension-two grazing bifurcation point.Selecting parameters in each area of Figure8 for numerical simulation, a series of phase diagrams can be obtained. Fixing ω=0.32, phase portraits obtained by perturbing d are shown in Figures 9(a)–9(i). When we decrease d, the phase portraits corresponding to the regions (a)–(i) in Figure 8 are the (1, 1, 3), (0, 2, 2), (0, 1, 1), (1, 0, 1), (1, 1, 1), (2, 2, 1), (3, 3, 1), (4, 4, 1), (5, 5, 1) orbits, respectively.Figure 9
The phase diagrams for periodic motions in regions (a)–(i) of Figure8. Poincare sections are marked by red dots on the phase planes.For fixedω=0.3, phase portraits obtained by perturbing the parameter are shown in Figures 10(a)–10(j) corresponding to (j)–(s) in Figure 8, respectively. When we decrease d, the phase portraits corresponding to the regions (j)–(s) in Figure 8 are the (0, 0, 1), (1, 1, 3), (2, 2, 3), (0, 1, 1), (1, 0, 1), (1, 1, 1), (2, 2, 1), (3, 3, 1), (4, 4, 1), (5, 5, 1) orbits, respectively.Figure 10
The phase diagrams in (a)-(j) correspond to the periodic motions of regions (j)-(s) in Figure8, respectively. Poincare sections are marked by red dots on the phase planes.


## 4.1. The Zero-Time-Discontinuity Mapping

Now suppose thatx∗1 and x∗2 are symmetrical grazing points of a state-space trajectory with the constraints D1 and D2, such that(10)hLCx∗1=0,hRCx∗2=0,v1∗1=hP1x∗1=∂hLCx∗1∂t=0,v1∗2=hP2x∗2=∂hRCx∗2∂t=0,a1∗1=∂2hLCx∗1∂t2<0,a2∗2=∂2hRCx∗2∂t2<0.We define two Poincare surfacesΠ1 and Π2: Π1=xεR5|θ=θ∗, Π2=xεR5|θ=θ∗+π, where the state vector x=x1,v1,x2,v2,θT, θ=ωt+τ, the modulus of θ being 2π, θ∗=τ∗, such that x∗1 and x∗2 correspond to a transversal intersection of a state-space trajectory with Poincare surfaces Π1 and Π2, respectively, and Φx∗1,t+T/2=Φx∗2,t, Φx∗2,t+T/2=Φx∗1,t. Thus, we can introduce two local Poincare maps P1,smooth:Π1⟶Π2 and P2,smooth:Π2⟶Π1, such that P1,smoothx∗1=x∗2 and P2,smoothx∗2=x∗1.According to the discontinuity mapping of the rigid constraint system introduced by Nordmark [4], the discontinuity mapping is extended to the elastic unilateral constraint system. As shown in Figure 5, there is a grazing trajectory with the discontinuous boundary. Consider the point x1 near the grazing point x∗; therefore, it can be seen that the trajectory x1,x2,x3,x4 passes through the discontinuous boundary. From the reverse direction of its trajectory, the time required to move from the point x3 to the point x0 is the same as that from the point x1 to the point x3. So, we assume that there is an instantaneous jump from x1 to x0; then, the zero-time-discontinuity mapping can be defined as the mapping from x1 to x0. For the specific process of discontinuity-mapping expression of the elastic constraint system, see [4].Figure 5
The sketch map of the zero-time-discontinuity mapping near the grazing trajectory.Based on the previous analysis, we extend the discontinuity mapping of elastic constraint systems from unilateral constraint to bilateral constraints. The two corresponding zero-time-discontinuity mappingsZDM1 and ZDM2 near two grazing points x∗1 and x∗2 can be written, respectively, as follows:(11)ZDM1=x,hminLCx≥0,x∗1+83−2a1∗1hminLC3/2Q,hminLCx<0,ZDM2=x,hmaxRCx<0,x∗2+83−2a1∗2hmaxRC3/2Q,hmaxRCx≥0,where Q=0,−μk0,0,0,0T,hmaxRCx represents the maximum displacement across the right collision surface, hminLCx represents the minimum displacement across the left collision surface.The two zero-time-discontinuity mappingsZDM1 and ZDM2 are defined in the vicinity of the grazing points x∗1 and x∗2, respectively, so the Poincare mapping P near the grazing trajectory has the following form:(12)P=P2,smooth∘ZDM2∘P1,smooth∘ZDM1.


## 4.2. The Poincare Mapping

The mappingP1,smooth can be expanded near the grazing point x∗1 as follows:(13)P1,smoothx=P1,smoothx∗1+P1,smooth,xx∗1x−x∗1+h.o.t=x∗2+P1,smooth,xx∗1x−x∗1+h.o.t.Similarly, we expand the event functionhLC near x∗1 and then hLC have the following form:(14)hLCx=hLCx∗1+hxLCx∗1x−x∗1+h.o.t.As the pointx∗1 is on the impact surface D1, it needs to satisfy hLCx∗1=0, such that (14) is expressed as(15)hLCx=hxLCx∗1x−x∗1+h.o.t,And, consequently,(16)P1=P1,smoothx°ZDM1x=x∗2+P1,smooth,xx∗1x−x∗1+h.o.t,hminLCx≥0,x∗2+P1,smooth,xx∗183−2a1∗1hminLC3/2Q+h.o.t,hminLCx<0.Similar to the above analysis, we have(17)P2=P2,smoothx°ZDM2x=x∗1+P2,smooth,xx∗2x−x∗2+h.o.t,hmaxRCx≤0,x∗1+P2,smooth,xx∗283−2a1∗2hmaxRC3/2Q+h.o.t,hmaxRCx>0.


## 4.3. Stability and Codimension-Two Grazing Bifurcation

### 4.3.1. Conditions of Codimension-Two Grazing Bifurcation

If the starting point near the grazing point, which starts from either the impact side or nonimpacting side, converges to the grazing point after the iteration of mappings (16) and (17), the grazing periodic orbit is stable. For an impact point x in the neighborhood of the grazing point x∗1, it satisfies hxLCx∗1x−x∗1<0. After the iteration of the mapping P1, it reaches the nonimpacting side, which satisfies(18)hxRCx∗2P1x−x∗2=hxRCx∗2P1,smooth,xx∗183−2a1∗1hminLC3/2Q<0;that is,(19)hxRCx∗2P1,smooth,xx∗1Q<0.And after the iteration of the mappingP2, the oscillator collides with the discontinuity surface D1 again, which satisfies(20)hxLCx∗1P2P1x−x∗1=hxLCx∗1P2,smooth,xx∗2P1,smooth,xx∗183−2a1∗1hminLC3/2Q<0;that is,(21)hxLCx∗1P2,smooth,xx∗2P1,smooth,xx∗1Q<0.This shows that if the oscillator collides with the impact surfaceD1 again after iteration, the impact will be persisted, which results in a large stretch of the orbit in a certain direction. Therefore, the grazing periodic orbit is unstable.Based on the above analysis, if(22)hxRCx∗2P1,smooth,xx∗1P2,smooth,xx∗2P1,smooth,xx∗1j−1Q<0,hxLCx∗1P2,smooth,xx∗2P1,smooth,xx∗1iQ<0,2≤j≤i,it means that the oscillator changes from the impact side to the nonimpacting side after several iterations but finally collides with the impact surfaceD1 again and then the impact will be perpetuated. Therefore, the grazing periodic orbit loses stability.Similarly, if(23)hxRCx∗2P1,smooth,xx∗1Q>0,hxRCx∗2P1,smooth,xx∗1P2,smooth,xx∗2Q>0,the oscillator collides with the impact surface D2 again after iteration, and the impact will be perpetuated; then, the grazing periodic orbit is unstable.According to the above analysis, if(24)hxRCx∗2P1,smooth,xx∗1P2,smooth,xx∗2P1,smooth,xx∗1iQ>0,hxRCx∗2P1,smooth,xx∗1P2,smooth,xx∗2j+1Q>0,1≤i≤j,that indicates the oscillator changes from the impact side to the nonimpacting side after several iterations but finally impacts the impact surface D2 again and the impact will be perpetuated. Thus, the grazing periodic orbit loses stability.In the same way, if(25)hxRCx∗2P1,smooth,xx∗1Q>0,hxLCx∗1P2,smooth,xx∗2Q<0,then the oscillator impacts the impact surface D1 and D2 after iteration. If(26)hxRCx∗2P1,smooth,xx∗1P2,smooth,xx∗2P1,smooth,xx∗1iQ>0,hxLCx∗1P2,smooth,xx∗2P1,smooth,xx∗1P2,smooth,xx∗2jQ<0,1≤i≤j,this means that the oscillator collides again with the discontinuity surfaces D1 and D2 after several iterations, and the impact will be perpetuated. Thus, the grazing periodic orbit loses stability.According to the above analysis, four special cases can be obtained as follows:(27)hxLCx∗1P2,smooth,xx∗2P1,smooth,xx∗1iQ=0,hxLCx∗1P2,smooth,xx∗2P1,smooth,xx∗1P2,smooth,xx∗2jQ=0,hxRCx∗2P1,smooth,xx∗1P2,smooth,xx∗2P1,smooth,xx∗1iQ=0,hxRCx∗2P1,smooth,xx∗1P2,smooth,xx∗2j+1Q=0,0≤i,0≤j.Based on the definition in the literature [24], the points corresponding to the above four cases are the codimension-two grazing bifurcation points. Let ξn=0n=1,2,3,4 express the condition that the codimension-two grazing bifurcation points satisfy; then, it can be written in four cases as follows:(H1)ξ1=hxLCx∗1P2,smooth,xx∗2P1,smooth,xx∗1iQ=0(H2)ξ2=hxLCx∗1P2,smooth,xx∗2P1,smooth,xx∗1P2,smooth,xx∗2jQ=0(H3)ξ3=hxRCx∗2P1,smooth,xx∗1P2,smooth,xx∗2P1,smooth,xx∗1iQ=0(H4)ξ4=hxRCx∗2P1,smooth,xx∗1P2,smooth,xx∗2j+1Q=00≤i,0≤jNext, we take the third case (H3) as an example to simplify the calculation formula:(28)P1,smooth,xx∗1P2,smooth,xx∗2P1,smooth,xx∗1n=P11n+1P12n+1P13n+1P14n+1P15n+1P21n+1P22n+1P23n+1P24n+1P25n+1P31n+1P32n+1P33n+1P34n+1P35n+1P41n+1P42n+1P43n+1P44n+1P45n+1P51n+1P52n+1P53n+1P54n+1P55n+1.Thus,(29)ξ3=10000P11n+1P12n+1P13n+1P14n+1P15n+1P21n+1P22n+1P23n+1P24n+1P25n+1P31n+1P32n+1P33n+1P34n+1P35n+1P41n+1P42n+1P43n+1P44n+1P45n+1P51n+1P52n+1P53n+1P54n+1P55n+10−μk0000=−μk0P12n+1.Similarly, we can simplify the calculation formulae of codimension-two bifurcation points for the other three cases (H1), (H2), and (H4). Since the system is symmetrically constrained, we can get thatξ1 and ξ4, ξ2 and ξ3 are the same results, respectively, after simplification. The corresponding results are as follows:(30)ξ1=ξ4=−μk0P12n,ξ2=ξ3=−μk0P12n+1,where P12n and P12n+1 denote the (1, 2) elements of matrices P2,smooth,xx∗2P1,smooth,xx∗1n and P1,smooth,xx∗1P2,smooth,xx∗2P1,smooth,xx∗1n, respectively.


### 4.3.2. Numerical Simulation

According to the theoretical results deduced above, the existence condition of codimension-two grazing bifurcation points isξn=0(n = 1, 2, 3, 4). If the first case (H1) is taken as an example, then the codimension-two grazing bifurcation points correspond to ξ1=−μk0P12n=0. Taking μk=0.25, μm=0.25, ζ=0.2, n = 1, when ω=0.1376,0.159,0.2091,0.3096,0.6189,…, we have ξ1=0. The grazing curve is shown in Figure 6. The points for which ξ1=0 are indicated by ‘∗’, which are the codimension-two grazing points. It can be seen from Figure 6 that the number of points of codimension-two grazing bifurcation increases with the decrease of frequency ω.Figure 6
The diagram of codimension-two grazing bifurcation points. The black solid line represents the grazing curve, and the points on the line represent the codimension-two grazing bifurcation points.Whenn=1, taking μk=0.25, μm=0.25, ζ=0.2, a codimension-two grazing bifurcation point corresponds to ω=0.3096, d=10.88, which satisfies the existence condition of grazing bifurcation. The grazing periodic trajectory is shown in Figure 7.Figure 7
A grazing periodic trajectory of oscillatorM1.In order to further study the dynamic behavior of the system in the region near a codimension-two grazing bifurcation point, firstly fixedμk=0.25, μm=0.25, ζ=0.2, μk0=100, and then ω=0.3096, d=10.88 corresponds to the codimension-two grazing bifurcation point. Changing the values of the parameters ω, d in the vicinity of ω=0.3096, d=10.88, different types of periodic-impact motions can be obtained by selecting a large number of parameters for numerical simulation. And finally, according to the dynamic behavior of the system, the parameter plane is divided into several regions. The unfolding diagram near the codimension-two grazing bifurcation point is drawn in Figure 8, where the dashed line in the diagram represents the grazing curve.Figure 8
Schematic of unfolding diagram near a codimension-two grazing bifurcation point.Selecting parameters in each area of Figure8 for numerical simulation, a series of phase diagrams can be obtained. Fixing ω=0.32, phase portraits obtained by perturbing d are shown in Figures 9(a)–9(i). When we decrease d, the phase portraits corresponding to the regions (a)–(i) in Figure 8 are the (1, 1, 3), (0, 2, 2), (0, 1, 1), (1, 0, 1), (1, 1, 1), (2, 2, 1), (3, 3, 1), (4, 4, 1), (5, 5, 1) orbits, respectively.Figure 9
The phase diagrams for periodic motions in regions (a)–(i) of Figure8. Poincare sections are marked by red dots on the phase planes.For fixedω=0.3, phase portraits obtained by perturbing the parameter are shown in Figures 10(a)–10(j) corresponding to (j)–(s) in Figure 8, respectively. When we decrease d, the phase portraits corresponding to the regions (j)–(s) in Figure 8 are the (0, 0, 1), (1, 1, 3), (2, 2, 3), (0, 1, 1), (1, 0, 1), (1, 1, 1), (2, 2, 1), (3, 3, 1), (4, 4, 1), (5, 5, 1) orbits, respectively.Figure 10
The phase diagrams in (a)-(j) correspond to the periodic motions of regions (j)-(s) in Figure8, respectively. Poincare sections are marked by red dots on the phase planes.


## 4.3.1. Conditions of Codimension-Two Grazing Bifurcation

If the starting point near the grazing point, which starts from either the impact side or nonimpacting side, converges to the grazing point after the iteration of mappings (16) and (17), the grazing periodic orbit is stable. For an impact point x in the neighborhood of the grazing point x∗1, it satisfies hxLCx∗1x−x∗1<0. After the iteration of the mapping P1, it reaches the nonimpacting side, which satisfies(18)hxRCx∗2P1x−x∗2=hxRCx∗2P1,smooth,xx∗183−2a1∗1hminLC3/2Q<0;that is,(19)hxRCx∗2P1,smooth,xx∗1Q<0.And after the iteration of the mappingP2, the oscillator collides with the discontinuity surface D1 again, which satisfies(20)hxLCx∗1P2P1x−x∗1=hxLCx∗1P2,smooth,xx∗2P1,smooth,xx∗183−2a1∗1hminLC3/2Q<0;that is,(21)hxLCx∗1P2,smooth,xx∗2P1,smooth,xx∗1Q<0.This shows that if the oscillator collides with the impact surfaceD1 again after iteration, the impact will be persisted, which results in a large stretch of the orbit in a certain direction. Therefore, the grazing periodic orbit is unstable.Based on the above analysis, if(22)hxRCx∗2P1,smooth,xx∗1P2,smooth,xx∗2P1,smooth,xx∗1j−1Q<0,hxLCx∗1P2,smooth,xx∗2P1,smooth,xx∗1iQ<0,2≤j≤i,it means that the oscillator changes from the impact side to the nonimpacting side after several iterations but finally collides with the impact surfaceD1 again and then the impact will be perpetuated. Therefore, the grazing periodic orbit loses stability.Similarly, if(23)hxRCx∗2P1,smooth,xx∗1Q>0,hxRCx∗2P1,smooth,xx∗1P2,smooth,xx∗2Q>0,the oscillator collides with the impact surface D2 again after iteration, and the impact will be perpetuated; then, the grazing periodic orbit is unstable.According to the above analysis, if(24)hxRCx∗2P1,smooth,xx∗1P2,smooth,xx∗2P1,smooth,xx∗1iQ>0,hxRCx∗2P1,smooth,xx∗1P2,smooth,xx∗2j+1Q>0,1≤i≤j,that indicates the oscillator changes from the impact side to the nonimpacting side after several iterations but finally impacts the impact surface D2 again and the impact will be perpetuated. Thus, the grazing periodic orbit loses stability.In the same way, if(25)hxRCx∗2P1,smooth,xx∗1Q>0,hxLCx∗1P2,smooth,xx∗2Q<0,then the oscillator impacts the impact surface D1 and D2 after iteration. If(26)hxRCx∗2P1,smooth,xx∗1P2,smooth,xx∗2P1,smooth,xx∗1iQ>0,hxLCx∗1P2,smooth,xx∗2P1,smooth,xx∗1P2,smooth,xx∗2jQ<0,1≤i≤j,this means that the oscillator collides again with the discontinuity surfaces D1 and D2 after several iterations, and the impact will be perpetuated. Thus, the grazing periodic orbit loses stability.According to the above analysis, four special cases can be obtained as follows:(27)hxLCx∗1P2,smooth,xx∗2P1,smooth,xx∗1iQ=0,hxLCx∗1P2,smooth,xx∗2P1,smooth,xx∗1P2,smooth,xx∗2jQ=0,hxRCx∗2P1,smooth,xx∗1P2,smooth,xx∗2P1,smooth,xx∗1iQ=0,hxRCx∗2P1,smooth,xx∗1P2,smooth,xx∗2j+1Q=0,0≤i,0≤j.Based on the definition in the literature [24], the points corresponding to the above four cases are the codimension-two grazing bifurcation points. Let ξn=0n=1,2,3,4 express the condition that the codimension-two grazing bifurcation points satisfy; then, it can be written in four cases as follows:(H1)ξ1=hxLCx∗1P2,smooth,xx∗2P1,smooth,xx∗1iQ=0(H2)ξ2=hxLCx∗1P2,smooth,xx∗2P1,smooth,xx∗1P2,smooth,xx∗2jQ=0(H3)ξ3=hxRCx∗2P1,smooth,xx∗1P2,smooth,xx∗2P1,smooth,xx∗1iQ=0(H4)ξ4=hxRCx∗2P1,smooth,xx∗1P2,smooth,xx∗2j+1Q=00≤i,0≤jNext, we take the third case (H3) as an example to simplify the calculation formula:(28)P1,smooth,xx∗1P2,smooth,xx∗2P1,smooth,xx∗1n=P11n+1P12n+1P13n+1P14n+1P15n+1P21n+1P22n+1P23n+1P24n+1P25n+1P31n+1P32n+1P33n+1P34n+1P35n+1P41n+1P42n+1P43n+1P44n+1P45n+1P51n+1P52n+1P53n+1P54n+1P55n+1.Thus,(29)ξ3=10000P11n+1P12n+1P13n+1P14n+1P15n+1P21n+1P22n+1P23n+1P24n+1P25n+1P31n+1P32n+1P33n+1P34n+1P35n+1P41n+1P42n+1P43n+1P44n+1P45n+1P51n+1P52n+1P53n+1P54n+1P55n+10−μk0000=−μk0P12n+1.Similarly, we can simplify the calculation formulae of codimension-two bifurcation points for the other three cases (H1), (H2), and (H4). Since the system is symmetrically constrained, we can get thatξ1 and ξ4, ξ2 and ξ3 are the same results, respectively, after simplification. The corresponding results are as follows:(30)ξ1=ξ4=−μk0P12n,ξ2=ξ3=−μk0P12n+1,where P12n and P12n+1 denote the (1, 2) elements of matrices P2,smooth,xx∗2P1,smooth,xx∗1n and P1,smooth,xx∗1P2,smooth,xx∗2P1,smooth,xx∗1n, respectively.


## 4.3.2. Numerical Simulation

According to the theoretical results deduced above, the existence condition of codimension-two grazing bifurcation points isξn=0(n = 1, 2, 3, 4). If the first case (H1) is taken as an example, then the codimension-two grazing bifurcation points correspond to ξ1=−μk0P12n=0. Taking μk=0.25, μm=0.25, ζ=0.2, n = 1, when ω=0.1376,0.159,0.2091,0.3096,0.6189,…, we have ξ1=0. The grazing curve is shown in Figure 6. The points for which ξ1=0 are indicated by ‘∗’, which are the codimension-two grazing points. It can be seen from Figure 6 that the number of points of codimension-two grazing bifurcation increases with the decrease of frequency ω.Figure 6
The diagram of codimension-two grazing bifurcation points. The black solid line represents the grazing curve, and the points on the line represent the codimension-two grazing bifurcation points.Whenn=1, taking μk=0.25, μm=0.25, ζ=0.2, a codimension-two grazing bifurcation point corresponds to ω=0.3096, d=10.88, which satisfies the existence condition of grazing bifurcation. The grazing periodic trajectory is shown in Figure 7.Figure 7
A grazing periodic trajectory of oscillatorM1.In order to further study the dynamic behavior of the system in the region near a codimension-two grazing bifurcation point, firstly fixedμk=0.25, μm=0.25, ζ=0.2, μk0=100, and then ω=0.3096, d=10.88 corresponds to the codimension-two grazing bifurcation point. Changing the values of the parameters ω, d in the vicinity of ω=0.3096, d=10.88, different types of periodic-impact motions can be obtained by selecting a large number of parameters for numerical simulation. And finally, according to the dynamic behavior of the system, the parameter plane is divided into several regions. The unfolding diagram near the codimension-two grazing bifurcation point is drawn in Figure 8, where the dashed line in the diagram represents the grazing curve.Figure 8
Schematic of unfolding diagram near a codimension-two grazing bifurcation point.Selecting parameters in each area of Figure8 for numerical simulation, a series of phase diagrams can be obtained. Fixing ω=0.32, phase portraits obtained by perturbing d are shown in Figures 9(a)–9(i). When we decrease d, the phase portraits corresponding to the regions (a)–(i) in Figure 8 are the (1, 1, 3), (0, 2, 2), (0, 1, 1), (1, 0, 1), (1, 1, 1), (2, 2, 1), (3, 3, 1), (4, 4, 1), (5, 5, 1) orbits, respectively.Figure 9
The phase diagrams for periodic motions in regions (a)–(i) of Figure8. Poincare sections are marked by red dots on the phase planes.For fixedω=0.3, phase portraits obtained by perturbing the parameter are shown in Figures 10(a)–10(j) corresponding to (j)–(s) in Figure 8, respectively. When we decrease d, the phase portraits corresponding to the regions (j)–(s) in Figure 8 are the (0, 0, 1), (1, 1, 3), (2, 2, 3), (0, 1, 1), (1, 0, 1), (1, 1, 1), (2, 2, 1), (3, 3, 1), (4, 4, 1), (5, 5, 1) orbits, respectively.Figure 10
The phase diagrams in (a)-(j) correspond to the periodic motions of regions (j)-(s) in Figure8, respectively. Poincare sections are marked by red dots on the phase planes.


## 5. Conclusions

For a two-degree-of-freedom system with symmetric elastic constraints, the codimension-one bifurcation of the system is discussed in detail by using the continuation method, and grazing bifurcation, saddle-node bifurcation, and period-doubling bifurcation are obtained. As far as we knew, there are little research on the codimension-two grazing bifurcations of the system with the symmetric elastic constraints. Then, based on the traditional method of discontinuity mapping under the single constraint, we extend it to the impact oscillator with two-sided elastic constraints, and the specific mathematical expressions of zero-time-discontinuity mapping of the symmetric elastic constraint system are deduced after a complicated calculation. By combining the zero-time-discontinuity mappings and local smooth Poincare maps, the two composite mappingsP1, P2 are obtained. According to the composite mappings, the criteria of codimension-two grazing bifurcation point under four different conditions are given out. The grazing curve and codimension-two grazing bifurcation points are shown by numerical simulation. And a series of complex dynamic behaviors in the vicinity of the codimension-two grazing bifurcation point are also presented, which reveal the rich dynamic behaviors of the impact system with two-sided elastic constraints.



---

*Source: 8854115-2021-04-01.xml*