# Periodic Motion and Transition of a Vibro‐Impact System with Multilevel Elastic Constraints

Authors: Jie Ding; Chao Wang; Wangcai Ding
Journal: Discrete Dynamics in Nature and Society
 (2021)

---

## Abstract

In this paper, a single-degree-of-freedom vibroimpact system with multilevel elastic constraints is taken as the research object. By constructing the Poincaré map of the system and calculating the Lyapunov exponent spectrum of the system, the stability of the system is determined. Using the multiparameter collaborative numerical simulation method, the parameter domains of various periodic motions are determined, and the diversity and transition characteristics of periodic motions are revealed. At the same time, combined with the cell mapping method, the coexistence of attractors induced due to grazing bifurcation, saddle-node bifurcation, and boundary crisis is studied. Finally, the influence of system parameters on periodic motion distribution is analyzed, which provides a scientific basis for system parameter optimization.


---

## Body

## 1. Introduction

The vibroimpact system with clearance, as a typical nonsmooth dynamic system, is widely used in engineering practice, such as vibrating sand shaker [1], vibration isolation system with actuators [2], impact damping device [3], railway train wheel-rail vibration [4], and vibroimpact capsule system [5]. This kind of system exhibits rich and complex nonsmooth dynamic behavior due to the sudden change of stiffness caused by clearance and elastic constraints.In the past few decades, many scholars have conducted extensive research on the segmented smooth impact vibration system with clearance, using analytical and numerical methods. Shaw and Holmes [6] first studied the piecewise linear oscillator under simple harmonic excitation with the view of modern dynamic system, analyzed the partial bifurcation of periodic motion by the central manifold theorem, and discussed the chaotic motion by the homoclinic phase intercept condition. Natsiavas [7, 8] proposed a joint method with a relatively small computational workload and studied the resonance response of symmetric and asymmetric piecewise linear oscillators excited by simple harmonics. Choi and Noah [9] introduced the fast Fourier transform in the harmonic balance method to calculate the harmonic coefficients of each order of the nonlinear restoring force and then studied the multiple resonance responses of the asymmetric elastic constrained oscillator. Ketema [10] applied the Melnikov method to the determination of the existence of the homoclinic point of the Poincaré mapping of the piecewise linear oscillator system. Hu [11] used a combination of segmented analytical solutions and shooting methods to improve the accuracy and efficiency of calculating the periodic motion of the elastic restraint system. Luo [12] used a nonsmooth analysis method to study various unstable and stable periodic motion patterns of the piecewise linear periodic excitation system.The instability, caused by the grazing impacts in the vibration shock system, is a hot spot of current research. As a pioneer in this field, Nordmark [13] studied the nonperiodic motion caused by the grazing impacts in the piecewise linear system. Chin et al. [14] analyzed a single-degree-of-freedom vibration system and proposed several different types of edge rubbing and bifurcation that caused changes in system dynamics. Chillingworth [15] studied the dynamic behavior of a single-degree-of-freedom shock oscillator near the edge rubbing point, including nondegenerate (secondary rubbing) and minimal degradation (three rubbing) phenomena. Humphries et al. [16] introduced a discontinuous geometry topology method to explain the cause of saddle-knot bifurcation near the edge rubbing bifurcation. Kryzhevich and Wiercigroch [17] studied the existence of nonhyperbolic homoclinic points of the periodic solutions of vibration shock systems and introduced a new type of robust chaotic dynamics for strongly nonlinear systems. Pavlovskaia et al. [18] used a large number of experiments to study the impact vibrator with unilateral elastic restraint and revealed the influence of different excitation frequencies near the edge rubbing point on the bifurcation. The result showed that there was attractors coexistence phenomenon when there was discontinuous transition among the two orbits via boundary crisis. Jiang et al. [19] used the path tracking method to study the dynamic behavior of the grazing bifurcation of the elastic and rigid impact vibrators. Under the constraints of the unilateral elastic impact vibrator, the smooth bifurcation (period-doubling, saddle-node) point was close to nonsmooth (grazing) bifurcation points, and there was no smooth bifurcation near the nonsmooth bifurcation point under the constraints of the unilateral rigid impact oscillator.The cell mapping method, created by the American mathematician Hsu [20] in the 1980s, is an effective tool for studying the global characteristics of nonsmooth dynamic systems. The advantage of this method is that it can efficiently and dynamically analyze the attractors, attraction domain, and chaos of the system. The behavior more comprehensively reveals the mechanism of the periodic movement of the system. Based on the idea of the cell mapping method, many scholars have proposed simple cell mapping [21], interpolated cell mapping [22], graph cell mapping [23], and other improved cell mapping methods and obtained a lot of results in the quantitative and qualitative analysis of the global dynamic characteristics of nonsmooth dynamic systems. Kong et al. [24] used the generalized cell mapping method to study the global dynamics of piecewise linear systems excited by a Gaussian white noise. Li et al. [25, 26] analyzed the diversity and evolution of basic periodic shock motions by combining the cell mapping method on the simple harmonic vibration dynamic model with gaps and studied the distribution of attractors and attractive domains. Chong et al. [27] studied in detail the influence of different parameters on system dynamics by combining the path following method and the cell mapping method for the harmonically excited linear oscillator with gaps.At present, most actual mechanical systems with clearances are multiparameter and multiconstraint systems, and the coordinated changes of multiple parameters show the dynamic behavior of the system more comprehensively. However, the existing research is mostly limited to the case of a single parameter. This paper not only considers multiple constraint factors but also conducts numerical simulation under the coordinated change of multiple parameters. At the same time, combined with the cell mapping method, the dynamic characteristics of a class of impact vibration systems with multilevel elastic constraints are studied. The periodic motion distribution and the transition process of periodic motion through various nonsmooth bifurcations are analyzed. The structure of this paper is as follows. In the second section, the physical model and motion equations of the collision system with two elastic constraints are established. In Section3, the Poincaré map and Lyapunov exponent of the system are derived, which can accurately determine the periodic motion and chaotic motion in the parameter plane. In Section 4, the distribution area of periodic motion and the transition law between adjacent periodic motions are studied through a multiparameter collaborative simulation, and the coexistence of periodic motion in the transition process of periodic motion is studied with the method of cell mapping. Section 5 studies the influence of system parameters on the distribution of periodic motion. Section 6 gives a summary of this article.


## 2. Mathematical Modelling

The physical model of a single-degree-of-freedom vibroimpact system with multilevel elastic constraints is shown in Figure1. A vibrator with a mass M is connected by a linear spring with a stiffness K and a linear damper with a damping coefficient C and moves along the horizontal direction under the action of harmonic exciting force F0sinΩT+τ. The coordinate system is established by taking the static equilibrium position of the system as the origin of the spatial coordinates, the displacement of the vibrator is X, and the elastic constraints with stiffnesses K1 and K2 are fixed at positions B1 and B1+B2 to the right of the vibrator. Under the action of the harmonic exciting force, when the displacement of the vibrator is small, the system is a single-degree-of-freedom linear vibration system; when the displacement of the vibrator is greater than B1 or B1+B2, the vibrator will interact with one or two elastic constraints on the right side. The system will exhibit strong nonlinearity and nonsmooth characteristics due to the existence of clearances and the occurrence of collisions.Figure 1
Physical model of a single-degree-of-freedom vibroimpact system with multilevel elastic constraints.The massM is subjected by the spring force, damping force, and external exciting force at the same time. According to Newton’s theorem, the differential equations of system motion can be expressed as(1)MX¨+CX˙+KX=F0sinΩT+τ,X≤B1,MX¨+CX˙+KX+K1X−B1=F0sinΩT+τ,B1<X<B2,MX¨+CX˙+KX+K1X−B1+K2X−B1−B2=F0sinΩT+τ,X≥B1+B2.Due to the different dimensions of these physical quantities, in order to thoroughly study the motion mechanism of the mechanical model and eliminate the influence of dimensions, it is necessary to perform nondimensional processing on the differential equation of motion and introduce the following nondimensional quantities:x=XK/F0, b1=B1K/F0, b2=B2K/F0, ζ=C/2MK, ω=ΩM/K, t=ΤK/M, μk1=K1/K, and μk2=K2/K.The nondimensional differential equations of the system motion are(2)x¨+2ζx˙+x=sinωt+τ,x≤b1,x¨+2ζx˙+1+μk1=sinωt+τ+μk1b1,b1<x<b1+b2,x¨+2ζx˙+1+μk1+μk2=sinωt+τ+μk1b1+μk2b1+b2,x≥b1+b2.The general solutions of the set of differential equations (2) are(3)xt=e−η1t−t0a1cosωd1t−t0+c1sinωd1t−t0+A1sinωt+τ0+C1cosωt+τ0,x≤b,xt=e−η2t−t1a2cosωd2t−t1+c2sinωd2t−t1+A2sinωt+τ0+C2cosωt+τ0+μk1b11+μk1,b1<x<b1+b2,xt=−e−η3t−t2a3cosωd3t−t2+c3sinωd3t−t2+A3sinωt+τ0+C3cosωt+τ0+μk1b1+μk2b1+b21+μk1+μk2,x≥b1+b2,where ωn12=1, ωn22=1+μk, ωn32=1+μk1+μk2, η1=ζ, η2=ζ, η3=ζ, ωd1=ωn12−η12, ωd2=ωn22−η22, and ωd3=ωn32−η32; a1, a2, a3, c1, c2, and c3 are the integral constants, respectively, which are determined by the initial conditions of the vibration system; A1, A2, A3, C1, C2, and C3 are the amplitude constants, respectively, which have(4)A1=1−ω21−ω22+2ξω2,C1=−2ξω1−ω22+2ξω2,A2=1+μk−ω21+μk−ω22+2ξω2,C2=−2ξω1+μk−ω22+2ξω2,A3=1+μk1+μk2−ω21+μk1+μk2−ω22+2ξω2,C3=−2ξω1+μk1+μk2−ω22+2ξω2.


## 3. Poincaré Mapping of the System

In order to study the stability, diversity, existence area, and local bifurcation of the periodic motion of the system, the Poincaré mapping of the system is established.Select the mapping sectionΣT=x,x˙,t∈R2×S,x=x−,modt=2π/ω to count the number of cycles of periodic motion of the system; use collision surfaces Σp1=x,x˙,θ∈R2×S,x=b1,x˙>0 and Σp2=x,x˙,θ∈R2×S,x=b1+b2,x˙>0 to count the number of collisions between the oscillator and the constrained surface.The Poincaré mapping diagram of the system is shown in Figure2. The collision surface Σp1 is used to establish the Poincaré mapping of the system: P:Σp1⟶Σp1; the eigenvalues of the linearization matrix of the Poincaré mapping at the fixed points are used to study the stability of the system’s periodic motions. According to the contact process, Poincaré mapping is divided into five stages: P1: the vibrator does not contact with the linear springs K1 and K2; P2: the vibrator contacts with spring K1 but not with spring K2; P3: the vibrator contacts with spring K2; P4: the vibrator changes direction and returns after contacting with spring K1; P5: the vibrator is detached from spring K2. Let DP1, DP2, DP3, DP4, and DP5 denote the mapping matrices of P1, P2, P3, P4, and P5; among them,(5)DPk=∂fki∂Yjx˙k−1,τk−1;i,j=1,2,k=1,2,3,DPk=∂fki∂Yjx˙m,τm;i,j=1,2,k=4,m=1ork=5,m=3.Figure 2
Schematic diagram of Poincaré mapping of the system.In the above formulas,fi represents the i component of Pk and Yj represents the j component of the coordinate. First, the integral constants al and bll=1,2,3 are obtained from the initial conditions of the initial stages of each movement, and then, according to the termination conditions, xt1−b1=0, xt2−b1−b2=0, xt3−b1−b2=0, and xt4−b1=0. Then the Poincaré mapping of the system’s n−p−q period movement is(6)fk1x˙k−1,τk−1=x˙t,fk2x˙k−1,τk−1=ωt+τk−1,k≥3,fk1x˙m,τm=x˙t,fk2x˙m,τm=ωt+τm,k=4,m=1ork=5,m=3.Define the function:(7)Gktk,x˙k−1,τk−1=xtk−b1,k≥3,Gktk,x˙m,τm=xtk−b1,k=4,m=1ork=5,m=3.According to the implicit function derivation theorem, each element inDPk can be found:(8)akij=∂fki∂Yj0−∂fki∂t∂Gk∂Yj0/∂Gk∂t,i,j=1,2.The Lyapunov exponential spectrum is an effective tool for judging the stability and chaos of the dynamic system. In this paper, the Poincaré mapping method is used to transform the continuous dynamic system into a discrete dynamic system. The calculation formula is as follows:(9)λi=1K∑k=1Klnzik,i=1,2,where K is the number of iterations, zik is the vector processed by Gram-Schmidt orthogonalization and norm normalization after each iteration, and ziki=1,2 is the norm of zik.


## 4. The Transition Law of Basic Periodic Motion

In order to study the distribution type of the periodic motion of the system and its transition law, the gapb1 and frequency ω are taken as the research objects. Take the system parameters ζ=0.1, μk1=20, μk2=50, ω∈0.1,2, b1∈0.1,2, and b2=0.1. The periodic motion distribution of the system in the parameter plane ω,b1 is shown in Figure 3 and combined with the Lyapunov exponent to determine the existence of stable periodic motion and chaotic motion. In Figure 3, different colors and corresponding symbols are used to indicate the existence area and distribution law of various periodic motions. In symbol n−p−q, n is the number of excitation cycles, p is the number of times the vibrator contacts with the elastic constraint of right stiffness K1 at forward speed, and q is the number of times the vibrator contacts with the elastic constraint of right stiffness K2 at forward speed. PF is used for period-doubling bifurcation, IPF for inverse period-doubling bifurcation, SN for saddle-node bifurcation, GR for edge rubbing bifurcation, and BC for the boundary crisis of the system. The areas not shown in the figure mainly include the following: ① the parameter field which is too small and the display which is not obvious enough; ② chaos area. In order to express the type of grazing bifurcation between the vibrator and different constraint surfaces, the following definitions are made:(a)
The grazing bifurcation between the vibrator and the elastic constraint surface with a stiffness coefficientK1 is defined as the first type of grazing bifurcation, represented by GR1, and the corresponding grazing motion is defined as the first type of grazing periodic motion.(b)
The grazing bifurcation between the vibrator and the elastic constraint surface with a stiffness coefficientK2 is defined as the second type of grazing bifurcation, represented by GR2, and the corresponding grazing motion is defined as the second type of grazing periodic motion.Figure 3
Distribution diagram of the periodic motion of the system in the parameter planeω,b1. (a) Overall view. (b) Partial enlarged views 1. (c) Transitional zone between 1−0−0 and 1−1−1 periodic motions. (d) Partial enlarged views 2 of Figure 3(a). (e) Partial enlarged views 3 of Figure 3(a). (f) Partial enlarged views 4 of Figure 3(a).
(a)(b)(c)(d)(e)(f)It can be seen from Figure3(a) that, with the increase or decrease of b1 and ω, the basic periodic motion 1−p−q mainly has the following forms: p, the number of collisions between the vibrator and the elastic restraint surface with a stiffness coefficient K1, is fixed, and q, the number of collisions with the elastic restraint surface with a stiffness coefficient K2, gradually increases or decreases by one. The number of collisions q between the vibrator and the elastic restraining surface with a stiffness coefficient K2 is fixed; the number of collisions p with the elastic restraint surface with a stiffness coefficient K1 gradually increases or decreases by one. When the number of collisions gradually increases, the system will flutter. There are two main types of flutter: the first one is that flutters occur only between the oscillator and the elastic constrained surface with a stiffness coefficient K1, which is expressed as 1−p˜−q; the second type is that flutters occur between the oscillator and two kinds of elastic constrained surfaces, which is expressed as 1−p˜−q˜. When both p and q are 0, the vibrator does not collide with the two types of elastic constraints on the right side.


### 4.1. Analysis of the Basic Transition Law between Adjacent Periodic Motions

In the low-frequency small gap area, as shown in Figure3(b), the 1−p−q periodic motion can be transferred to 1−p+1−q and 1−p−q+1 periodic motion through the first and second types of grazing bifurcations, respectively.Under theb1 or ω change, the transition between 1−p−q and 1−p+1−q periodic motion conforms to the following rules:(10)b1ω↓:1−p+∞−q←GR1⋯1−p+2−q←GR11−p+1−q←GR11−p−q.Under theb1 or ω change, the transition between 1−p−q and 1−p−q+1 periodic motion conforms to the following rules:(11)b1ω↓:1−p−q+∞←GR2⋯1−p−q+2←GR21−p−q+1←GR21−p−q.Selectω=0.1722 and b1∈0.537,1.035 to analyze the transition process between 1−p−q and 1−p+1−q periodic motions in detail. It can be seen from Figure 4(a) that, in the process of large decrease of b1, the first type of grazing bifurcation occurs in the system in turn, which causes the number of collisions p between the vibrator and the elastic restraint surface with a stiffness coefficient K1 to increase one by one. For example, the first type of grazing bifurcation occurs in the system at b1=0.80143, and 1−4−0 periodic motion transitions to 1−5−0 periodic motion. Figure 5(a) is the phase diagram of the grazing periodic motion, and Figure 5(b) is a partially enlarged view of the collision surface.Figure 4
Single-parameter bifurcation diagram of the system. (a)ω=0.1722. (b) ω=0.2822.
(a)(b)Figure 5
Phase diagram of periodic motion. (a)1−4−0 phase diagram of grazing periodic motion. (b) Partially enlarged view of Figure 5(a). (c) 1−6−2 phase diagram of grazing periodic motion. (d) Partially enlarged view of Figure 5(d).
(a)(b)(c)(d)Selectω=0.2822 and b1∈0.1,0.2232 to specifically analyze the transition process between the 1−p−q and 1−p−q+1 periodic motions. It can be seen from Figure 4(b) that when b1 reduces, the second type of grazing bifurcation occurs in the system in turn, which causes the number of collisions q between the vibrator and the elastic restraint surface with a stiffness coefficient K2 to increase one by one. For example, the first type of grazing bifurcation occurs in the system at b1=0.12254, and 1−6−2 periodic motion transits to 1−6−3 periodic motion. Figure 5(c) is a phase diagram of the periodic motion of the grazing, and Figure 5(d) is a partially enlarged view of the collision surface.When adjacent periodic motions transfer to each other, because the parameter points of different types of bifurcation are not in the same position, the transfer process is irreversible, so that a polymorphic coexistence region composed of two adjacent periodic motions will be formed in the process of the transition of the two adjacent periodic motions. Take the system parametersb1=0.5826 and ω=0.6928,0.74 to analyze the polymorphic coexistence region when the 1−1−1 and 1−2−2 periodic motions transit mutually. Figure 6 is the bidirectional bifurcation diagram of the ω change, and the red curve shows the single-parameter bifurcation when ω is greatly reduced, and the black curve represents the single-parameter bifurcation diagram when ω increases gradually; it is convenient to observe the coexistence of different attractors. When ω decreases, the system transits from 1−1−1 periodic motion at ω=0.7092 to 1−2−2 periodic motion via the first type of grazing bifurcation; when ω increases, the system transits from 1−2−2 periodic motion to 1−1−1 via saddle-node bifurcation to 1−2−2 cycles of movement. In this way, an area where 1−1−1 and 1−2−2 periodic motions coexist stably is formed between the two bifurcation points. Combined with the cell mapping method, select the initial state domain: Ω=x,x˙|−2<x<−0.5,0<x˙<2, divide the initial state domain into 1000×1000 state cells, and the distribution diagram of the attractor of ω=0.715 is shown in Figure 7.Figure 6
Bidirectional single-parameter bifurcation diagram ofb1=0.5826.Figure 7
The attraction domain on the phase-fixing surfaceΣn at ω=0.715.


### 4.2. Analysis of the Transition Law with Transition Zone between Adjacent Periodic Motions

In the parameter area shown in Figure3(a), in addition to the above-mentioned basic transition methods between adjacent periodic motions, the transition process between adjacent basic periodic motions will be more complicated due to the emergence of the transition zone.The first type of transition zone shown in Figure3(e) generally exists at the peak position of the boundary line of two adjacent periodic motions. This type of transition zone has self-similar characteristics. The transition process of the first type of transition zone is as follows: with the decrease of b1, the system will have two grazing bifurcations of the first type in sequence, resulting in the basic periodic motion 1−p−q transiting to 2−2p+2−2q periodic motion and then to 1−p+1−q periodic motion via the inverse period-doubling bifurcation. The transitional process can be expressed as(12)b1↓:1−p+1−q←IPB2−2p+2−2q←GR12−2p+1−2q←GR11−p−q.Selectω=0.2706 and b1∈0.5261,0.5889 to analyze the process of 1−p−q periodic motion transiting to 1−p+1−q periodic motion via the first transition zone. It can be seen from Figure 8 that, in the process of b1 decreasing from large, 1−4−0 periodic motion transits to 2−9−0 periodic motion through the first type of grazing bifurcation at b1=0.5782; the 2−9−0 periodic motion transits to the 2−10−0 periodic motion through the first type of grazing bifurcation at b1=0.5467. Figure 9(a) is the phase diagram of the 2−9−0 grazing periodic motion, and Figure 9(b) is the partially enlarged view of the collision surface. Afterwards, the 2−10−0 periodic motion transits to the 1−5−0 periodic motion through the reverse period-doubling bifurcation.Figure 8
Single-parameter bifurcation diagram atω=0.2706.Figure 9
Phase diagram of the periodic motion. (a)2−9−0 phase diagram of grazing periodic motion. (b) Partially enlarged view of Figure 9(a).
(a)(b)The second type of transition zone shown in Figure3(f) mainly exists at the intersection of multiple adjacent periodic motions. Due to the appearance of period-doubling bifurcation, 1−p+1−q and 1−p−q+1 periodic motions will transit to the 2i−2ip+1−2iq and 2i−2ip−2iq+1 periodic motions, respectively; 1−p−q will also transit directly to the chaotic motion via the saddle-node bifurcation. In this way, a transition zone with similar periodic motion distribution will be formed at the intersection of multiple adjacent periodic motions.The parametersb1=0.2501 and w∈0.382,0.396 are selected to specifically analyze the transitional process between various periodic motions in the transition zone. The single-parameter bifurcation diagram of the system and the Lyapunov exponent spectrum of the system are shown in Figures 10 and 11. When ω decreases gradually, the system transits from 1−4−2 periodic motion to chaos through saddle-node bifurcation and then transits to stable periodic motion through a series of inverse period-doubling bifurcations and then to 2−9−0 periodic motion in the same way. After two first-type grazing bifurcations and two inverse period-doubling bifurcations, the system eventually transits to 1−5−0 periodic motion. When ω increases gradually, the system transits from 1−5−0 periodic motion to 2−9−0 periodic motion through period-doubling bifurcation and saddle-node bifurcation and then to stable periodic motion through a series of period-doubling bifurcations and boundary crises and then enters chaos state, followed by transiting to 1−4−2 periodic motion via a series of inverse period-doubling bifurcations. It can be seen from Figure 11 that when ω increases gradually, the system enters the chaos state via saddle-node bifurcation at ω=0.3924. The maximum Lyapunov exponent of the system is shown with blue line; at this point it changes from a negative number to a positive number; when ω decreases gradually, the system transits from chaos motion to 1−4−2 periodic motion at ω=0.3935 via a series of inverse period-doubling bifurcations; the maximum Lyapunov exponent of the system is shown with red line. After a series of swings at the zero-scale mark, the system stabilizes to a negative number, and the system is stable. The specific transitional process during the two-way change of ω is as follows:(13)ω↓:1−5−0⟵IPB2−10−0⟵IPB4−20−0⟵GR14−19−0⟵GR12−9−0⟵IPB⋯⟵IPBchaos⟵SN3−13⟵0⟵IPB⋯⟵IPBchaos⟵SN1−4−0,ω↑:1−5−0⟶PB2−10−0⟶PB4−20−0⟶SN4−19−0⟶SN2−9−0⟶PB⋯⟶PBchaos⟶BC3−13−0⟶PB⋯⟶PBchaos⟶IPB⋯⟶IPB1−4−2.Figure 10
Bidirectional single parameter.Figure 11
Maximum Lyapunov exponent ofb1=0.2501.It can be seen from Figure8 that since the saddle-node bifurcation and the inverse period-doubling bifurcation are not in the same position, there will be a polymorphic coexistence zone where the 1−4−2 periodic motion and the chaotic motion coexist between the two parameter points. Combined with the cell mapping method, select the initial state domain: Ω=x,x˙|−0.3<x<0.3,0<x˙<0.6, divide the initial state domain into 1000×1000 state cells, and the system periodic attractor distribution diagram is shown in Figure 12.Figure 12
The attraction domain on the phase-fixed surfaceΣn. (a) ω=0.39347. (b) ω=0.39357. (c) ω=0.39259
(a)(b)(c)


## 4.1. Analysis of the Basic Transition Law between Adjacent Periodic Motions

In the low-frequency small gap area, as shown in Figure3(b), the 1−p−q periodic motion can be transferred to 1−p+1−q and 1−p−q+1 periodic motion through the first and second types of grazing bifurcations, respectively.Under theb1 or ω change, the transition between 1−p−q and 1−p+1−q periodic motion conforms to the following rules:(10)b1ω↓:1−p+∞−q←GR1⋯1−p+2−q←GR11−p+1−q←GR11−p−q.Under theb1 or ω change, the transition between 1−p−q and 1−p−q+1 periodic motion conforms to the following rules:(11)b1ω↓:1−p−q+∞←GR2⋯1−p−q+2←GR21−p−q+1←GR21−p−q.Selectω=0.1722 and b1∈0.537,1.035 to analyze the transition process between 1−p−q and 1−p+1−q periodic motions in detail. It can be seen from Figure 4(a) that, in the process of large decrease of b1, the first type of grazing bifurcation occurs in the system in turn, which causes the number of collisions p between the vibrator and the elastic restraint surface with a stiffness coefficient K1 to increase one by one. For example, the first type of grazing bifurcation occurs in the system at b1=0.80143, and 1−4−0 periodic motion transitions to 1−5−0 periodic motion. Figure 5(a) is the phase diagram of the grazing periodic motion, and Figure 5(b) is a partially enlarged view of the collision surface.Figure 4
Single-parameter bifurcation diagram of the system. (a)ω=0.1722. (b) ω=0.2822.
(a)(b)Figure 5
Phase diagram of periodic motion. (a)1−4−0 phase diagram of grazing periodic motion. (b) Partially enlarged view of Figure 5(a). (c) 1−6−2 phase diagram of grazing periodic motion. (d) Partially enlarged view of Figure 5(d).
(a)(b)(c)(d)Selectω=0.2822 and b1∈0.1,0.2232 to specifically analyze the transition process between the 1−p−q and 1−p−q+1 periodic motions. It can be seen from Figure 4(b) that when b1 reduces, the second type of grazing bifurcation occurs in the system in turn, which causes the number of collisions q between the vibrator and the elastic restraint surface with a stiffness coefficient K2 to increase one by one. For example, the first type of grazing bifurcation occurs in the system at b1=0.12254, and 1−6−2 periodic motion transits to 1−6−3 periodic motion. Figure 5(c) is a phase diagram of the periodic motion of the grazing, and Figure 5(d) is a partially enlarged view of the collision surface.When adjacent periodic motions transfer to each other, because the parameter points of different types of bifurcation are not in the same position, the transfer process is irreversible, so that a polymorphic coexistence region composed of two adjacent periodic motions will be formed in the process of the transition of the two adjacent periodic motions. Take the system parametersb1=0.5826 and ω=0.6928,0.74 to analyze the polymorphic coexistence region when the 1−1−1 and 1−2−2 periodic motions transit mutually. Figure 6 is the bidirectional bifurcation diagram of the ω change, and the red curve shows the single-parameter bifurcation when ω is greatly reduced, and the black curve represents the single-parameter bifurcation diagram when ω increases gradually; it is convenient to observe the coexistence of different attractors. When ω decreases, the system transits from 1−1−1 periodic motion at ω=0.7092 to 1−2−2 periodic motion via the first type of grazing bifurcation; when ω increases, the system transits from 1−2−2 periodic motion to 1−1−1 via saddle-node bifurcation to 1−2−2 cycles of movement. In this way, an area where 1−1−1 and 1−2−2 periodic motions coexist stably is formed between the two bifurcation points. Combined with the cell mapping method, select the initial state domain: Ω=x,x˙|−2<x<−0.5,0<x˙<2, divide the initial state domain into 1000×1000 state cells, and the distribution diagram of the attractor of ω=0.715 is shown in Figure 7.Figure 6
Bidirectional single-parameter bifurcation diagram ofb1=0.5826.Figure 7
The attraction domain on the phase-fixing surfaceΣn at ω=0.715.


## 4.2. Analysis of the Transition Law with Transition Zone between Adjacent Periodic Motions

In the parameter area shown in Figure3(a), in addition to the above-mentioned basic transition methods between adjacent periodic motions, the transition process between adjacent basic periodic motions will be more complicated due to the emergence of the transition zone.The first type of transition zone shown in Figure3(e) generally exists at the peak position of the boundary line of two adjacent periodic motions. This type of transition zone has self-similar characteristics. The transition process of the first type of transition zone is as follows: with the decrease of b1, the system will have two grazing bifurcations of the first type in sequence, resulting in the basic periodic motion 1−p−q transiting to 2−2p+2−2q periodic motion and then to 1−p+1−q periodic motion via the inverse period-doubling bifurcation. The transitional process can be expressed as(12)b1↓:1−p+1−q←IPB2−2p+2−2q←GR12−2p+1−2q←GR11−p−q.Selectω=0.2706 and b1∈0.5261,0.5889 to analyze the process of 1−p−q periodic motion transiting to 1−p+1−q periodic motion via the first transition zone. It can be seen from Figure 8 that, in the process of b1 decreasing from large, 1−4−0 periodic motion transits to 2−9−0 periodic motion through the first type of grazing bifurcation at b1=0.5782; the 2−9−0 periodic motion transits to the 2−10−0 periodic motion through the first type of grazing bifurcation at b1=0.5467. Figure 9(a) is the phase diagram of the 2−9−0 grazing periodic motion, and Figure 9(b) is the partially enlarged view of the collision surface. Afterwards, the 2−10−0 periodic motion transits to the 1−5−0 periodic motion through the reverse period-doubling bifurcation.Figure 8
Single-parameter bifurcation diagram atω=0.2706.Figure 9
Phase diagram of the periodic motion. (a)2−9−0 phase diagram of grazing periodic motion. (b) Partially enlarged view of Figure 9(a).
(a)(b)The second type of transition zone shown in Figure3(f) mainly exists at the intersection of multiple adjacent periodic motions. Due to the appearance of period-doubling bifurcation, 1−p+1−q and 1−p−q+1 periodic motions will transit to the 2i−2ip+1−2iq and 2i−2ip−2iq+1 periodic motions, respectively; 1−p−q will also transit directly to the chaotic motion via the saddle-node bifurcation. In this way, a transition zone with similar periodic motion distribution will be formed at the intersection of multiple adjacent periodic motions.The parametersb1=0.2501 and w∈0.382,0.396 are selected to specifically analyze the transitional process between various periodic motions in the transition zone. The single-parameter bifurcation diagram of the system and the Lyapunov exponent spectrum of the system are shown in Figures 10 and 11. When ω decreases gradually, the system transits from 1−4−2 periodic motion to chaos through saddle-node bifurcation and then transits to stable periodic motion through a series of inverse period-doubling bifurcations and then to 2−9−0 periodic motion in the same way. After two first-type grazing bifurcations and two inverse period-doubling bifurcations, the system eventually transits to 1−5−0 periodic motion. When ω increases gradually, the system transits from 1−5−0 periodic motion to 2−9−0 periodic motion through period-doubling bifurcation and saddle-node bifurcation and then to stable periodic motion through a series of period-doubling bifurcations and boundary crises and then enters chaos state, followed by transiting to 1−4−2 periodic motion via a series of inverse period-doubling bifurcations. It can be seen from Figure 11 that when ω increases gradually, the system enters the chaos state via saddle-node bifurcation at ω=0.3924. The maximum Lyapunov exponent of the system is shown with blue line; at this point it changes from a negative number to a positive number; when ω decreases gradually, the system transits from chaos motion to 1−4−2 periodic motion at ω=0.3935 via a series of inverse period-doubling bifurcations; the maximum Lyapunov exponent of the system is shown with red line. After a series of swings at the zero-scale mark, the system stabilizes to a negative number, and the system is stable. The specific transitional process during the two-way change of ω is as follows:(13)ω↓:1−5−0⟵IPB2−10−0⟵IPB4−20−0⟵GR14−19−0⟵GR12−9−0⟵IPB⋯⟵IPBchaos⟵SN3−13⟵0⟵IPB⋯⟵IPBchaos⟵SN1−4−0,ω↑:1−5−0⟶PB2−10−0⟶PB4−20−0⟶SN4−19−0⟶SN2−9−0⟶PB⋯⟶PBchaos⟶BC3−13−0⟶PB⋯⟶PBchaos⟶IPB⋯⟶IPB1−4−2.Figure 10
Bidirectional single parameter.Figure 11
Maximum Lyapunov exponent ofb1=0.2501.It can be seen from Figure8 that since the saddle-node bifurcation and the inverse period-doubling bifurcation are not in the same position, there will be a polymorphic coexistence zone where the 1−4−2 periodic motion and the chaotic motion coexist between the two parameter points. Combined with the cell mapping method, select the initial state domain: Ω=x,x˙|−0.3<x<0.3,0<x˙<0.6, divide the initial state domain into 1000×1000 state cells, and the system periodic attractor distribution diagram is shown in Figure 12.Figure 12
The attraction domain on the phase-fixed surfaceΣn. (a) ω=0.39347. (b) ω=0.39357. (c) ω=0.39259
(a)(b)(c)


## 5. The Influence of System Parameters on Periodic Motion Distribution

This section uses the multiparameter collaborative simulation method to focus on the impact of the vibration system parameters, shown in Figure1, on the periodic motion distribution and provides a theoretical basis for the reasonable design of matching parameters. The dynamic characteristics of the dimensionless equation (1) are mainly determined by the six parameters of ζ, μk1, μk2, b1, b2, and ω. Take the system parameters ζ=0.1, μk1=20, μk2=50, ω∈0.1,2, b1∈0.1,2, and b2=0.1 as benchmark parameters. The influence of a certain parameter on the distribution of periodic motion of the system is analyzed and studied, taking discrete values within a reasonable range.


### 5.1. The Influence of Clearanceb2 on the Dynamic Characteristics of the System

The reference parameter of the clearanceb2 is 0.1. For discrete values within a reasonable range, take b2=0.5, b2=1.0, and b2=1.5 to calculate the distribution of periodic motion on the parameter plane ω,b1 under different b2 values. Figures 13(a)–13(c) correspond to periodic motion distribution diagram under different values of b2. By comparison with the reference parameter in Figure 3(a), it can be seen that when the clearance b2 gradually increases, the system periodic motion parameter domain and its distribution law remain roughly unchanged, and the transition mode in the transition zone is basically the same. Comparing Figures 13(a)–13(c) with Figure 3(a), when b2 increases to 0.5, the 1−1−1 periodic motion distribution area is significantly reduced, and it shifts to a larger place of b1, there is no impact between the vibrator and the elastic constraint surface with a stiffness coefficient K2, and the 1−p−0p≥1 periodic motion distribution area is slightly widened. When b2 increases to 1, the 1−1−1 periodic motion distribution area basically disappears and only exists in a part of the area at b1=1.6. When b2 increases to 1.5, the distribution area of the 1−1−1 periodic motion disappears completely, and the vibrator will have no impact on the elastic restraint surface with a stiffness coefficient K2. At this time, the system can be regarded as only a kind of elastic restraint with a stiffness coefficient K1.Figure 13
Distribution of periodic motion withb2 change. (a) b2=0.5. (b) b2=1.0. (c) b2=1.5.
(a)(b)(c)


### 5.2. The Influence of Damping Coefficientζ on the Dynamic Characteristics of the System

The datum parameter of the damping coefficientζ is 0.1. Take discrete values for ζ within a reasonable range. Take ζ=0.05, ζ=0.3, and ζ=0.5, respectively, to calculate the distribution of periodic motion on the ω,b1 parameter. Figures 14(a)–14(c) correspond to the periodic motion distribution diagram under different ζ values. Compared with the reference parameter in Figure 3(a), the change of damping coefficient ζ has a greater impact on the periodic motion distribution area. As shown in Figure 14(a), the peak value of each periodic motion generally increases, and the 1−2−2 periodic motion distribution area is expanded to a larger clearance interval and the reduction of the damping coefficient of ζ causes the impact vibration system to present a very complex and fractal distribution law in the low-frequency range. The number of transition regions has increased significantly, the transition method is more complicated, and flutter and chaos occur in the system in many places. When the damping coefficient is ζ=0.3, the dynamic characteristics of the system are no longer complicated. The transition between adjacent periodic motions is only through the basic transition method, the transition zone disappears, and the distribution area of 1−0−0 periodic motion is slightly expanded. When the damping coefficient is ζ=0.5, the types of basic periodic motion are further reduced, the peak value is generally reduced, and the distribution area of 1−0−0 periodic motion is further expanded.Figure 14
Distribution of periodic motion withζ changes. (a) kc=0.05. (b) kc=0.3. (c) kc=0.5.
(a)(b)(c)


## 5.1. The Influence of Clearanceb2 on the Dynamic Characteristics of the System

The reference parameter of the clearanceb2 is 0.1. For discrete values within a reasonable range, take b2=0.5, b2=1.0, and b2=1.5 to calculate the distribution of periodic motion on the parameter plane ω,b1 under different b2 values. Figures 13(a)–13(c) correspond to periodic motion distribution diagram under different values of b2. By comparison with the reference parameter in Figure 3(a), it can be seen that when the clearance b2 gradually increases, the system periodic motion parameter domain and its distribution law remain roughly unchanged, and the transition mode in the transition zone is basically the same. Comparing Figures 13(a)–13(c) with Figure 3(a), when b2 increases to 0.5, the 1−1−1 periodic motion distribution area is significantly reduced, and it shifts to a larger place of b1, there is no impact between the vibrator and the elastic constraint surface with a stiffness coefficient K2, and the 1−p−0p≥1 periodic motion distribution area is slightly widened. When b2 increases to 1, the 1−1−1 periodic motion distribution area basically disappears and only exists in a part of the area at b1=1.6. When b2 increases to 1.5, the distribution area of the 1−1−1 periodic motion disappears completely, and the vibrator will have no impact on the elastic restraint surface with a stiffness coefficient K2. At this time, the system can be regarded as only a kind of elastic restraint with a stiffness coefficient K1.Figure 13
Distribution of periodic motion withb2 change. (a) b2=0.5. (b) b2=1.0. (c) b2=1.5.
(a)(b)(c)


## 5.2. The Influence of Damping Coefficientζ on the Dynamic Characteristics of the System

The datum parameter of the damping coefficientζ is 0.1. Take discrete values for ζ within a reasonable range. Take ζ=0.05, ζ=0.3, and ζ=0.5, respectively, to calculate the distribution of periodic motion on the ω,b1 parameter. Figures 14(a)–14(c) correspond to the periodic motion distribution diagram under different ζ values. Compared with the reference parameter in Figure 3(a), the change of damping coefficient ζ has a greater impact on the periodic motion distribution area. As shown in Figure 14(a), the peak value of each periodic motion generally increases, and the 1−2−2 periodic motion distribution area is expanded to a larger clearance interval and the reduction of the damping coefficient of ζ causes the impact vibration system to present a very complex and fractal distribution law in the low-frequency range. The number of transition regions has increased significantly, the transition method is more complicated, and flutter and chaos occur in the system in many places. When the damping coefficient is ζ=0.3, the dynamic characteristics of the system are no longer complicated. The transition between adjacent periodic motions is only through the basic transition method, the transition zone disappears, and the distribution area of 1−0−0 periodic motion is slightly expanded. When the damping coefficient is ζ=0.5, the types of basic periodic motion are further reduced, the peak value is generally reduced, and the distribution area of 1−0−0 periodic motion is further expanded.Figure 14
Distribution of periodic motion withζ changes. (a) kc=0.05. (b) kc=0.3. (c) kc=0.5.
(a)(b)(c)


## 6. Conclusions

This paper takes a single-degree-of-freedom vibroimpact system with multilevel elastic constraints as the research object. By constructing the Poincaré map, calculating the Lyapunov exponent spectrum to determine the mode type of the periodic motion of the system, the numerical simulation is carried out under the change of the system parameterω,b1, and, combined with the cell mapping method, the law of transition between adjacent periodic motions of the system is analyzed. The influence of system parameters on the periodic motion distribution is analyzed.(1)
When directly transferring between adjacent periodic motions, in the process ofb1 and ω decreasing, 1−p−q periodic motion can be transferred to 1−p+1−q and 1−p−q+1 periodic motions through the first and second types of grazing bifurcations, respectively. When b1 and ω increase, 1−p+1−q and 1−p−q+1 periodic motions can be transferred to 1−p−q periodic motion through saddle-node bifurcation. When the adjacent periodic motion transits through the transition zone, the 1−p−q periodic motion transits to 2i−2ip−2iq via period-doubling bifurcation, to 2−2p−2q+1 via the grazing bifurcation, or to the chaotic motion via the saddle-node bifurcation. Due to the different positions of different types of bifurcations in the system, the transition process is irreversible. In this way, a polymorphic coexistence zone, composed of two adjacent periodic motions, will be formed during the transition of adjacent periodic motions.(2)
By selecting different parameters and comparing them with benchmark parameters, the influence of each parameter on the distribution of system periodic motion is analyzed. The change of gapb2 will change the periodic motion type of the system. With the increase of b2, the area where the vibrator collides with the elastic restraint surface with a stiffness coefficient K2 gradually decreases. When the collision area is b2=2, the collision area completely disappears. It is equivalent to the fact that there is only one elastic constraint surface with stiffness coefficient K1 in the system. The change of the damping coefficient ζ will have a greater impact on the periodic motion distribution area, and a smaller damping coefficient makes the basic periodic motion extremely complex and fractal at low frequencies, larger damping will cause the system to have only basic periodic motion, and there is no transition zone transition between adjacent periodic motions.A future study in this field may include chaos control, bifurcation control, and coexistence attractor transition control of the systems with multiple parameters and multiple constraints. The study on the global dynamic characteristics of high-dimensional vibroimpact system is of great importance to the solution of practical engineering problems.



---

*Source: 6687887-2021-02-09.xml*