# Normal Form for High-Dimensional Nonlinear System and Its Application to a Viscoelastic Moving Belt

Authors: S. P. Chen; Y. H. Qian
Journal: Abstract and Applied Analysis
 (2014)

---

## Abstract

This paper is concerned with the computation of the normal form and its application to a viscoelastic moving belt. First, a new computation method is proposed for significantly refining the normal forms for high-dimensional nonlinear systems. The improved method is described in detail by analyzing the four-dimensional nonlinear dynamical systems whose Jacobian matrices evaluated at an equilibrium point contain three different cases, that are, (i) two pairs of pure imaginary eigenvalues, (ii) one nonsemisimple double zero and a pair of pure imaginary eigenvalues, and (iii) two nonsemisimple double zero eigenvalues. Then, three explicit formulae are derived, herein, which can be used to compute the coefficients of the normal form and the associated nonlinear transformation. Finally, employing the present method, we study the nonlinear oscillation of the viscoelastic moving belt under parametric excitations. The stability and bifurcation of the nonlinear vibration system are studied. Through the illustrative example, the feasibility and merit of this novel method are also demonstrated and discussed.


---

## Body

## 1. Introduction

Bifurcation and stability analysis of nonlinear differential equations is one of the challenging problems of mathematicians and engineers. Normal form theory is one of the most important tools for such analysis [1–3]. The normal form theory for differential equations can be dated back to the pioneer works of the renowned mathematician Poincaré [4]. He tried to use some change of variables to alter nonlinear systems into linear ones. The idea of the method is to simplify the system such that the topological behavior of the system in the vicinity of a singularity point remains unchanged.Some recent developments of the theory of normal form can be found in [5–15]. Chen et al. [5] presented the renormalization group theory, which was used for the search of normal forms for large classes of finite-dimensional vector fields [6–8]. Stróżyna and Żołądek [9, 10] made use of the time rescaling to achieve the results on the orbital equivalence of vector fields. In addition, Dullin and Meiss [11] and Murdock [12] dealt with the problems for the normal forms of nilpotent systems, whose linear part about the origin is a nilpotent matrix. Benderesky and Churchill [13, 14] and Sanders [15] studied the spectral sequences for the normal forms of vectors. Such spectral sequences can provide valuable information on the normal forms.On the other hand, Kuznetsov [16] considered the normal form theory in an application-oriented way for computation of high-dimensional nonlinear systems. Zhang and his coworkers [17, 18] employed the adjoint operator method to obtain the higher-order normal forms of high-dimensional nonlinear dynamical systems and the associated nonlinear transformations. Zhang and Leung [19] considered a general four-dimensional normal form of a double Hopf bifurcation. Yu and his associates [20–22] developed efficient computing methods for parametric normal forms. They also applied the new method to consider controlling bifurcations of the nonlinear dynamical systems. Chen and Dora [23] were devoted to the development of effective methods for further reductions of the classical normal forms for complex dynamical systems. From a practical point of view, a better understanding and knowledge of the normal forms of various complex nonlinear systems will further promote the potential interest for the analysis of real engineering problems.The normal form theory plays an important role in the study of bifurcation behavior of differential dynamical systems. Itovich and Moiola [24] made use of the frequency domain and the normal form methodologies to analyze the unfolding of a nonresonant double Hopf singularity. Zhang et al. [25] employed the center manifold reduction and normal form method to obtain the singular bifurcation of a ring of three coupled advertising oscillators with delay. Jiang and Yuan [26] studied the classical Van Der Pol equation with Bogdanov-Takens singularity and bifurcation. Gattulli et al. [27] analyzed the postcritical behavior of a single degree of freedom system equipped with a Tuned Mass Damper for double Hopf bifurcation in the neighbourhood of 1 : 1 resonance. Using the normal form method and the center manifold theory, Li and his associates [28] investigated the double Hopf bifurcation of the trivial equilibrium for delay-coupled limit cycle oscillators. Buono and Belair [29] studied the normal form of a vector field, which is generated by a scalar delay-differential equations at nonresonant double Hopf bifurcation points.Some of the above works use the Jordan canonical form of the leading matrixA. However, it is well known that handling the eigenvalues and Jordan canonical forms is very difficult in computer algebra system. In this paper, a new computation method by direct computation is developed to refine the normal forms for high-dimensional nonlinear systems. We do not need to compute the Jordan canonical form of A nor its eigenvalues. Our method is applicable in both the nilpotent and the nonnilpotent cases.In this paper, we will develop an efficient method for computing the normal forms directly for general four-dimension systems and apply the method to consider controlling bifurcations. The approach is efficient since it does not require the computation of the Jordan canonical form ofA or its eigenvalues. Besides, the proposed method is applied to investigate the nonlinear oscillations of a viscoelastic moving belt under parametric excitations. The rest of the paper is organized as follows. In Section 2, the essential idea behind the method in [23] is briefly introduced. The new computation method is described in detail by analyzing the four-dimensional nonlinear dynamical systems in Section 3. The applications to stability and bifurcation analysis on the viscoelastic moving belt are presented in Sections 4 and 5 to show the efficiency of the method. Finally, conclusions are drawn in Section 6.


## 2. Normal Forms for Nonlinear System

Consider a dynamical system described by the following differential equation:(1)
x
˙
=
F
(
x
)
=
A
x
+
∑
k
=
2
∞
f
k
(
x
)
,
x
∈
R
n
,
where A
x represents the linear part, A is the Jacobian matrix, and f
k
(
x
) denotes the kth-order vector homogeneous polynomials of x.Without loss of generality,A is expressed in terms of the standard Jordan canonical form. Note that system (1) is assumed to have an equilibrium at the origin x
=
0.We take the coordinate transformation as follows:(2)
x
=
y
+
φ
k
(
y
)
.
Substituting (2) into (1) gives
(3)
y
˙
=
(
I
+
∂
y
φ
k
)
-
1
F
(
y
+
φ
k
(
y
)
)
,
(4)
(
I
+
∂
y
φ
k
(
y
)
)
-
1
=
I
-
∂
y
φ
k
(
y
)
+
O
(
∥
y
∥
2
k
-
2
)
,
where ∂
y
φ
k is the Jacobian matrix of φ
k with respect to y.Then, (4) is substituted back into (3) to form
(5)
y
˙
=
A
y
+
f
2
(
y
)
+
⋯
+
f
k
-
1
(
y
)
+
{
f
k
(
y
)
-
[
∂
y
φ
k
(
y
)
A
y
-
A
φ
k
(
y
)
]
}
+
O
(
∥
y
∥
k
+
1
)
.Employing the normal form theory, we introduce a linear operator as follows:(6)
L
A
k
:
H
k
n
⟶
H
k
n
,
(
L
A
k
φ
k
)
(
y
)
=
∂
y
φ
k
(
y
)
A
y
-
A
φ
k
(
y
)
,
φ
k
∈
H
k
n
.
Hence, we can write
(7)
H
k
n
=
R
k
+
C
k
,
where R
k represents the range of L
A
k and C
k is the complementary space to R
k.Hence, thekth-order terms can be simplified to
(8)
g
k
(
y
)
∈
C
k
.The rationale for the classical normal form theory can be explained by the following theorem (see [30]).Theorem 1.
Let the notations be the same as above. Suppose that the decomposition (7) is given for k
=
2
,
…
,
N. Then, there exists a sequence of near identity changes of variables x
=
y
+
φ
k
(
y
), in which φ
k
(
y
)
∈
H
k
n. Therefore, the dynamical system (1) is transformed into
(9)
y
˙
=
A
y
+
g
2
(
y
)
+
⋯
+
g
N
(
y
)
+
O
(
∥
y
∥
N
+
1
)
,
where g
k
(
y
)
∈
C
k for k
=
2
,
…
,
N.By applying the Takens normal form theory [30], one arrives at the kth-order normal form g
k
(
y
), while those parts belonging to R
k can be removed by properly choosing the coefficients of the nonlinear transformation φ
k
(
y
).


## 3. Computation of Normal Forms and Their Coefficients

Consider a four-dimensional generalized averaged system withZ
2
⊕
Z
2-symmetry governed by
(10)
x
˙
=
F
(
x
)
=
A
x
+
f
3
(
x
)
,
x
∈
R
4
,
where f
3
(
x
)
∈
H
4
3; that is
(11)
f
3
(
x
)
=
(
f
1
3
(
x
)
,
f
2
3
(
x
)
,
f
3
3
(
x
)
,
f
4
3
(
x
)
)
T
=
(
∑
|
m
|
=
3
a
m
1
m
2
m
3
m
4
x
1
m
1
x
2
m
2
x
3
m
3
x
4
m
4
,
∑
|
m
|
=
3
b
m
1
m
2
m
3
m
4
x
1
m
1
x
2
m
2
x
3
m
3
x
4
m
4
,
∑
|
m
|
=
3
c
m
1
m
2
m
3
m
4
x
1
m
1
x
2
m
2
x
3
m
3
x
4
m
4
,
∑
|
m
|
=
3
d
m
1
m
2
m
3
m
4
x
1
m
1
x
2
m
2
x
3
m
3
x
4
m
4
)
T
,
and |
m
|
=
m
1
+
m
2
+
m
3
+
m
4.At the same time, (6) becomes
(12)
∂
y
φ
3
(
y
)
A
y
-
A
φ
3
(
y
)
=
f
3
-
g
3
,
The problem at hand is to determine φ, so that g
3 contains the smallest possible number of monomials.Consider the following three cases of the Jordan matrixA in four-dimensional nonlinear systems.(i)
The Jordan matrixA has two pairs of pure imaginary eigenvalues.
(ii)
The Jordan matrixA has one nonsemisimple double zero and a pair of pure imaginary eigenvalues.
(iii)
The Jordan matrixA has two nonsemisimple double zero eigenvalues.The forms of the Jordan matrixA in these cases (i)–(iii) can be, respectively, represented by (13)–(15) as follows:
(13)
A
=
[
0
-
w
1
0
0
w
1
0
0
0
0
0
0
-
w
2
0
0
w
2
0
]
=
[
J
1
0
0
J
2
]
,
(14)
A
=
[
0
1
0
0
0
0
0
0
0
0
0
-
w
0
0
w
0
]
=
[
J
1
0
0
J
2
]
,
(15)
A
=
[
0
1
0
0
0
0
0
0
0
0
0
1
0
0
0
0
]
=
[
J
1
0
0
J
2
]
.It is easy to see that the three-order polynomial solutions in four variables can be obtained from (12). To achieve this, we write
(16)
f
3
(
x
)
=
(
∑
|
q
|
=
3
α
q
x
q
∑
|
q
|
=
3
β
q
x
q
∑
|
q
|
=
3
γ
q
x
q
∑
|
q
|
=
3
η
q
x
q
)
,
g
3
(
y
)
=
(
∑
|
q
|
=
3
α
q
′
y
q
∑
|
q
|
=
3
β
q
′
y
q
∑
|
q
|
=
3
γ
q
′
y
q
∑
|
q
|
=
3
η
q
′
y
q
)
,
φ
(
y
)
=
(
∑
|
q
|
=
3
a
q
y
q
∑
|
q
|
=
3
b
q
y
q
∑
|
q
|
=
3
c
q
y
q
∑
|
q
|
=
3
d
q
y
q
)
,
with a
q
,
b
q
,
c
q
,
d
q
,
α
q
′
,
β
q
′
,
γ
q
′ and η
q
′ to be determined later.For the case (i), (12) can be expressed as(17a)
∑
|
q
|
=
3
[
-
w
1
(
q
1
+
1
)
a
q
+
e
1
-
e
2
+
w
1
(
q
2
+
1
)
a
q
-
e
1
+
e
2
-
w
2
(
q
3
+
1
)
a
q
+
e
3
-
e
4
+
w
2
(
q
4
+
1
)
a
q
-
e
3
+
e
4
+
w
1
b
q
]
y
q
=
∑
|
q
|
=
3
(
α
q
-
α
q
′
)
y
q
,
(17b)
∑
|
q
|
=
3
[
-
w
1
(
q
1
+
1
)
b
q
+
e
1
-
e
2
+
w
1
(
q
2
+
1
)
b
q
-
e
1
+
e
2
-
w
2
(
q
3
+
1
)
b
q
+
e
3
-
e
4
+
w
2
(
q
4
+
1
)
b
q
-
e
3
+
e
4
-
w
1
a
q
]
y
q
=
∑
|
q
|
=
3
(
β
q
-
β
q
′
)
y
q
,
(17c)
∑
|
q
|
=
3
[
-
w
1
(
q
1
+
1
)
c
q
+
e
1
-
e
2
+
w
1
(
q
2
+
1
)
c
q
-
e
1
+
e
2
-
w
2
(
q
3
+
1
)
c
q
+
e
3
-
e
4
+
w
2
(
q
4
+
1
)
c
q
-
e
3
+
e
4
+
w
2
d
q
]
y
q
=
∑
|
q
|
=
3
(
γ
q
-
γ
q
′
)
y
q
,
(17d)
∑
|
q
|
=
3
[
-
w
1
(
q
1
+
1
)
d
q
+
e
1
-
e
2
+
w
1
(
q
2
+
1
)
d
q
-
e
1
+
e
2
-
w
2
(
q
3
+
1
)
d
q
+
e
3
-
e
4
+
w
2
(
q
4
+
1
)
d
q
-
e
3
+
e
4
-
w
2
c
q
]
y
q
=
∑
|
q
|
=
3
(
η
q
-
η
q
′
)
y
q
,where e
1
=
(
1
,
0
,
0
,
0
), e
2
=
(
0
,
1
,
0
,
0
), e
3
=
(
0
,
0
,
1
,
0
), and e
4
=
(
0
,
0
,
0
,
1
).Then, we solve the following equations:(18a)
-
w
1
(
q
1
+
1
)
a
q
+
e
1
-
e
2
+
w
1
(
q
2
+
1
)
a
q
-
e
1
+
e
2
-
w
2
(
q
3
+
1
)
a
q
+
e
3
-
e
4
+
w
2
(
q
4
+
1
)
a
q
-
e
3
+
e
4
+
w
1
b
q
=
α
q
-
α
q
′
,
(18b)
-
w
1
(
q
1
+
1
)
b
q
+
e
1
-
e
2
+
w
1
(
q
2
+
1
)
b
q
-
e
1
+
e
2
-
w
2
(
q
3
+
1
)
b
q
+
e
3
-
e
4
+
w
2
(
q
4
+
1
)
b
q
-
e
3
+
e
4
-
w
1
a
q
=
β
q
-
β
q
′
,
(18c)
-
w
1
(
q
1
+
1
)
c
q
+
e
1
-
e
2
+
w
1
(
q
2
+
1
)
c
q
-
e
1
+
e
2
-
w
2
(
q
3
+
1
)
c
q
+
e
3
-
e
4
+
w
2
(
q
4
+
1
)
c
q
-
e
3
+
e
4
+
w
2
d
q
=
γ
q
-
γ
q
′
,
(18d)
-
w
1
(
q
1
+
1
)
d
q
+
e
1
-
e
2
+
w
1
(
q
2
+
1
)
d
q
-
e
1
+
e
2
-
w
2
(
q
3
+
1
)
d
q
+
e
3
-
e
4
+
w
2
(
q
4
+
1
)
d
q
-
e
3
+
e
4
-
w
2
c
q
=
η
q
-
η
q
′
.For the case (ii), we obtain(19a)
∑
|
q
|
=
3
[
(
q
1
+
1
)
a
q
+
e
1
-
e
2
-
w
(
q
3
+
1
)
a
q
+
e
3
-
e
4
+
w
(
q
4
+
1
)
a
q
-
e
3
+
e
4
-
b
q
]
y
q
=
∑
|
q
|
=
3
(
α
q
-
α
q
′
)
y
q
,
(19b)
∑
|
q
|
=
3
[
(
q
1
+
1
)
b
q
+
e
1
-
e
2
-
w
(
q
3
+
1
)
b
q
+
e
3
-
e
4
+
w
(
q
4
+
1
)
b
q
-
e
3
+
e
4
]
y
q
=
∑
|
q
|
=
3
(
β
q
-
β
q
′
)
y
q
,
(19c)
∑
|
q
|
=
3
[
(
q
1
+
1
)
c
q
+
e
1
-
e
2
-
w
(
q
3
+
1
)
c
q
+
e
3
-
e
4
+
w
(
q
4
+
1
)
c
q
-
e
3
+
e
4
+
w
d
q
]
y
q
=
∑
|
q
|
=
3
(
γ
q
-
γ
q
′
)
y
q
,
(19d)
∑
|
q
|
=
3
[
(
q
1
+
1
)
d
q
+
e
1
-
e
2
-
w
(
q
3
+
1
)
d
q
+
e
3
-
e
4
+
w
(
q
4
+
1
)
d
q
-
e
3
+
e
4
-
w
c
q
]
y
q
=
∑
|
q
|
=
3
(
η
q
-
η
q
′
)
y
q
,and we need to solve the following equations:(20a)
(
q
1
+
1
)
a
q
+
e
1
-
e
2
-
w
(
q
3
+
1
)
a
q
+
e
3
-
e
4
+
w
(
q
4
+
1
)
a
q
-
e
3
+
e
4
-
b
q
=
α
q
-
α
q
′
,
(20b)
(
q
1
+
1
)
b
q
+
e
1
-
e
2
-
w
(
q
3
+
1
)
b
q
+
e
3
-
e
4
+
w
(
q
4
+
1
)
b
q
-
e
3
+
e
4
=
β
q
-
β
q
′
,
(20c)
(
q
1
+
1
)
c
q
+
e
1
-
e
2
-
w
(
q
3
+
1
)
c
q
+
e
3
-
e
4
+
w
(
q
4
+
1
)
c
q
-
e
3
+
e
4
+
w
d
q
=
γ
q
-
γ
q
′
,
(20d)
(
q
1
+
1
)
d
q
+
e
1
-
e
2
-
w
(
q
3
+
1
)
d
q
+
e
3
-
e
4
+
w
(
q
4
+
1
)
d
q
-
e
3
+
e
4
-
w
c
q
=
η
q
-
η
q
′
.Finally, the case (iii) leads to(21a)
∑
|
q
|
=
3
[
(
q
1
+
1
)
a
q
+
e
1
-
e
2
+
(
q
3
+
1
)
a
q
+
e
3
-
e
4
-
b
q
]
y
q
=
∑
|
q
|
=
3
(
α
q
-
α
q
′
)
y
q
,
(21b)
∑
|
q
|
=
3
[
(
q
1
+
1
)
b
q
+
e
1
-
e
2
+
(
q
3
+
1
)
b
q
+
e
3
-
e
4
]
y
q
=
∑
|
q
|
=
3
(
β
q
-
β
q
′
)
y
q
,
(21c)
∑
|
q
|
=
3
[
(
q
1
+
1
)
c
q
+
e
1
-
e
2
+
(
q
3
+
1
)
c
q
+
e
3
-
e
4
-
d
q
]
y
q
=
∑
|
q
|
=
3
(
γ
q
-
γ
q
′
)
y
q
,
(21d)
∑
|
q
|
=
3
[
(
q
1
+
1
)
d
q
+
e
1
-
e
2
+
(
q
3
+
1
)
d
q
+
e
3
-
e
4
]
y
q
=
∑
|
q
|
=
3
(
η
q
-
η
q
′
)
y
q
.Similarly, the following equations must be solved:(22a)
(
q
1
+
1
)
a
q
+
e
1
-
e
2
+
(
q
3
+
1
)
a
q
+
e
3
-
e
4
-
b
q
=
α
q
-
α
q
′
,
(22b)
(
q
1
+
1
)
b
q
+
e
1
-
e
2
+
(
q
3
+
1
)
b
q
+
e
3
-
e
4
=
β
q
-
β
q
′
,
(22c)
(
q
1
+
1
)
c
q
+
e
1
-
e
2
+
(
q
3
+
1
)
c
q
+
e
3
-
e
4
-
d
q
=
γ
q
-
γ
q
′
,
(22d)
(
q
1
+
1
)
d
q
+
e
1
-
e
2
+
(
q
3
+
1
)
a
q
+
e
3
-
e
4
=
η
q
-
η
q
′
.Let(23)
M
q
=
(
a
q
b
q
c
q
d
q
)
,
Λ
q
=
(
α
q
β
q
γ
q
η
q
)
,
Λ
q
′
=
(
α
q
′
β
q
′
γ
q
′
η
q
′
)
.
Substituting (23) into (18a)–(18d), (20a)–(20d), and (22a)–(22d), separately, the following three sets of 3-order nonlinear algebraic equations are resulted.(i)
For the case of two pairs of pure imaginary eigenvalues, we arrive at(24)
-
w
1
(
q
1
+
1
)
M
q
+
e
1
-
e
2
+
w
1
(
q
2
+
1
)
M
q
-
e
1
+
e
2
-
w
2
(
q
3
+
1
)
M
q
+
e
3
-
e
4
+
w
2
(
q
4
+
1
)
M
q
-
e
3
+
e
4
-
A
M
q
=
Λ
q
-
Λ
q
′
.
(ii)
For the case of one nonsemisimple double zero and a pair of pure imaginary eigenvalues, we get(25)
(
q
1
+
1
)
M
q
+
e
1
-
e
2
-
w
(
q
3
+
1
)
M
q
+
e
3
-
e
4
+
w
(
q
4
+
1
)
M
q
-
e
3
+
e
4
-
A
M
q
=
Λ
q
-
Λ
q
′
.
(iii)
For the case of two nonsemisimple double zero eigenvalues, we derive(26)
(
q
1
+
1
)
M
q
+
e
1
-
e
2
+
(
q
3
+
1
)
M
q
+
e
3
-
e
4
-
A
M
q
=
Λ
q
-
Λ
q
′
.It should be clear thatM
q
=
0 if one of the components of q is completely negative. We choose the lexicographical order on the set {
q
=
(
q
1
,
q
2
,
q
3
,
q
4
)
∈
N
4
:
|
q
|
=
3
}.Letq
1
=
i, q
2
=
j, q
3
=
l, and q
4
=
3
-
i
-
j
-
l. Then, we have
(27)
Z
i
,
j
,
l
=
M
(
i
,
j
,
l
,
3
-
i
-
j
-
l
)
,
Λ
i
,
j
,
l
=
f
3
,
Λ
i
,
j
,
l
′
=
g
3
.
Substituting (27) into (24)–(26) gives the following three conclusions.(i)
For the case of two pairs of pure imaginary eigenvalues, we have(28)
-
w
1
(
i
+
1
)
Z
i
+
1
,
j
-
1
,
l
+
w
1
(
j
+
1
)
Z
i
-
1
,
j
+
1
,
l
-
w
2
(
l
+
1
)
Z
i
,
j
,
l
+
1
+
w
2
(
4
-
i
-
j
-
l
)
Z
i
,
j
,
l
-
1
-
A
Z
i
,
j
,
l
=
Λ
i
,
j
,
l
-
Λ
i
,
j
,
l
′
.
Let Z
0
,
j
,
l
(
0
≤
j
≤
3
,
0
≤
l
≤
3
) be given arbitrarily. We can determine Z
i
+
1
,
j
-
1
,
l by the above formulae to achieve Λ
i
,
j
,
l
′
=
0 for 0
≤
i
≤
k
-
1, 1
≤
j
≤
k, and 0
≤
l
≤
k. In fact, we have
(29)
w
1
(
i
+
1
)
Z
i
+
1
,
j
-
1
,
l
=
w
1
(
j
+
1
)
Z
i
-
1
,
j
+
1
,
l
-
w
2
(
l
+
1
)
Z
i
,
j
,
l
+
1
+
w
2
(
4
-
i
-
j
-
l
)
Z
i
,
j
,
l
-
1
-
A
Z
i
,
j
,
l
-
Λ
i
,
j
,
l
.
The remaining equation is (for j
=
0)
(30)
Λ
i
,
0
,
l
′
=
Λ
i
,
0
,
l
-
w
1
Z
i
-
1
,
1
,
l
+
w
2
(
l
+
1
)
Z
i
,
0
,
l
+
1
-
w
2
(
4
-
i
-
l
)
Z
i
,
0
,
l
-
1
+
A
Z
i
,
0
,
l
.
(ii)
For the case of one nonsemisimple double zero and a pair of pure imaginary eigenvalues, (25) becomes
(31)
(
i
+
1
)
Z
i
+
1
,
j
-
1
,
l
-
w
(
l
+
1
)
Z
i
,
j
,
l
+
1
+
w
(
4
-
i
-
j
-
l
)
Z
i
,
j
,
l
-
1
-
A
Z
i
,
j
,
l
=
Λ
i
,
j
,
l
-
Λ
i
,
j
,
l
′
.
In view of (31), we have
(32)
(
i
+
1
)
Z
i
+
1
,
j
-
1
,
l
=
w
(
l
+
1
)
Z
i
,
j
,
l
+
1
-
w
(
4
-
i
-
j
-
l
)
Z
i
,
j
,
l
-
1
+
A
Z
i
,
j
,
l
+
Λ
i
,
j
,
l
.
Hence, the remaining equation is given by (for j
=
0)
(33)
Λ
i
,
0
,
l
′
=
Λ
i
,
0
,
l
+
w
(
l
+
1
)
Z
i
,
0
,
l
+
1
-
w
(
4
-
i
-
l
)
Z
i
,
0
,
l
-
1
+
A
Z
i
,
0
,
l
.
(iii)
For the case of two nonsemisimple double zero eigenvalues, we deduce(34)
(
i

+

1
)
Z
i
+
1
,
j
-
1
,
l
+
(
l

+

1
)
Z
i
,
j
,
l
+
1
-
A
Z
i
,
j
,
l

=

Λ
i
,
j
,
l
-
Λ
i
,
j
,
l
′
.Similarly, making use of (34) leads to
(35)
(
i
+
1
)
Z
i
+
1
,
j
-
1
,
l
=
Λ
i
,
j
,
l
-
(
l
+
1
)
Z
i
,
j
,
l
+
1
+
A
Z
i
,
j
,
l
.From (34) and (35), the remaining equation can be written as (for j
=
0)
(36)
Λ
i
,
0
,
l
′
=
Λ
i
,
0
,
l
-
(
l
+
1
)
Z
i
,
0
,
l
+
1
+
A
Z
i
,
0
,
l
.We determine allZ
i
,
j
,
l for i
≥
1 as functions of Z
0
,
j
,
l, which may be solved in some cases by making some of the terms Λ
i
,
0
,
l
′ equal zeroes.


## 4. Application to a Viscoelastic Moving Belt Model

In this section, we apply the proposed method in Section3 to a parametrically excited viscoelastic moving belt with the external damping. Consider the viscoelastic moving belt model [31] with cross-sectional area A, length L between two end supports, axial velocity v, and viscous damping coefficient c as shown in Figure 1. A Cartesian coordinate system (Oxyz) is adopted, which is located in the plane of the viscoelastic moving belt. Another coordinate system is a moving coordinate fixed on the belt. The u and w denote the displacements in the x and y directions, respectively. It is assumed that the tension T is characterized as a small periodic perturbation T
1
cos
⁡
Ω
t on the steady-state tension T
0; that is, T
=
T
0
+
T
1
cos
⁡
Ω
t. Since the belt tension is assumed to dominate the transverse stiffness, the bending stiffness of the viscoelastic moving belt is neglected. The equations of motion for the transverse vibration of the belt are based on an axially moving string model. The nondimensional nonlinear governing equation of motion for the viscoelastic moving belt under parametric excitations can be written as follows [31–33]:
(37)
w
,
t
t

+
2
γ
w
,
t
x

+
(
γ
2
-
1
-

a
cos
⁡
ω
t
)
w
,
x
x
+

2
μ
w
,
t

-
N
(
w
)
=
0
,
where the comma subscript denotes the partial differentiation, and
(38)
a

=
T
1
T
0
,
γ
=

v
ρ
A
T
0
,
ω
=

Ω
ρ
A
L
2
T
0
,
E
e
=

E
0
A
T
0
,
E
v
=
η
A
ρ
T
0
L
2
,
μ
=
1
2
c
A
T
0
ρ
,
(39)
N
(
w
)
=
3
2
E
e
(
w
,
x
)
2
w
,
x
x
+
E
v
∂
∂
t
(
1
2
(
w
,
x
)
2
)
w
,
x
x
+
E
v
w
,
x
∂
∂
t
(
w
,
x
w
,
x
x
)
.
The boundary conditions are imposed by
(40)
w
(
0
,
t
)
=
0
,
w
(
L
,
t
)
=
0
.Figure 1
The model of a viscoelastic moving belt and the coordinate system.In the subsequent analysis, we use the method of multiple scales and Galerkin’s approach in the partial differential governing equations of the viscoelastic moving belt. We introduce the mass, gyroscopic, and linear stiffness operators as follows:(41)
M
=
I
,
G
=
2
γ
∂
∂
x
,
K
=
(
γ
2
-
1
)
∂
2
∂
x
2
.
Substituting (41) into (37) leads to the standard symbolic form
(42)
M
w
,
t
t
+
G
w
,
t
+
K
w
=
N
(
w
)
+
a
cos
⁡
ω
t
∂
2
w
∂
x
2
-
2
μ
∂
w
∂
t
.To obtain a system that is suitable for the application of the method of multiple scales, we introduce the scale transformations(43)
a
⟶
ε
a
,
N
⟶
ε
N
,
μ
⟶
ε
μ
,
where ε is the small perturbation parameter.Substituting (43) into (42), we obtain the following dimensionless nonlinear system under parametric excitations
(44)
M
w
,
t
t
+
G
w
,
t
+
K
w
=
ε
N
(
w
)
+
ε
a
cos
⁡
ω
t
∂
2
w
∂
x
2
-
2
ε
μ
∂
w
∂
t
.
The method of multiple scales can now be applied to search for the uniform solutions of (44) in the following form:
(45)
w
(
x
,
t
,
ε
)
=
w
0
(
x
,
T
0
,
T
1
)
+
ε
w
1
(
x
,
T
0
,
T
1
)
+
⋯
,
where T
0
=
t
,
T
1
=
ε
t.The differential operators of the method of multiple scales can be defined as(46a)
d
d
t
=
∂
∂
T
0
∂
T
0
∂
t
+
∂
∂
T
1
∂
T
1
∂
t
+
⋯
=
D
0
+
ε
D
1
+
⋯
,
(46b)
d
2
d
t
2
=
(
D
0
+
ε
D
1
+
⋯
)
2
=
D
0
2
+
2
ε
D
0
D
1
+
⋯
,where D
n
=
∂
/
∂
T
n, n
=
0,1
,
….In this paper, we investigate the case of primary parametric resonance for thenth andlth order modes of (44), and we introduce the following linear transformation:
(47)
x
1
⟶
x
2
,
x
2
⟶
x
1
,
x
3
⟶
x
4
,
x
4
⟶
x
3
.
Then, the four-dimensional averaged equations in the Cartesian form can be obtained as(48a)
x
˙
1
=
-
μ
x
1
-
(
σ
1
-
f
)
x
2
+
a
e
2
v
1
x
1
(
x
1
2
+
x
2
2
)
+
a
e
1
v
2
x
2
(
x
1
2
+
x
2
2
)
+
b
e
2
v
1
x
1
(
x
3
2
+
x
4
2
)
+
b
e
1
v
2
x
2
(
x
3
2
+
x
4
2
)
+
c
e
2
v
1
x
3
(
x
1
2
-
x
2
2
)
-
c
e
1
v
2
x
4
(
x
1
2
-
x
2
2
)
+
2
c
e
1
v
2
x
1
x
2
x
3
+
2
c
e
2
v
1
x
1
x
2
x
4
,
(48b)
x
˙
2
=
(
σ
1
+
f
)
x
1
-
μ
x
2
-
a
e
1
v
2
x
1
(
x
1
2
+
x
2
2
)
+
a
e
2
v
1
x
2
(
x
1
2
+
x
2
2
)
-
b
e
1
v
2
x
1
(
x
3
2
+
x
4
2
)
+
b
e
2
v
1
x
2
(
x
3
2
+
x
4
2
)
-
c
e
1
v
2
x
3
(
x
1
2
-
x
2
2
)
-
c
e
2
v
1
x
4
(
x
1
2
-
x
2
2
)
+
2
c
e
2
v
1
x
1
x
2
x
3
-
2
c
e
1
v
2
x
1
x
2
x
4
,
(48c)
x
˙
3
=
-
μ
x
3
+
(
D
-
σ
2
)
x
4
+
A
e
2
v
1
x
3
(
x
3
2
+
x
4
2
)
+
A
e
1
v
2
x
4
(
x
3
2
+
x
4
2
)
+
B
e
2
v
1
x
3
(
x
1
2
+
x
2
2
)
+
B
e
1
v
2
x
4
(
x
1
2
+
x
2
2
)
+
C
e
2
v
1
x
1
(
x
3
2
-
x
4
2
)
+
C
e
2
v
1
x
1
(
x
3
2
-
x
4
2
)
-
C
e
1
v
2
x
2
(
x
3
2
-
x
4
2
)
+
2
C
e
1
v
2
x
1
x
3
x
4
+
2
C
e
2
v
1
x
2
x
3
x
4
,
(48d)
x
˙
4
=
(
σ
2
+
D
)
x
3
-
μ
x
4
-
A
e
1
v
2
x
3
(
x
3
2
+
x
4
2
)
+
A
e
2
v
1
x
4
(
x
3
2
+
x
4
2
)
-
B
e
1
v
2
x
3
(
x
1
2
+
x
2
2
)
+
B
e
2
v
1
x
4
(
x
1
2
+
x
2
2
)
-
C
e
1
v
2
x
1
(
x
3
2
-
x
4
2
)
-
C
e
2
v
1
x
2
(
x
3
2
-
x
4
2
)
+
2
C
e
2
v
1
x
1
x
3
x
4
-
2
C
e
1
v
2
x
2
x
3
x
4
,where the coefficients presented in (48a)–(48d) are listed in the appendix.The above algorithm applied to the system (48a)–(48d) leads to the following normal form(49a)
x
˙
1
=
-
μ
x
1
-
(
σ
1
-
f
)
x
2
+
a
e
2
v
1
σ
1
σ
1
-
f
x
1
3
+
a
e
2
v
1
σ
1
σ
1
+
f
x
1
x
2
2
+
a
e
1
v
2
(
2
σ
1
2
+
f
2
)
2
(
σ
1
2
-
f
2
)
x
2
x
1
2
+
a
e
1
v
2
(
2
σ
1
2
+
f
2
)
2
(
σ
1
+
f
)
2
x
2
3
+
b
e
2
v
1
σ
2
σ
2
-
D
x
1
x
3
2
+
b
e
2
v
1
σ
2
σ
2
+
D
x
1
x
4
2
+
b
e
1
v
2
σ
1
σ
2
(
σ
2
-
D
)
(
σ
1
+
f
)
x
2
x
3
2
+
b
e
1
v
2
σ
1
σ
2
(
σ
2
+
D
)
(
σ
1
+
f
)
x
2
x
4
2
,
(49b)
x
˙
2
=
(
σ
1
+
f
)
x
1
-
μ
x
2
+
a
e
2
v
1
σ
1
σ
1
-
f
x
2
x
1
2
+
a
e
2
v
1
σ
1
σ
1
+
f
x
2
3
-
a
e
1
v
2
(
2
σ
1
2
+
f
2
)
2
(
σ
1
-
f
)
2
x
1
3
-
a
e
1
v
2
(
2
σ
1
2
+
f
2
)
2
(
σ
1
2
-
f
2
)
x
1
x
2
2
+
b
e
2
v
1
σ
2
σ
2
-
D
x
2
x
3
2
+
b
e
2
v
1
σ
2
σ
2
+
D
x
2
x
4
2
-
b
e
1
v
2
σ
1
σ
2
(
σ
2
-
D
)
(
σ
1
-
f
)
x
1
x
3
2
-
b
e
1
v
2
σ
1
σ
2
(
σ
2
+
D
)
(
σ
1
-
f
)
x
1
x
4
2
,
(49c)
x
˙
3
=
-
μ
x
3
+
(
D
-
σ
2
)
x
4
+
A
e
2
v
1
σ
2
σ
2
-
D
x
3
2
+
A
e
2
v
1
σ
2
σ
2
+
D
x
3
x
4
2
+
A
e
1
v
2
(
2
σ
2
2
+
D
2
)
2
(
σ
2
2
-
D
2
)
x
4
x
3
2
+
A
e
1
v
2
(
2
σ
2
2
+
D
2
)
2
(
σ
2
+
D
)
2
x
4
3
+
B
e
2
v
1
σ
1
σ
1
-
f
x
3
x
1
2
+
B
e
2
v
1
σ
1
σ
1
+
f
x
3
x
2
2
+
B
e
1
v
2
σ
1
σ
2
(
σ
1
-
f
)
(
σ
2
+
D
)
x
4
x
1
2
+
B
e
1
v
2
σ
1
σ
2
(
σ
1
+
f
)
(
σ
2
+
D
)
x
4
x
2
2
,
(49d)
x
˙
4
=
(
σ
2
+
D
)
x
3
-
μ
x
4
+
A
e
2
v
1
σ
2
σ
2
-
D
x
4
x
3
2
+
A
e
2
v
1
σ
2
σ
2
+
D
x
4
3
-
A
e
1
v
2
(
2
σ
2
2
+
D
2
)
2
(
σ
2
-
D
)
2
x
3
3
-
A
e
1
v
2
(
2
σ
2
2
+
D
2
)
2
(
σ
2
2
-
D
2
)
x
3
x
4
2

+

B
e
2
v
1
σ
1
σ
1
-
f
x
4
x
1
2

+

B
e
2
v
1
σ
1
σ
1
+
f
x
4
x
2
2
-
B
e
1
v
2
σ
1
σ
2
(
σ
1
-
f
)
(
σ
2
-
D
)
x
3
x
1
2
-
B
e
1
v
2
σ
1
σ
2
(
σ
1
+
f
)
(
σ
2
-
D
)
x
3
x
2
2
.Comparing the method developed here with other methods given in [17, 31], it is observed that normal form of the averaged (49a)–(49d) is simpler than normal forms obtained in [17, 31].


## 5. Stability and Bifurcation Analysis on the Viscoelastic Moving Belt

It is known that (49a)–(49d) has a trivial zero solution (
x
1
,
x
2
,
x
3
,
x
4
)
=
(
0
,
0
,
0
,
0
) in which the Jacobi matrix can be written as
(50)
J
=
D
x
X
=
[
-
μ
-
(
σ
1
-
f
)
0
0
σ
1
+
f
-
μ
0
0
0
0
-
μ
D
-
σ
2
0
0
σ
2
+
D
-
μ
]
.The characteristic equation corresponding to the trivial zero solution is(51)
(
λ
2
+
2
μ
λ
+
μ
2
+
σ
1
2
-
f
2
)
(
λ
2
+
2
μ
λ
+
μ
2
+
σ
2
2
-
D
2
)
=
0
.The eigenvalues of the above equations are(52)
λ
1
=
-
μ
+
i
σ
1
2
-
f
2
,
λ
2
=
-
μ
-
i
σ
1
2
-
f
2
,
λ
3
=
-
μ
+
i
σ
2
2
-
D
2
,
λ
4
=
-
μ
-
i
σ
2
2
-
D
2
.Take a linear transformation of coordinate as follows:(53)
{
X
}
=
[
D
]

{
Z
}
,
where [
D
] is the matrix of the eigenvectors of the linear part in (49a)–(49d),
(54)
{
X
}
=
[
x
1
,
x
2
,
x
3
,
x
4
]
T
,
{
Z
}
=
[
z
1
,
z
-
1
,
z
2
,
z
-
2
]
T
.
After substituting the above transformation into (49a)–(49d), one obtains the equations in new complex coordinate as follows:
(55a)
z
˙
1
=
(
-
μ
+
i
σ
1
2
-
f
2
)
z
1
+
(
4
a
e
2
v
1
σ
1
σ
1
+
f
-
i
2
a
e
1
v
2
(
2
σ
1
2
+
f
2
)
(
σ
1
+
f
)
σ
1
2
-
f
2
)
z
1
|
z
1
|
2
+
(
4
b
e
2
v
1
σ
2
σ
2
+
D
-
i
4
b
e
1
v
2
σ
1
σ
2
(
σ
2
+
D
)
σ
1
2
-
f
2
)
z
1
|
z
2
|
2
,
(55b)
z
˙
2
=
(
-
u
+
i
σ
2
2
-
D
2
)
z
2
+
(
4
B
e
2
v
1
σ
1
σ
1
+
f
-
i
4
B
e
1
v
2
σ
1
σ
2
(
σ
1
+
f
)
σ
2
2
-
D
2
)
z
2
|
z
1
|
2
+
(
4
A
e
2
v
1
σ
2
σ
2
+
D
-
i
2
A
e
1
v
2
(
2
σ
2
2
+
D
2
)
(
σ
2
+
D
)
σ
2
2
-
D
2
)
z
2
|
z
2
|
2
.Let(56)
z
1
=
r
1
e
i
φ
1
,
z
2
=
r
2
e
i
φ
2
.
In polar coordinates (
r
1
,
r
2
,
φ
1
,
φ
2
), systems (55a) and (55b) can be written as(57a)
r
˙
1
=
r
1
(
-
u
+
4
a
e
2
v
1
σ
1
σ
1
+
f
r
1
2
+
4
b
e
2
v
1
σ
2
σ
2
+
D
r
2
2
)
,
(57b)
r
˙
2
=
r
2
(
-
u
+
4
B
e
2
v
1
σ
1
σ
1
+
f
r
1
2
+
4
A
e
2
v
1
σ
2
σ
2
+
D
r
2
2
)
,
(57c)
φ
˙
1
=
σ
1
2
-
f
2
-
2
a
e
1
v
2
(
2
σ
1
2
+
f
2
)
(
σ
1
+
f
)
σ
1
2
-
f
2
r
1
2
-
4
b
e
1
v
2
σ
1
σ
2
(
σ
2
+
D
)
σ
1
2
-
f
2
r
2
2
,
(57d)
φ
˙
2
=
σ
2
2
-
D
2
-
4
B
e
1
v
2
σ
1
σ
2
(
σ
1
+
f
)
σ
2
2
-
D
2
r
1
2
-
2
A
e
1
v
2
(
2
σ
2
2
+
D
2
)
(
σ
2
+
D
)
σ
2
2
-
D
2
r
2
2
.Since our analysis is local, we can truncate the higher order and consider the system(58a)
r
˙
1
=
r
1
(
-
u
+
4
a
e
2
v
1
σ
1
σ
1
+
f
r
1
2
+
4
b
e
2
v
1
σ
2
σ
2
+
D
r
2
2
)
,
(58b)
r
˙
2
=
r
2
(
-
u
+
4
B
e
2
v
1
σ
1
σ
1
+
f
r
1
2
+
4
A
e
2
v
1
σ
2
σ
2
+
D
r
2
2
)
,
(58c)
φ
˙
1
=
σ
1
2
-
f
2
,
(58d)
φ
˙
2
=
σ
2
2
-
D
2
.The first two equations are independent of the last two. The last two equations describe rotations in the planesr
2
=
0 and r
1
=
0 with angular velocities σ
1
2
-
f
2 and σ
2
2
-
D
2, respectively. This does not change the bifurcation diagrams. So, it is enough to study the first two equations. Therefore, we can study the original four-dimensional system (49a)–(49d) by analyzing the planar system(59a)
r
˙
1
=
r
1
(
-
u
+
4
a
e
2
v
1
σ
1
σ
1
+
f
r
1
2
+
4
b
e
2
v
1
σ
2
σ
2
+
D
r
2
2
)
,
(59b)
r
˙
2
=
r
2
(
-
u
+
4
B
e
2
v
1
σ
1
σ
1
+
f
r
1
2
+
4
A
e
2
v
1
σ
2
σ
2
+
D
r
2
2
)
.This system is called the amplitude system. The trivial equilibrium,r
1
=
0, r
2
=
0, corresponds to the trivial equilibrium of the original system. The study of the amplitude system is simplified if we use squares ρ
1,2 of the amplitudes as follows:
(60)
ρ
k
=
r
k
2
,
k
=
1
,
2
.
The equations for ρ
1,2 read(61a)
ρ
˙
1
=
2
ρ
1
(
-
u
+
4
a
e
2
v
1
σ
1
σ
1
+
f
ρ
1
+
4
b
e
2
v
1
σ
2
σ
2
+
D
ρ
2
)
,
(61b)
ρ
˙
2
=
2
ρ
2
(
-
u
+
4
B
e
2
v
1
σ
1
σ
1
+
f
ρ
1
+
4
A
e
2
v
1
σ
2
σ
2
+
D
ρ
2
)
.The behavior of systems (61a) and (61b) depends on the coefficients of ρ
1 and ρ
2. We start our analysis with the case when 4
a
e
2
v
1
σ
1
/
(
σ
1
+
f
)
<
0 and 4
A
e
2
v
1
σ
2
/
(
σ
2
+
D
)
<
0. The case 4
a
e
2
v
1
σ
1
/
(
σ
1
+
f
)
>
0 and 4
A
e
2
v
1
σ
2
/
(
σ
2
+
D
)
>
0 can be reduced to the previous one by time reversal. The choices 4
a
e
2
v
1
σ
1
/
(
σ
1
+
f
)
<
0 and 4
A
e
2
v
1
σ
2
/
(
σ
2
+
D
)
<
0 imply that both of the primary Hopf bifurcations are supercritical and stable.For the case of16
A
e
2
v
1
a
e
2
v
1
σ
1
σ
2
/
(
σ
1
+
f
)
(
σ
2
+
D
)
<
0, first we can reduce the number of the coefficients in (61a) and (61b) by rescaling. Let
(62)
ξ
1
=
4
a
e
2
v
1
σ
1
σ
1
+
f
ρ
1
,
ξ
2
=
-
4
A
e
2
v
1
σ
2
σ
2
+
D
ρ
2
,
τ
=
2
t
.
We obtain the system(63a)
ξ
˙
1
=
ξ
1
(
-
u
+
ξ
1
-
b
e
2
v
1
A
e
2
v
1
ξ
2
)
,
(63b)
ξ
˙
2
=
ξ
2
(
-
u
+
B
e
2
v
1
a
e
2
v
1
ξ
1
-
ξ
2
)
.The trivial and the nontrivial equilibria of (63a) and (63b) have the representations
(64)
E
1
=
(
u
,
0
)
,
E
2
=
(
0
,
u
)
,
E
3
=
(
a
e
2
v
1
u
(
A
e
2
v
1
-
b
e
2
v
1
)
a
e
2
v
1
A
e
2
v
1
-
b
e
2
v
1
B
e
2
v
1
,
A
e
2
v
1
u
(
B
e
2
v
1
-
a
e
2
v
1
)
a
e
2
v
1
A
e
2
v
1
-
b
e
2
v
1
B
e
2
v
1
)
.The nontrivial equilibriumE
3 can bifurcate. Suppose that 4
a
e
2
v
1
σ
1
/
(
σ
1
+
f
)
>
0 and 4
A
e
2
v
1
σ
2
/
(
σ
2
+
D
)
<
0. The opposite case can be treated similarly. The Hopf bifurcation and consequent existence of cycles are only possible in following three cases: (
1
)
b
e
2
v
1
/
A
e
2
v
1
>
1, B
e
2
v
1
/
a
e
2
v
1
>
1; (
2
)
b
e
2
v
1
/
A
e
2
v
1
>
1, B
e
2
v
1
/
a
e
2
v
1
<
1, b
e
2
v
1
B
e
2
v
1
/
A
e
2
v
1
a
e
2
v
1
>
1; (
3
)
b
e
2
v
1
/
A
e
2
v
1
<
0, B
e
2
v
1
/
a
e
2
v
1
<
0, b
e
2
v
1
B
e
2
v
1
/
A
e
2
v
1
a
e
2
v
1
>
1.The numerical simulation result is given in Figure2, which is the trajectory from the 4-dimensional space (
x
1
,
x
2
,
x
3
,
x
4
) into the 3-dimensional space (
x
1
,
x
2
,
x
3
). Figure 2 shows the trajectory of the systems (61a) and (61b) for μ
=
0.0001, σ
1
=
2, f
=
-
124, a
e
2
v
1
=
56, a
e
1
v
2
=
2.3, b
e
2
v
1
=
15, b
e
1
v
2
=
19.6, σ
2
=
4.5, A
e
2
v
1
=
-
127, A
e
1
v
2
=
12, B
e
2
v
1
=
1.82, B
e
1
v
2
=
12, and D
=
190.Figure 2
Periodic solution of the viscoelastic moving belt system.


## 6. Conclusions

An efficient method for computing the normal form of high-dimensional nonlinear systems is presented in this paper. This computation method is applied to obtain the normal form of the averaged equation for the viscoelastic moving belt under parametric excitations. Based on the current studies, it is found that the newly developed computation method improves the classical normal form. Therefore, it is a further reduction of the classical normal form. Meanwhile, the normal form derived herein is used to explore the bifurcation and stability analysis of axially viscoelastic moving belts under parametric excitations.In contrast to earlier works, this article argues that the normal forms of high-dimensional nonlinear systems may always be achieved without computing either the Jordan canonical form ofA or its eigenvalues. This significantly reduces the difficulties which other computing methods may face in obtaining normal forms. Therefore, it is more convenient to utilize the approach developed here to compute the normal forms of averaged equations for different resonant cases. It is also found that the normal form of the averaged equation by using the new method is simpler than those obtained in [17, 31]. This is very useful for the complex behavior patterns analysis of the nonlinear dynamical systems in an abstract sense. Therefore, this has opened the area for more research works.



---

*Source: 879564-2014-08-27.xml*