# Double Grazing Periodic Motions and Bifurcations in a Vibroimpact System with Bilateral Stops

Authors: Qunhong Li; Limei Wei; Jieyan Tan; Jiezhen Xi
Journal: Abstract and Applied Analysis
 (2014)

---

## Abstract

The double grazing periodic motions and bifurcations are investigated for a two-degree-of-freedom vibroimpact system with symmetrical rigid stops in this paper. From the initial condition and periodicity, existence of the double grazing periodic motion of the system is discussed. Using the existence condition derived, a set of parameter values is found that generates a double grazing periodic motion in the considered system. By extending the discontinuity mapping of one constraint surface to that of two constraint surfaces, the Poincaré map of the vibroimpact system is constructed in the proximity of the grazing point of a double grazing periodic orbit, which has a more complex form than that of the single grazing periodic orbit. The grazing bifurcation of the system is analyzed through the Poincaré map with clearance as a bifurcation parameter. Numerical simulations show that there is a continuous transition from the chaotic band to a period-1 periodic motion, which is confirmed by the numerical simulation of the original system.


---

## Body

## 1. Introduction

Nonsmooth dynamical systems have some special bifurcations, such as grazing bifurcation, sliding bifurcation, and corner-collision bifurcation, besides the bifurcations occurring in the smooth dynamical systems. These bifurcations exhibit the complicated behavior of dynamical systems [1]. The research of nonsmooth dynamical systems is more difficult than that of smooth systems, so it attracts much attention of scholars from all over the world.In the 1980s, Shaw and his coworkers [2, 3] studied the impact oscillators with the theory of modern dynamical systems. They considered the motion of a single-degree-of-freedom periodically forced oscillator subjected to a rigid amplitude constraint and found that grazing impact leads to the singularity of the Poincaré map, which makes a great effect on global dynamical behavior of systems. Nordmark [4] considered a single-degree-of-freedom periodically forced oscillator subjected to a rigid amplitude constraint. By the analytical methods, a nonconventional bifurcation caused by grazing impact, that is, grazing bifurcation is studied. Chin et al. [5, 6] investigated in detail the dynamics of a vibroimpact system near the grazing impact by the Nordmark map deduced in [4]. Some phenomena that appeared only in nonsmooth systems were found, for example, the bifurcation from periodic motion into chaos and period-adding bifurcation scenario and so forth. Ivanov [7] studied how to obtain a stable periodic impact motion from the nonimpact periodic motion. Virgin and Begley [8] explored some interesting global dynamic behavior in the response of a double-sided, harmonically forced, impact oscillator including the influence of Coulomb damping. Both basins of attraction and grazing bifurcations were studied. Using discontinuity mapping, Fredriksson and Nordmark [9, 10] developed a normal-form calculation for nonsmooth systems with several degrees of freedom which are useful for analyzing the dynamics close to bifurcations. Bernardo et al. [11] performed local analysis of grazing bifurcations in n-dimensional piecewise-smooth systems of ordinary differential equations. Under quite general circumstances, they showed that this leads to a normal-form mapping containing in lowest order either a square-root or a (3/2)-type singularity. Bernardo et al. [12] analyzed sliding bifurcations in n-dimensional piecewise-smooth dynamical systems with discontinuous vector field and derived the normal form of the discontinuity mapping. Foale and Bishop [13] and Hu [14] also studied the grazing bifurcations of impact oscillators and found that it is the grazing bifurcation that leads to the change of stability of periodic motion and, as a result, gives rise to chaos. Dankowicz et al. [15–17] investigated codimension-one and codimension-two grazing bifurcation in impact microactuators and the single-degree-of-freedom impact oscillators, respectively. Unfoldings of the degenerate grazing dynamics in impact microactuators were analyzed. Kowalczyk et al. [18] proposed a strategy for the classification of codimension-two discontinuity-induced bifurcations of limit cycles in piecewise-smooth systems. Li and Tan [19] presented a method for Lyapunov exponent calculation of a two-degree-of-freedom vibroimpact system with symmetrical rigid stops, which can be used for chaotic motions of the system. Shen et al. [20] analyzed subharmonic and grazing bifurcations for a simple bilinear oscillator via a combination of analytical and numerical methods. They found that the dynamics of the system for the case of large dissipation is quite different from that for the case of small dissipation.In this paper, existence of double grazing period-n orbit is derived analytically for a two-degree-of-freedom vibroimpact system with symmetrical rigid stops, while most of the literature mainly focused on the grazing motions with only unilateral constraint in the past decades. Adopting the idea of discontinuity mapping, the Poincaré map is constructed in the vicinity of the grazing point for the system with two-sided constraints, which is more complicated than that of unilateral constraint. The grazing bifurcation of periodic orbit is explored with the map derived and a continuous transition from the chaotic band to a period-1 periodic motion is found.The rest of the paper is outlined as follows. Section2 describes the mechanical model and gives the equations of motion. The existence of grazing periodic motion is discussed by the analytical method in Section 3. In Section 4, the discontinuity mapping near the double grazing periodic motion is deduced and the bifurcations are performed by the numerical simulations. Finally, we give a brief conclusion in Section 5.


## 2. The Mechanical Model

A two-degree-of-freedom system having symmetrically placed rigid stops and subjected to periodic excitation is shown in Figure1 [19, 21]. Displacements of the masses M
1 and M
2 are represented by X
1 and X
2, respectively. The masses are connected to linear springs with stiffnesses K
1 and K
2 and linear viscous dashpots with damping constants C
1 and C
2. The masses move only in the horizontal direction and the excitations on both masses are harmonic, which take the forms of P
i
sin
(
Ω
T
+
τ
)
(
i
=
1,2
), respectively. M
1 moves between the rigid stops A and C. When the displacement X
1 of the mass M
1 is B(or -
B), the mass M
1 will hit the rigid stop A (or C), and the velocity of the mass M
1 will change its value and direction. Then M
1 hits the stop C (or A) again under some condition and so on.Figure 1
Schematic representation of a two-degree-of-freedom vibroimpact system with symmetrical rigid stops.Damping in the mechanical model is assumed as proportional damping of the Rayleigh type, which in this case implies thatC
1
/
K
1
=
C
2
/
K
2. The impact is described by a coefficient of restitution R, and the duration of impact is negligible compared to the period of the force.Between any two consecutive impacts(
|
X
1
|
<
B
), the differential equations of motion are given by
(1)
[
M
1
0
0
M
2
]
{
X
¨
1
X
¨
2
}
+
[
C
1
-
C
1
-
C
1
C
1
+
C
2
]
{
X
˙
1
X
˙
2
}
+
[
K
1
-
K
1
-
K
1
K
1
+
K
2
]
{
X
1
X
2
}
=
[
P
1
P
2
]
sin
(
Ω
T
+
τ
)
,
|
X
1
|
<
B
.
The impacting equations of mass M
1 are as follows:
(2)
X
˙
1
A
+
=
-
R
X
˙
1
A
-
(
X
1
=
B
)
,
X
˙
1
C
+
=
-
R
X
˙
1
C
-
(
X
1
=
-
B
)
,
where the subscripts -
,
+ denote the values just before and after impact, respectively.Introduce the nondimensional quantitiesμ
m
=
M
2
/
M
1, μ
k
=
K
2
/
K
1, μ
k
=
μ
c, f
2
=
P
2
/
(
P
1
+
P
2
), ω
=
Ω
M
1
/
K
1, t
=
T
K
1
/
M
1, ζ
=
C
1
/
2
K
1
M
1,  b
=
B
K
1
/
(
P
1
+
P
2
), and x
i
=
X
i
K
1
/
(
P
1
+
P
2
). Thus nondimensional differential equations of motion without impacting will now have the form
(3)
[
1
0
0
μ
m
]
{
x
¨
1
x
¨
2
}
+
[
2
ζ
-
2
ζ
-
2
ζ
2
ζ
(
1
+
μ
c
)
]
{
x
˙
1
x
˙
2
}
+
[
1
-
1
-
1
1
+
μ
k
]
{
x
1
x
2
}
=
{
1
-
f
2
f
2
}
sin
(
ω
t
+
τ
)
,
|
x
1
|
<
b
,
and the impacting equations of mass M
1 become
(4)
x
˙
1
A
+
=
-
R
x
˙
1
A
-
,
(
x
1
=
b
)
,
x
˙
1
C
+
=
-
R
x
˙
1
C
-
,
(
x
1
=
-
b
)
.LetΨ stand for the canonical model matrix of (3) and ω
1 and ω
2 denote the eigenfrequencies of the system as no impact occurs. By taking Ψ as a transition matrix, under the transformation X
=
Ψ
ξ, (3) becomes
(5)
I
ξ
¨
+
C
ξ
˙
+
Λ
ξ
=
F
-
sin
(
ω
t
+
τ
)
,
where X
=
(
x
1
,
x
2
)
T, ξ
=
(
ξ
1
,
ξ
2
)
T, I is a unit matrix of degree 2
×
2, and C and Λ are diagonal matrices with C
=
diag
⁡
[
2
ζ
ω
1
2
,
2
ζ
ω
2
2
] and Λ
=
diag
⁡
[
ω
1
2
,
ω
2
2
]. In addition, F
-
=
(
f
-
1
,
f
-
2
)
T
=
Ψ
T
P
k and P
k
=
(
1
-
f
2
,
f
2
)
T. Equation (3) can be solved by using the modal coordinates and the modal matrix approach. The general solution is
(6)
x
i
(
t
)
=
∑
j
=
1
2
ψ
i
j
(
e
-
η
j
(
t
-
t
0
)
(
a
j
cos
⁡
ω
d
j
(
t
-
t
0
)
+
b
j
sin
ω
d
j
(
t
-
t
0
)
)
+
A
j
sin
(
ω
t
+
τ
)
+
B
j
cos
⁡
(
ω
t
+
τ
)
)
,
where t
0 denotes the time when the mass M
1 collides with the constraint A or C, ψ
i
j
(
i
,
j
=
1,2
) are the elements of the canonical modal matrix Ψ, η
j
=
ζ
ω
j
2, ω
d
j
=
ω
j
2
-
η
j
2, and a
j,  b
j are the constants of integration which are determined by the initial condition and modal parameters of the system. A
j,  B
j are the amplitude parameters given by
(7)
A
j
=
1
2
ω
d
j
(
ω
+
ω
d
j
(
ω
+
ω
d
j
)
2
+
η
j
2
-
ω
-
ω
d
j
(
ω
-
ω
d
j
)
2
+
η
j
2
)
f
-
j
,
B
j
=
η
j
2
ω
d
j
(
1
(
ω
-
ω
d
j
)
2
+
η
j
2
-
1
(
ω
+
ω
d
j
)
2
+
η
j
2
)
f
-
j
.Letx
˙
1
=
v
1, x
˙
2
=
v
2, and φ
=
ω
t
(
mod
⁡
2
π
); then (5) can be changed into one-order autonomous dynamical system (8)
(
x
˙
1
v
˙
1
x
˙
2
v
˙
2
φ
˙
)
=
(
v
1
-
2
ζ
v
1
+
2
ζ
v
2
-
x
1
+
x
2
+
(
1
-
f
2
)
sin
(
φ
+
τ
)
v
2
2
ζ
v
1
-
2
ζ
(
1
+
μ
c
)
v
2
+
x
1
-
(
1
+
μ
k
)
x
2
+
f
2
sin
(
φ
+
τ
)
μ
m
ω
)
,where the analytical expressions of x
i and v
i
(
i
=
1,2
) can be obtained from (6).


## 3. Existence of Grazing Periodic Motion

If the oscillatorM
1 impacts each rigid constraint with zero velocity and the direction of the acceleration is opposite to the motion, then we say that the system is undergoing grazing motion. A grazing periodic motion may be denoted by n
-
p which means that the oscillator M
1 grazes with each constraint for p times in n periodic external excitation force. In the following, we will derive an existence condition of grazing motion with period n
T, where T is the period of external excitation. Assume that the grazing periodic motion begins from the grazing point on the constraint A. Inserting the initial conditions
(9)
x
1
(
0
)
=
x
10
=
b
,
x
2
(
0
)
=
x
20
,
v
1
(
0
)
=
v
10
=
0
,
v
2
(
0
)
=
v
20
and the periodic conditions
(10)
x
1
(
2
n
π
ω
)
=
x
10
=
b
,
x
2
(
2
n
π
ω
)
=
x
20
,
v
1
(
2
n
π
ω
)
=
v
10
=
0
,
v
2
(
2
n
π
ω
)
=
v
20
into the general solution of system (3), we have
(11)
a
j
=
e
j
b
j
s
j
1
-
e
j
c
j
,
where e
j
=
e
-
η
j
(
2
n
π
/
ω
), c
j
=
cos
⁡
(
ω
d
j
(
2
n
π
/
ω
)
), and s
j
=
sin
(
ω
d
j
(
2
n
π
/
ω
)
), j
=
1
,
2.If the grazing periodic motion begins from the grazing point on the constraintC, similar to the case for (9) and (10), the initial conditions and the periodic conditions are
(12)
x
1
(
0
)
=
x
10
=
-
b
,
x
2
(
0
)
=
x
20
,
v
1
(
0
)
=
v
10
=
0
,
v
2
(
0
)
=
v
20
,
(13)
x
1
(
2
n
π
ω
)
=
x
10
=
-
b
,
x
2
(
2
n
π
ω
)
=
x
20
,
v
1
(
2
n
π
ω
)
=
v
10
=
0
,
v
2
(
2
n
π
ω
)
=
v
20
,
respectively. Inserting (12) and (13) into (6), we can also obtain the expression of a
j, which is the same as (11). Substituting (11) into (9) and (10) yields
(14)
ψ
11
l
-
1
b
1
+
ψ
12
l
-
2
b
2
=
0
,
ψ
21
l
-
1
b
1
+
ψ
22
l
-
2
b
2
=
0
,
where l
-
j
=
ω
d
j
(
1
-
(
e
j
(
c
j
-
e
j
)
/
(
1
-
e
j
c
j
)
)
), j
=
1
,
2. Thus, if
(15)
|
ψ
11
l
-
1
ψ
12
l
-
2
ψ
21
l
-
1
ψ
22
l
-
2
|
=
(
ψ
11
ψ
22
-
ψ
12
ψ
21
)
l
-
1
l
-
2
≠
0
,
we have b
1
=
0 and b
2
=
0. Hence, a
1
=
0 and a
2
=
0.For simplicity, assume that the parameters are chosen such that the integral constantsa
1, a
2, b
1, and b
2 are vanishing. Inserting a
1, a
2, b
1, and b
2 into (6) gives
(16)
d
1
sin
τ
+
d
2
cos
⁡
τ
=
b
,
d
1
cos
⁡
τ
-
d
2
sin
τ
=
0
,
as the grazing periodic motion sets off from the grazing point on the constraint A or
(17)
d
1
sin
τ
+
d
2
cos
⁡
τ
=
-
b
,
d
1
cos
⁡
τ
-
d
2
sin
τ
=
0
,
as the grazing periodic motion leaves from the grazing point on the constraint C. Then it follows that
(18)
b
2
=
d
1
2
+
d
2
2
,
τ
=
arctan
(
d
1
d
2
)
,
where d
1
=
ψ
11
A
1
+
ψ
12
A
2 and d
2
=
ψ
11
B
1
+
ψ
12
B
2. Denote the acceleration of the oscillator M
1 as A
g
1 (or A
g
2) for the case in which the periodic grazing motion begins from the grazing point on the constraint A (or C) with A
g
1
=
x
¨
(
t
)
|
t
=
0
=
-
ω
2
b and A
g
2
=
x
¨
(
t
)
|
t
=
0
=
ω
2
b. If b
>
0, we have A
g
1
<
0 and A
g
2
>
0.Based on the analysis above, we have the conclusion as follows.Theorem 1.
If there exists a double grazing periodic orbit in the system (3)-(4) with initial conditions (9), (12) and periodic conditions (10), (13), then system parameters must satisfy the following condition:
(19)
(
ψ
11
ψ
22
-
ψ
12
ψ
21
)
l
-
1
l
-
2
≠
0
,
b
=
d
1
2
+
d
2
2
,
τ
=
arctan
(
d
1
d
2
)
.In order to verify the existence condition obtained, numerical simulation of the original system will be given in the following.When the parameters have the valuesR
=
0.8, μ
k
=
5, μ
m
=
10, ζ
=
0.05, f
2
=
0, ω
=
0.63, τ
=
2.229311, and b
=
6.7107356, it can be verified that the existence condition is satisfied and the period-1 grazing motion occurs. Phase portrait and time history of the oscillator M
1 are shown in Figures 2 and 3, respectively. Figure 2 shows that the oscillator M
1 collides with the constraints A and C with zero velocity. Figure 4 is the phase portrait of the oscillator M
1 in (
x
1
,
x
2
) plane.Figure 2
Phase portrait of the oscillatorM
1 in (
x
1
,
v
1
) plane.Figure 3
Time history of the displacement of the oscillatorM
1.Figure 4
Phase portrait of the oscillatorM
1 in (
x
1
,
x
2
) plane.


## 4. Grazing Bifurcation

The grazing periodic motion is likely to change qualitatively for a tiny variation in system parameters; that is, grazing bifurcation will take place. As a result, impact periodic motion or nonimpact periodic motion will occur. Obviously, the method for impact periodic motion cannot be applied directly to grazing periodic motion. By extending the discontinuity mapping introduced in [4, 9], we derive a new discontinuity mapping near the double grazing orbit. And then we will analyze the grazing bifurcation with the help of the discontinuity mapping obtained.


### 4.1. Derivation of Discontinuity Mapping of Double Grazing Orbit

Assume that there is a double grazing orbit, denoted asγ, in the system (8). And the grazing point of constraint A is marked as X
1
* with X
1
*
=
(
x
1
A
*
,
v
1
A
*
,
x
2
A
*
,
v
2
A
*
). Define a function h
1
(
X
) in a neighborhood B
1
(
X
1
*
) of X
1
* by
(20)
h
1
(
X
)
=
x
1
-
b
.
Therefore, h
1
(
X
)
<
0 implies that the mass M
1 has no impact with the constraint plane A, while h
1
(
X
)
>
0 corresponds to a penetration of the constraint plane A. Thus the impact surface Σ
1 of constraint A can be defined as
(21)
Σ
1
=
{
X
∈
B
1
(
X
1
*
)
∣
h
1
(
X
)
=
0
}
.
A new function H
1
(
X
,
t
)
=
h
1
(
Φ
t
(
X
)
) may be obtained by inserting the solution of the system, X
=
Φ
t
(
X
), into h
1
(
X
), which describes the minimum distance to the constraint A after a flight of time t from a starting point X.According to the relationship between the orbit and the constraintC, we can define some function similarly as before. The grazing point of constraint C is written as X
2
* with X
2
*
=
(
x
1
C
*
,
v
1
C
*
,
x
2
C
*
,
v
2
C
*
). Define a function in a neighbourhood B
2
(
X
2
*
) of X
2
* as follows:
(22)
h
2
(
X
)
=
x
1
+
b
.
Hence, h
2
(
X
)
>
0 implicates that the mass M
1 has no impact with the constraint plane C, while h
2
(
X
)
<
0 corresponds to a penetration of the constraint plane C. The notation Σ
2 is used to denote the impact surface of constraint C; namely,
(23)
Σ
2
=
{
X
∈
B
2
(
X
2
*
)
∣
h
2
(
X
)
=
0
}
.
Another new function H
2
(
X
,
t
)
=
h
2
(
Φ
t
(
X
)
) can be obtained if we insert the solution of the system X
=
Φ
t
(
X
) into h
2
(
X
), which indicates the distance to the constraint C after a flight of time t with X as a starting point.For a double grazing periodic motion, we have(24)
H
1
(
X
1
*
,
0
)
=
0
,
∂
H
1
∂
t
(
X
1
*
,
0
)
=
0
,
A
g
1
=
∂
2
H
1
∂
t
2
(
X
1
*
,
0
)
<
0
,
H
2
(
X
2
*
,
0
)
=
0
,
∂
H
2
∂
t
(
X
2
*
,
0
)
=
0
,
A
g
2
=
∂
2
H
2
∂
t
2
(
X
2
*
,
0
)
>
0
.
If X is restricted on the impact surfaces Σ
1 and Σ
2, respectively, then from the signs of ν
1
=
(
∂
H
1
/
∂
t
)
(
X
,
0
) and ν
2
=
(
∂
H
2
/
∂
t
)
(
X
,
0
), the impact surfaces can be divided into the following subsets:
(25)
Σ
1
+
=
{
X
∈
Σ
1
∣
ν
1
>
0
}
,
Σ
1
0
=
{
X
∈
Σ
1
∣
ν
1
=
0
}
,
Σ
1
-
=
{
X
∈
Σ
1
∣
ν
1
<
0
}
,
Σ
2
+
=
{
X
∈
Σ
2
∣
ν
2
>
0
}
,
Σ
2
0
=
{
X
∈
Σ
2
∣
ν
2
=
0
}
,
Σ
2
-
=
{
X
∈
Σ
2
∣
ν
2
<
0
}
.
So with the definitions above, the impact process on the impact surfaces Σ
1 and Σ
2 can be expressed as two mappings; that is,
(26)
G
1
:
Σ
1
+
⟶
Σ
1
-
,
G
2
:
Σ
2
-
⟶
Σ
2
+
.Select constant phase planes just before the impact as the Poincaré sections; namely,(27)
Π
1
=
{
(
X
,
θ
)
∣
θ
=
θ
1
}
,
Π
2
=
{
(
X
,
θ
)
∣
θ
=
θ
2
}
,
where θ
1 and θ
2 are the phase angles of the system just before the oscillator M
1 contacts with the constraints A and C, respectively. Let X
-
1
=
γ
∩
Π
1 and X
-
2
=
γ
∩
Π
2. Thus based on the assumptions of periodicity and symmetry of the system, we have Φ
T
/
2
(
X
-
1
)
=
X
-
2 and Φ
T
/
2
(
X
-
2
)
=
X
-
1, where Φ
T
/
2
(
X
) is the smooth flow mapping of the system with T
=
2
n
π
/
ω. As in Figure 5, let X
1
* be the image of X
-
1 after the time of flight t
1, X
-
2 the image of X
1
* after the time of flight t
2, X
2
* the image of X
-
2 after the time of flight t
3, and X
-
1 the image of X
2
* after the time of flight t
4. Then the total time in a period is T
=
t
1
+
t
2
+
t
3
+
t
4. If the points in B
1
(
X
-
1
)
∩
Π
1 are mapped by the flow mapping Φ
t
1
(
X
) without consideration of impact effect, the image of an impacting point would penetrate the surface Σ
1, which obviously leads to an error. This problem would happen for points inB
2
(
X
-
2
)
∩
Π
2 if impact effect is ignored. To deal with the impact incident, two mappings C
1 and C
2 must be constructed for points near X
1
* and X
2
*, which will enable us to write a valid Poincaré map. Consider
(28)
P
(
X
)
=
Φ
t
4
∘
C
2
∘
Φ
t
3
∘
Φ
t
2
∘
C
1
∘
Φ
t
1
(
X
)
.Figure 5
The schematic figure of the discontinuity mappingsC
1 and C
2.In what follows, the mappingsC
1 and C
2 will be derived. Starting at a point X
1 in the neighborhood B
1
(
X
1
*
) of X
1
*, H
1
(
X
1
,
t
) indicates the distance from the oscillator M
1 to the constraint A as the variation of time t. Since X
1
∈
B
1
(
X
1
*
), H
1
(
X
1
,
t
) will reach a local maximum through a small time interval τ
1
=
τ
1
(
X
1
). Introduce a function ψ
1
(
X
1
)
=
H
1
(
X
1
,
τ
1
(
X
1
)
) such that ψ
1
(
X
1
)
≥
0 for impacting points on the constraint A and ψ
1
(
X
1
)
<
0 for nonimpacting points. Similarly, the distance function H
2
(
X
2
,
t
) will get to a local minimum through a small time interval τ
2
=
τ
2
(
X
2
), where X
2 is a starting point in a neighborhood B
2
(
X
2
*
) of X
2
*. A function ψ
2
(
X
2
)
=
H
2
(
X
2
,
τ
2
(
X
2
)
) is introduced such that ψ
2
(
X
2
)
≤
0 for impacting points, while ψ
2
(
X
2
)
>
0 for nonimpacting points.For impacting points on constraint planeA, we define a mapping Φ
t
c
1, where t
c
1 is the time of flight from the impacting point to the crossing point of Σ
1
+ along the flow. Hence, the mapping Φ
t
c
1 takes an impacting point near the constraint A to the impact surface and then under the impact mapping G
1 a point is obtained on Σ
1
- after rigid impact. Since the mapping G
1
∘
Φ
t
c
1 does not take place in zero time [10], it must be composited with Φ
-
t
c
1to reach this point. So for X
1
∈
B
1
(
X
1
*
), the full discontinuity mapping C
1 may be written as follows:
(29)
C
1
=
{
Φ
-
t
c
1
∘
G
1
∘
Φ
t
c
1
ψ
1
≥
0
,
I
ψ
1
<
0
,
where G
1
(
X
1
)
=
(
x
1
A
+
,
v
1
A
+
,
x
2
A
+
,
v
2
A
+
)
T
=
(
x
1
A
-
,
-
R
v
1
A
-
,
x
2
A
-
,
v
2
A
-
)
T. The full discontinuity mapping C
2 can be written similarly as follows:
(30)
C
2
=
{
Φ
-
t
c
2
∘
G
2
∘
Φ
t
c
2
ψ
2
≤
0
,
I
ψ
2
>
0
,
where I is the identity map, t
c
2 is the time of flight from the impacting point of constraint C to crossing point of Σ
2
-, and G
2
(
X
2
)
=
(
x
1
C
+
,
v
1
C
+
,
x
2
C
+
,
v
2
C
+
)
T
=
(
x
1
C
-
,
-
R
v
1
C
-
,
x
2
C
-
,
v
2
C
-
)
T. We now consider the series expansions of the discontinuity mappings C
1 and C
2.Forψ
1
(
X
1
)
≥
0, an impacting point of constraint A, denoted as X
1, can be mapped to the impact surface Σ
1 along the flow. For simplicity, we will still use the notation X
1 to denote the image of the impacting point, which is on Σ
1. By the impact law, it follows that
(31)
Δ
X
1
=
G
1
(
X
1
)
-
X
1
=
(
x
1
A
+
v
1
A
+
x
2
A
+
v
2
A
+
)
-
(
x
1
A
-
v
1
A
-
x
2
A
-
v
2
A
-
)
=
-
v
1
A
-
(
0
R
+
1
0
0
)
≜
-
v
1
A
-
δ
,
where δ
=
(
0
,
R
+
1,0
,
0
)
T.To determinet
c
1, we start by expanding the function H
1 with respect to t. For X
1
∈
B
1
(
X
-
1
), near the time τ
1
(
X
1
), H
1
(
X
1
,
t
) has the form
(32)
H
1
(
X
1
,
t
)
=
ψ
1
(
X
1
)
+
(
t
-
τ
1
(
X
1
)
)
2
A
g
1
2
+
o
(
ψ
1
(
X
1
)
,
t
)
.
From (8), it is clear that A
g
1
=
(
∂
2
H
1
/
∂
t
2
)
(
X
1
*
,
0
)
=
-
x
1
A
*
+
x
2
A
*
-
2
ζ
v
1
A
*
+
2
ζ
v
2
A
*
+
(
1
-
f
2
)
sin
τ. Since H
1
(
X
1
,
t
c
1
)
=
0, we have
(33)
t
c
1
=
τ
1
(
X
1
)
-
-
2
ψ
1
(
X
1
)
A
g
1
+
o
(
X
1
)
.
Expressing ν
1 as a function of X
1 yields
(34)
ν
1
=
∂
H
1
∂
t
(
X
1
,
t
c
1
)
=
-
2
A
g
1
ψ
1
(
X
1
)
+
o
(
X
1
,
t
)
.
Hence,
(35)
Δ
X
1
=
-
v
1
A
-
δ
=
-
ν
1
δ
=
-
-
2
A
g
1
ψ
1
(
X
1
)
δ
+
o
(
X
1
)
.
For X
1
∈
B
1
(
X
1
*
),
(36)
Φ
t
c
1
(
X
1
)
=
X
1
+
F
(
X
1
*
)
t
c
1
+
o
(
X
1
,
t
c
1
)
;
then
(37)
Φ
-
t
c
1
(
X
1
)
∘
G
1
∘
Φ
t
c
1
(
X
1
)
=
X
1
-
-
2
A
g
1
ψ
1
(
X
1
)
δ
+
o
(
X
1
)
.
Finally, we obtain
(38)
C
1
(
X
1
)
=
{
X
1
-
ψ
1
(
X
1
)
-
2
A
g
1
δ
+
o
(
X
1
2
)
,
ψ
1
≥
0
,
X
1
,
ψ
1
<
0
.
The series expression of C
2 can be given in the same way; that is,
(39)
C
2
(
X
2
)
=
{
X
2
+
-
ψ
2
(
X
2
)
2
A
g
2
δ
+
o
(
X
2
2
)
,
ψ
2
≤
0
,
X
2
,
ψ
2
>
0
.
Therefore, the Poincaré map can be written in the following form:
(40)
P
=
Φ
t
4
∘
C
2
∘
Φ
t
3
∘
Φ
t
2
∘
C
1
∘
Φ
t
1
.
Let t
1
→
0 and t
3
→
0; the Poincaré map can be simplified as
(41)
P
=
Φ
T
/
2
∘
C
2
∘
Φ
T
/
2
∘
C
1
.
Let L
1
=
(
∂
Φ
T
/
2
/
∂
X
)
(
X
1
*
), L
2
=
(
∂
Φ
T
/
2
/
∂
X
)
(
X
2
*
), η
1
=
(
∂
h
1
/
∂
X
)
(
X
1
*
), and η
2
=
(
∂
h
2
/
∂
X
)
(
X
2
*
). Thus Φ
T
/
2, ψ
1, ψ
2 can be linearized as the expressions of L
1, L
2, η
1, and η
2. Here L
1 and L
2 are both 4
×
4 matrices, and the elements can be computed by the chain rule of the composite function. For X
1
∈
Π
1, we will discuss the expression of the Poincaré map in the following cases.(
1
)  For the case ψ
1
(
X
1
)
<
0, namely, η
1
(
X
1
-
X
1
*
)
<
0, it yields
(42)
X
2
=
Φ
T
/
2
∘
C
1
(
X
1
)
=
X
2
*
+
L
1
(
C
1
(
X
1
)
-
X
1
*
)
+
o
(
X
1
2
)
=
X
2
*
+
L
1
(
X
1
-
X
1
*
)
+
o
(
X
1
2
)
.
(i)
Ifψ
2
(
X
2
)
>
0, namely, η
2
(
X
2
-
X
2
*
)
>
0, then
(43)
Φ
T
/
2
∘
C
2
∘
Φ
T
/
2
∘
C
1
(
X
1
)
=
Φ
T
/
2
∘
C
2
(
X
2
)
=
X
1
*
+
L
2
(
X
2
-
X
2
*
)
+
o
(
X
2
2
)
.
(ii)
Ifψ
2
(
X
2
)
≤
0, namely, η
2
(
X
2
-
X
2
*
)
≤
0, then
(44)
Φ
T
/
2
∘
C
2
∘
Φ
T
/
2
∘
C
1
(
X
1
)
=
Φ
T
/
2
∘
C
2
(
X
2
)
=
X
1
*
+
L
2
(
C
2
(
X
2
)
-
X
2
*
)
+
o
(
X
2
2
)
=
X
1
*
+
L
2
(
X
2
+
-
2
A
g
2
ψ
2
(
X
2
)
δ
-
X
2
*
)
+
o
(
X
2
2
)
.(
2
) For the case ψ
1
(
X
1
)
≥
0, namely, η
1
(
X
1
-
X
1
*
)
≥
0, it follows that
(45)
X
2
=
Φ
T
/
2
∘
C
1
(
X
1
)
=
X
2
*
+
L
1
(
C
1
(
X
1
)
-
X
1
*
)
+
o
(
X
1
2
)
=
X
2
*
+
L
1
(
X
1
-
-
2
A
g
1
ψ
1
(
X
1
)
δ
-
X
1
*
)
+
o
(
X
1
2
)
.
(i)
Ifψ
2
(
X
2
)
>
0, namely, η
2
(
X
2
-
X
2
*
)
>
0, then
(46)
Φ
T
/
2
∘
C
2
∘
Φ
T
/
2
∘
C
1
(
X
1
)
=
Φ
T
/
2
∘
C
2
(
X
2
)
=
X
1
*
+
L
2
(
X
2
-
X
2
*
)
+
o
(
X
2
2
)
.
(ii)
Ifψ
2
(
X
2
)
≤
0, namely, η
2
(
X
2
-
X
2
*
)
≤
0, then
(47)
Φ
T
/
2
∘
C
2
∘
Φ
T
/
2
∘
C
1
(
X
1
)
=
Φ
T
/
2
∘
C
2
(
X
2
)
=
X
1
*
+
L
2
(
C
2
(
X
2
)
-
X
2
*
)
+
o
(
X
2
2
)
=
X
1
*
+
L
2
(
X
2
+
-
2
A
g
2
ψ
2
(
X
2
)
δ
-
X
2
*
)
+
o
(
X
2
2
)
.Combining expressions (43)-(44) and (46)-(47), for arbitrary point X
1
∈
B
1
(
X
1
*
), we have (48)
P
=
Φ
T
/
2
∘
C
2
∘
Φ
T
/
2
∘
C
1
=
{
X
1
*
+
L
2
(
X
2
-
X
2
*
)
+
o
(
X
2
2
)
,
η
1
(
X
1
-
X
1
*
)
<
0
,
η
2
(
X
2
-
X
2
*
)
>
0
;
X
1
*
+
L
2
(
X
2
+
-
2
A
g
2
ψ
2
(
X
2
)
δ
-
X
2
*
)
+
o
(
X
2
2
)
,
η
1
(
X
1
-
X
1
*
)
<
0
,
η
2
(
X
2
-
X
2
*
)
≤
0
;
X
1
*
+
L
2
(
X
2
-
X
2
*
)
+
o
(
X
2
2
)
,
η
1
(
X
1
-
X
1
*
)
≥
0
,
η
2
(
X
2
-
X
2
*
)
>
0
;
X
1
*
+
L
2
(
X
2
+
-
2
A
g
2
ψ
2
(
X
2
)
δ
-
X
2
*
)
+
o
(
X
2
2
)
,
η
1
(
X
1
-
X
1
*
)
≥
0
,
η
2
(
X
2
-
X
2
*
)
≤
0
.


### 4.2. Numerical Simulation

In this subsection, we will use a set of parameter values that generates the grazing periodic motion as given in Section3 to investigate the grazing bifurcation of the system.TakingR
=
0.8, μ
k
=
5, μ
m
=
10, ζ
=
0.05, f
2
=
0, ω
=
0.63, and τ
=
2.229311, let the parameter b gradually increase from 6.5 to 6.9 by a step of 0.0002. For each b, we take 1500 iterations and the first 1300 of which are omitted as the transient process. Figure 6 is the bifurcation diagram of map (48) which shows that the displacement x
1 of the first oscillator M
1 varies with the parameter b.Figure 6
Bifurcation diagram of map (48).From Figure6, we see that for 6.5
<
b
<
6.7107356 there is a chaotic band in the system and the motion states of the system are complex. When b is equal to 6.69 and the right-hand side of impact surface is selected as the Poincaré section, the velocity of the oscillator M
1 changing with the number of iteration is shown in Figure 7 which shows that the system has the periodic motion. When b is 6.68, the variation of velocity of the oscillator M
1 with the number of iteration is shown in Figure 8. By means of the Lyapunov exponents, we say that chaos appears in the system. When b
>
6.7107365, there is only one fixed point on the Poincaré section for each b. As expected, the grazing bifurcation occurs at b
=
6.7107365. Thus structural stability of the system begins to change such that transition from chaotic motion to a stable period-1 periodic motion occurs.Figure 7
The velocity of the oscillatorM
1 versus the number of iteration for b
=
6.69.Figure 8
The velocity of the oscillatorM
1 versus the number of iteration for b
=
6.68.To compare the numerical results of the map (48) with that of the original system (3)-(4) under the same parameter values, now we begin to do numerical simulation for the original system. At first, we take bifurcation parameter b
=
6.8 and plot the phase portrait of the oscillator M
1 as shown in Figure 9. It can be seen that there is a nonimpacting periodic orbit located at the right-hand side of the grazing bifurcation point of the original system. When b
=
6.6, Figure 10 indicates that there exists the chaotic motion located at the left-hand side of the grazing bifurcation point of the original system. These numerical results are in agreement with that of the map (48) (see Figure 6).Figure 9
Phase portrait of the oscillatorM
1 with b
=
6.8.Figure 10
Phase portrait of the oscillatorM
1 with b
=
6.6.


## 4.1. Derivation of Discontinuity Mapping of Double Grazing Orbit

Assume that there is a double grazing orbit, denoted asγ, in the system (8). And the grazing point of constraint A is marked as X
1
* with X
1
*
=
(
x
1
A
*
,
v
1
A
*
,
x
2
A
*
,
v
2
A
*
). Define a function h
1
(
X
) in a neighborhood B
1
(
X
1
*
) of X
1
* by
(20)
h
1
(
X
)
=
x
1
-
b
.
Therefore, h
1
(
X
)
<
0 implies that the mass M
1 has no impact with the constraint plane A, while h
1
(
X
)
>
0 corresponds to a penetration of the constraint plane A. Thus the impact surface Σ
1 of constraint A can be defined as
(21)
Σ
1
=
{
X
∈
B
1
(
X
1
*
)
∣
h
1
(
X
)
=
0
}
.
A new function H
1
(
X
,
t
)
=
h
1
(
Φ
t
(
X
)
) may be obtained by inserting the solution of the system, X
=
Φ
t
(
X
), into h
1
(
X
), which describes the minimum distance to the constraint A after a flight of time t from a starting point X.According to the relationship between the orbit and the constraintC, we can define some function similarly as before. The grazing point of constraint C is written as X
2
* with X
2
*
=
(
x
1
C
*
,
v
1
C
*
,
x
2
C
*
,
v
2
C
*
). Define a function in a neighbourhood B
2
(
X
2
*
) of X
2
* as follows:
(22)
h
2
(
X
)
=
x
1
+
b
.
Hence, h
2
(
X
)
>
0 implicates that the mass M
1 has no impact with the constraint plane C, while h
2
(
X
)
<
0 corresponds to a penetration of the constraint plane C. The notation Σ
2 is used to denote the impact surface of constraint C; namely,
(23)
Σ
2
=
{
X
∈
B
2
(
X
2
*
)
∣
h
2
(
X
)
=
0
}
.
Another new function H
2
(
X
,
t
)
=
h
2
(
Φ
t
(
X
)
) can be obtained if we insert the solution of the system X
=
Φ
t
(
X
) into h
2
(
X
), which indicates the distance to the constraint C after a flight of time t with X as a starting point.For a double grazing periodic motion, we have(24)
H
1
(
X
1
*
,
0
)
=
0
,
∂
H
1
∂
t
(
X
1
*
,
0
)
=
0
,
A
g
1
=
∂
2
H
1
∂
t
2
(
X
1
*
,
0
)
<
0
,
H
2
(
X
2
*
,
0
)
=
0
,
∂
H
2
∂
t
(
X
2
*
,
0
)
=
0
,
A
g
2
=
∂
2
H
2
∂
t
2
(
X
2
*
,
0
)
>
0
.
If X is restricted on the impact surfaces Σ
1 and Σ
2, respectively, then from the signs of ν
1
=
(
∂
H
1
/
∂
t
)
(
X
,
0
) and ν
2
=
(
∂
H
2
/
∂
t
)
(
X
,
0
), the impact surfaces can be divided into the following subsets:
(25)
Σ
1
+
=
{
X
∈
Σ
1
∣
ν
1
>
0
}
,
Σ
1
0
=
{
X
∈
Σ
1
∣
ν
1
=
0
}
,
Σ
1
-
=
{
X
∈
Σ
1
∣
ν
1
<
0
}
,
Σ
2
+
=
{
X
∈
Σ
2
∣
ν
2
>
0
}
,
Σ
2
0
=
{
X
∈
Σ
2
∣
ν
2
=
0
}
,
Σ
2
-
=
{
X
∈
Σ
2
∣
ν
2
<
0
}
.
So with the definitions above, the impact process on the impact surfaces Σ
1 and Σ
2 can be expressed as two mappings; that is,
(26)
G
1
:
Σ
1
+
⟶
Σ
1
-
,
G
2
:
Σ
2
-
⟶
Σ
2
+
.Select constant phase planes just before the impact as the Poincaré sections; namely,(27)
Π
1
=
{
(
X
,
θ
)
∣
θ
=
θ
1
}
,
Π
2
=
{
(
X
,
θ
)
∣
θ
=
θ
2
}
,
where θ
1 and θ
2 are the phase angles of the system just before the oscillator M
1 contacts with the constraints A and C, respectively. Let X
-
1
=
γ
∩
Π
1 and X
-
2
=
γ
∩
Π
2. Thus based on the assumptions of periodicity and symmetry of the system, we have Φ
T
/
2
(
X
-
1
)
=
X
-
2 and Φ
T
/
2
(
X
-
2
)
=
X
-
1, where Φ
T
/
2
(
X
) is the smooth flow mapping of the system with T
=
2
n
π
/
ω. As in Figure 5, let X
1
* be the image of X
-
1 after the time of flight t
1, X
-
2 the image of X
1
* after the time of flight t
2, X
2
* the image of X
-
2 after the time of flight t
3, and X
-
1 the image of X
2
* after the time of flight t
4. Then the total time in a period is T
=
t
1
+
t
2
+
t
3
+
t
4. If the points in B
1
(
X
-
1
)
∩
Π
1 are mapped by the flow mapping Φ
t
1
(
X
) without consideration of impact effect, the image of an impacting point would penetrate the surface Σ
1, which obviously leads to an error. This problem would happen for points inB
2
(
X
-
2
)
∩
Π
2 if impact effect is ignored. To deal with the impact incident, two mappings C
1 and C
2 must be constructed for points near X
1
* and X
2
*, which will enable us to write a valid Poincaré map. Consider
(28)
P
(
X
)
=
Φ
t
4
∘
C
2
∘
Φ
t
3
∘
Φ
t
2
∘
C
1
∘
Φ
t
1
(
X
)
.Figure 5
The schematic figure of the discontinuity mappingsC
1 and C
2.In what follows, the mappingsC
1 and C
2 will be derived. Starting at a point X
1 in the neighborhood B
1
(
X
1
*
) of X
1
*, H
1
(
X
1
,
t
) indicates the distance from the oscillator M
1 to the constraint A as the variation of time t. Since X
1
∈
B
1
(
X
1
*
), H
1
(
X
1
,
t
) will reach a local maximum through a small time interval τ
1
=
τ
1
(
X
1
). Introduce a function ψ
1
(
X
1
)
=
H
1
(
X
1
,
τ
1
(
X
1
)
) such that ψ
1
(
X
1
)
≥
0 for impacting points on the constraint A and ψ
1
(
X
1
)
<
0 for nonimpacting points. Similarly, the distance function H
2
(
X
2
,
t
) will get to a local minimum through a small time interval τ
2
=
τ
2
(
X
2
), where X
2 is a starting point in a neighborhood B
2
(
X
2
*
) of X
2
*. A function ψ
2
(
X
2
)
=
H
2
(
X
2
,
τ
2
(
X
2
)
) is introduced such that ψ
2
(
X
2
)
≤
0 for impacting points, while ψ
2
(
X
2
)
>
0 for nonimpacting points.For impacting points on constraint planeA, we define a mapping Φ
t
c
1, where t
c
1 is the time of flight from the impacting point to the crossing point of Σ
1
+ along the flow. Hence, the mapping Φ
t
c
1 takes an impacting point near the constraint A to the impact surface and then under the impact mapping G
1 a point is obtained on Σ
1
- after rigid impact. Since the mapping G
1
∘
Φ
t
c
1 does not take place in zero time [10], it must be composited with Φ
-
t
c
1to reach this point. So for X
1
∈
B
1
(
X
1
*
), the full discontinuity mapping C
1 may be written as follows:
(29)
C
1
=
{
Φ
-
t
c
1
∘
G
1
∘
Φ
t
c
1
ψ
1
≥
0
,
I
ψ
1
<
0
,
where G
1
(
X
1
)
=
(
x
1
A
+
,
v
1
A
+
,
x
2
A
+
,
v
2
A
+
)
T
=
(
x
1
A
-
,
-
R
v
1
A
-
,
x
2
A
-
,
v
2
A
-
)
T. The full discontinuity mapping C
2 can be written similarly as follows:
(30)
C
2
=
{
Φ
-
t
c
2
∘
G
2
∘
Φ
t
c
2
ψ
2
≤
0
,
I
ψ
2
>
0
,
where I is the identity map, t
c
2 is the time of flight from the impacting point of constraint C to crossing point of Σ
2
-, and G
2
(
X
2
)
=
(
x
1
C
+
,
v
1
C
+
,
x
2
C
+
,
v
2
C
+
)
T
=
(
x
1
C
-
,
-
R
v
1
C
-
,
x
2
C
-
,
v
2
C
-
)
T. We now consider the series expansions of the discontinuity mappings C
1 and C
2.Forψ
1
(
X
1
)
≥
0, an impacting point of constraint A, denoted as X
1, can be mapped to the impact surface Σ
1 along the flow. For simplicity, we will still use the notation X
1 to denote the image of the impacting point, which is on Σ
1. By the impact law, it follows that
(31)
Δ
X
1
=
G
1
(
X
1
)
-
X
1
=
(
x
1
A
+
v
1
A
+
x
2
A
+
v
2
A
+
)
-
(
x
1
A
-
v
1
A
-
x
2
A
-
v
2
A
-
)
=
-
v
1
A
-
(
0
R
+
1
0
0
)
≜
-
v
1
A
-
δ
,
where δ
=
(
0
,
R
+
1,0
,
0
)
T.To determinet
c
1, we start by expanding the function H
1 with respect to t. For X
1
∈
B
1
(
X
-
1
), near the time τ
1
(
X
1
), H
1
(
X
1
,
t
) has the form
(32)
H
1
(
X
1
,
t
)
=
ψ
1
(
X
1
)
+
(
t
-
τ
1
(
X
1
)
)
2
A
g
1
2
+
o
(
ψ
1
(
X
1
)
,
t
)
.
From (8), it is clear that A
g
1
=
(
∂
2
H
1
/
∂
t
2
)
(
X
1
*
,
0
)
=
-
x
1
A
*
+
x
2
A
*
-
2
ζ
v
1
A
*
+
2
ζ
v
2
A
*
+
(
1
-
f
2
)
sin
τ. Since H
1
(
X
1
,
t
c
1
)
=
0, we have
(33)
t
c
1
=
τ
1
(
X
1
)
-
-
2
ψ
1
(
X
1
)
A
g
1
+
o
(
X
1
)
.
Expressing ν
1 as a function of X
1 yields
(34)
ν
1
=
∂
H
1
∂
t
(
X
1
,
t
c
1
)
=
-
2
A
g
1
ψ
1
(
X
1
)
+
o
(
X
1
,
t
)
.
Hence,
(35)
Δ
X
1
=
-
v
1
A
-
δ
=
-
ν
1
δ
=
-
-
2
A
g
1
ψ
1
(
X
1
)
δ
+
o
(
X
1
)
.
For X
1
∈
B
1
(
X
1
*
),
(36)
Φ
t
c
1
(
X
1
)
=
X
1
+
F
(
X
1
*
)
t
c
1
+
o
(
X
1
,
t
c
1
)
;
then
(37)
Φ
-
t
c
1
(
X
1
)
∘
G
1
∘
Φ
t
c
1
(
X
1
)
=
X
1
-
-
2
A
g
1
ψ
1
(
X
1
)
δ
+
o
(
X
1
)
.
Finally, we obtain
(38)
C
1
(
X
1
)
=
{
X
1
-
ψ
1
(
X
1
)
-
2
A
g
1
δ
+
o
(
X
1
2
)
,
ψ
1
≥
0
,
X
1
,
ψ
1
<
0
.
The series expression of C
2 can be given in the same way; that is,
(39)
C
2
(
X
2
)
=
{
X
2
+
-
ψ
2
(
X
2
)
2
A
g
2
δ
+
o
(
X
2
2
)
,
ψ
2
≤
0
,
X
2
,
ψ
2
>
0
.
Therefore, the Poincaré map can be written in the following form:
(40)
P
=
Φ
t
4
∘
C
2
∘
Φ
t
3
∘
Φ
t
2
∘
C
1
∘
Φ
t
1
.
Let t
1
→
0 and t
3
→
0; the Poincaré map can be simplified as
(41)
P
=
Φ
T
/
2
∘
C
2
∘
Φ
T
/
2
∘
C
1
.
Let L
1
=
(
∂
Φ
T
/
2
/
∂
X
)
(
X
1
*
), L
2
=
(
∂
Φ
T
/
2
/
∂
X
)
(
X
2
*
), η
1
=
(
∂
h
1
/
∂
X
)
(
X
1
*
), and η
2
=
(
∂
h
2
/
∂
X
)
(
X
2
*
). Thus Φ
T
/
2, ψ
1, ψ
2 can be linearized as the expressions of L
1, L
2, η
1, and η
2. Here L
1 and L
2 are both 4
×
4 matrices, and the elements can be computed by the chain rule of the composite function. For X
1
∈
Π
1, we will discuss the expression of the Poincaré map in the following cases.(
1
)  For the case ψ
1
(
X
1
)
<
0, namely, η
1
(
X
1
-
X
1
*
)
<
0, it yields
(42)
X
2
=
Φ
T
/
2
∘
C
1
(
X
1
)
=
X
2
*
+
L
1
(
C
1
(
X
1
)
-
X
1
*
)
+
o
(
X
1
2
)
=
X
2
*
+
L
1
(
X
1
-
X
1
*
)
+
o
(
X
1
2
)
.
(i)
Ifψ
2
(
X
2
)
>
0, namely, η
2
(
X
2
-
X
2
*
)
>
0, then
(43)
Φ
T
/
2
∘
C
2
∘
Φ
T
/
2
∘
C
1
(
X
1
)
=
Φ
T
/
2
∘
C
2
(
X
2
)
=
X
1
*
+
L
2
(
X
2
-
X
2
*
)
+
o
(
X
2
2
)
.
(ii)
Ifψ
2
(
X
2
)
≤
0, namely, η
2
(
X
2
-
X
2
*
)
≤
0, then
(44)
Φ
T
/
2
∘
C
2
∘
Φ
T
/
2
∘
C
1
(
X
1
)
=
Φ
T
/
2
∘
C
2
(
X
2
)
=
X
1
*
+
L
2
(
C
2
(
X
2
)
-
X
2
*
)
+
o
(
X
2
2
)
=
X
1
*
+
L
2
(
X
2
+
-
2
A
g
2
ψ
2
(
X
2
)
δ
-
X
2
*
)
+
o
(
X
2
2
)
.(
2
) For the case ψ
1
(
X
1
)
≥
0, namely, η
1
(
X
1
-
X
1
*
)
≥
0, it follows that
(45)
X
2
=
Φ
T
/
2
∘
C
1
(
X
1
)
=
X
2
*
+
L
1
(
C
1
(
X
1
)
-
X
1
*
)
+
o
(
X
1
2
)
=
X
2
*
+
L
1
(
X
1
-
-
2
A
g
1
ψ
1
(
X
1
)
δ
-
X
1
*
)
+
o
(
X
1
2
)
.
(i)
Ifψ
2
(
X
2
)
>
0, namely, η
2
(
X
2
-
X
2
*
)
>
0, then
(46)
Φ
T
/
2
∘
C
2
∘
Φ
T
/
2
∘
C
1
(
X
1
)
=
Φ
T
/
2
∘
C
2
(
X
2
)
=
X
1
*
+
L
2
(
X
2
-
X
2
*
)
+
o
(
X
2
2
)
.
(ii)
Ifψ
2
(
X
2
)
≤
0, namely, η
2
(
X
2
-
X
2
*
)
≤
0, then
(47)
Φ
T
/
2
∘
C
2
∘
Φ
T
/
2
∘
C
1
(
X
1
)
=
Φ
T
/
2
∘
C
2
(
X
2
)
=
X
1
*
+
L
2
(
C
2
(
X
2
)
-
X
2
*
)
+
o
(
X
2
2
)
=
X
1
*
+
L
2
(
X
2
+
-
2
A
g
2
ψ
2
(
X
2
)
δ
-
X
2
*
)
+
o
(
X
2
2
)
.Combining expressions (43)-(44) and (46)-(47), for arbitrary point X
1
∈
B
1
(
X
1
*
), we have (48)
P
=
Φ
T
/
2
∘
C
2
∘
Φ
T
/
2
∘
C
1
=
{
X
1
*
+
L
2
(
X
2
-
X
2
*
)
+
o
(
X
2
2
)
,
η
1
(
X
1
-
X
1
*
)
<
0
,
η
2
(
X
2
-
X
2
*
)
>
0
;
X
1
*
+
L
2
(
X
2
+
-
2
A
g
2
ψ
2
(
X
2
)
δ
-
X
2
*
)
+
o
(
X
2
2
)
,
η
1
(
X
1
-
X
1
*
)
<
0
,
η
2
(
X
2
-
X
2
*
)
≤
0
;
X
1
*
+
L
2
(
X
2
-
X
2
*
)
+
o
(
X
2
2
)
,
η
1
(
X
1
-
X
1
*
)
≥
0
,
η
2
(
X
2
-
X
2
*
)
>
0
;
X
1
*
+
L
2
(
X
2
+
-
2
A
g
2
ψ
2
(
X
2
)
δ
-
X
2
*
)
+
o
(
X
2
2
)
,
η
1
(
X
1
-
X
1
*
)
≥
0
,
η
2
(
X
2
-
X
2
*
)
≤
0
.


## 4.2. Numerical Simulation

In this subsection, we will use a set of parameter values that generates the grazing periodic motion as given in Section3 to investigate the grazing bifurcation of the system.TakingR
=
0.8, μ
k
=
5, μ
m
=
10, ζ
=
0.05, f
2
=
0, ω
=
0.63, and τ
=
2.229311, let the parameter b gradually increase from 6.5 to 6.9 by a step of 0.0002. For each b, we take 1500 iterations and the first 1300 of which are omitted as the transient process. Figure 6 is the bifurcation diagram of map (48) which shows that the displacement x
1 of the first oscillator M
1 varies with the parameter b.Figure 6
Bifurcation diagram of map (48).From Figure6, we see that for 6.5
<
b
<
6.7107356 there is a chaotic band in the system and the motion states of the system are complex. When b is equal to 6.69 and the right-hand side of impact surface is selected as the Poincaré section, the velocity of the oscillator M
1 changing with the number of iteration is shown in Figure 7 which shows that the system has the periodic motion. When b is 6.68, the variation of velocity of the oscillator M
1 with the number of iteration is shown in Figure 8. By means of the Lyapunov exponents, we say that chaos appears in the system. When b
>
6.7107365, there is only one fixed point on the Poincaré section for each b. As expected, the grazing bifurcation occurs at b
=
6.7107365. Thus structural stability of the system begins to change such that transition from chaotic motion to a stable period-1 periodic motion occurs.Figure 7
The velocity of the oscillatorM
1 versus the number of iteration for b
=
6.69.Figure 8
The velocity of the oscillatorM
1 versus the number of iteration for b
=
6.68.To compare the numerical results of the map (48) with that of the original system (3)-(4) under the same parameter values, now we begin to do numerical simulation for the original system. At first, we take bifurcation parameter b
=
6.8 and plot the phase portrait of the oscillator M
1 as shown in Figure 9. It can be seen that there is a nonimpacting periodic orbit located at the right-hand side of the grazing bifurcation point of the original system. When b
=
6.6, Figure 10 indicates that there exists the chaotic motion located at the left-hand side of the grazing bifurcation point of the original system. These numerical results are in agreement with that of the map (48) (see Figure 6).Figure 9
Phase portrait of the oscillatorM
1 with b
=
6.8.Figure 10
Phase portrait of the oscillatorM
1 with b
=
6.6.


## 5. Conclusions

For a two-degree-of-freedom vibroimpact system with symmetrical rigid stops, using the initial condition and the periodicity condition of grazing periodic motion, the parameter condition is derived that may generate a double grazing period-n motion. According to the derived condition, a group of parameters are obtained for the original vibroimpact system, which leads to a double grazing periodic motion.By the improved discontinuity mapping presented in this paper, the Poincaré map is constructed near the known double grazing periodic orbit. Thus grazing bifurcation of the system is analyzed using the Poincaré map with clearanceb as a bifurcation parameter. Numerical simulation shows that there exists a transition from the chaotic band to a period-1 periodic motion via grazing bifurcation. This numerical result of the Poincaré map is in agreement with that of the original vibroimpact system. In the later study, the discontinuity mapping and the Poincaré map for the nonsmooth dynamical systems with many constraints can be set up in a similar way.



---

*Source: 642589-2014-06-25.xml*