An Algorithm for Higher Order Hopf Normal Forms
Shock and Vibration, Vol.
An Algorithm for Higher Order Hopf Normal Forms
0 Department of Civil and Structural Engineering University of Hong Kong Hong Kong

INTRODUCTION
The invariant manifold of a nonlinear oscillator
near an equilibrium or a periodic orbit is deter
mined by the structure of its vector field. Two
oftenused mathematical tools to simplify the
original system are center manifold and normal
forms. The normal form theory is a technique of
transforming the original nonlinear differential
equation to a simpler standard form by appropri
ate changes of coordinates, so that the essential
features of the manifold become more evident.
Basic references on normal forms and their appli
cations may be found in
Poincare (1889)
,
Birkhoff (1927)
,
Arnold (1983)
,
Chow and Hole
(1982)
,
Guckenheimer and Holmes (1983)
, Iooss
and Joseph (1980),
Sethna and Sell (1978)
,
Van
der Beek (1989
),
Vakakis and Rand (1992)
, and
Leung and Zhang (1994)
. In this article a compu
tational approach based on the classical normal
form theory of Poincare and Birkhoff is intro
duced. The relationship of the coefficients
between the original equations and the normal form
equations is explicitly constructed. The tech
nique presented in our approach follows the idea
of
Takens (1973)
. A linear operator and its adja
cent operator are defined with an inner product
on the space of homogeneous polynomials. The
resultant normal form keeps only the resonant
terms, which cannot be eliminated by a nonlinear
polynomial changes of variables, in the kernel of
the adjacent linear operator. The normal form
simplifies the original systems so that the dy
namic stability and bifurcation can be studied in a
standard manner and the classification of mani
fold in the neighborhood of an equilibrium or a
periodic orbit can be achieved with relatively lit
tle efforts. This article is a further development
of previous work
(Leung and Zhang, 1994)
. The
higher order normal forms of several typical non
linear oscillators: cubic system, quadratic sys
tem, and DuffingVan der Pol system are pro
vided and the steadystate solutions of the
method are compared with existing results.
TRANSFORMATION TO NORMAL FORM
Substituting Eq. (7) into Eq. (6) gives
Consider the nonlinear ordinary differential
equations
i = (/ + Dyhk(Y))Y = J(y + hk(y))
+ ... + Fk(y + hk(y)) + O(lylk+I),
(8)
wherefis an nvector function of u, differentia
ble up to order r. Suppose Eq. (1) has a fixed
point at u = Uo. We first perform a few linear
transformations to simplify eq. (1). By the vari
able change v = u  uo, we eliminate the constant
terms and shift the fixed point to the origin under
which Eq. (1) becomes
v = f(v + uo) == H(v),
where H(v) is at least linear in v. We next split
the linear part of the ordinary differential equa
tion and write (2) as follows
v = D"H(O)v + H(v),
where H(v) == H(v)  D"H(O)v and, H(v) =
O(lv F) is at least quadratic in v and D" denotes
differentiation with respect to v. We further
transform D"H(O) into Jordan canonical form by
the canonical matrix T, i.e., v = Tx, and obtain
i = lID"H(O)Tx + TIH(Tx).
which can be written alternatively as
i = Jx + F(x).
where J is the Jordan Canonical form of D"H(O)
and F(x) is the nonlinear part of the equation. In
Poincare's normal form theory, the nonlinear
function F(x) is expanded by a series of homoge
neous polynomials,
i = Jx + F 2(x) + F3(X) + ... + F,(x) (6)
+ O(lxl,+I)
where Fk E Hko which is the set of homogeneous
polynomials or order k. To transform Eq. (6) into
its normal form, Poincare introduces a nearly
identity nonlinear coordinate transformation of
the form
(2)
(3)
(4)
(5)
where / denotes the n x n identity matrix and the
term (/ + Dyhk(y)) is invertible for sufficient
small Y so that
Substituting Eq. (9) into Eq. (8) and applying
similar nearly identity nonlinear transformations
up to the kth order, we obtain
Y = Jy + F 2(y) + ... + FkI(y) + Fk(y)
+ (Jhk(y)  Dyhk(y)Jy) + O(lylk+l)
= Jy + F 2(y) + ... + FkI(y) + Fk(y)
+ adJhk(y) + O(lylk+l)
where the overbar is used to represent the origi
nal polynomial matrices and adJhk(y) denotes an
adjacent operator equivalent to the function of
the Lie Bracket,
To simplify the terms of order k as much as pos
sible, we choose a specific form for hk(y) such
that
(10)
(12)
If Fk(y) is the range of L J(Hk)' then all terms of
order k can be eliminated completely from Eq.
(10). Otherwise we must find a complementary
space Gk to LJ(Hk) and let Hk = LAHk) E9 Gk so
that only terms of order k which are in Gk remain
in the resultant expression. It is interesting to
note that simplifying terms at order k would not
affect the coefficients of any lower terms. How
ever, terms of order higher than k will be
changed. Therefore, it is only necessary to keep
track of the way that the higher order terms are
modified by the successive coordinate transfor
mations. This will be discussed in the next sec
tion.
(7)
NORMAL FORMS OF OSCILLATING
SYSTEMS
where hk(y) is kth order in y, hk(y) E Hb 2 :S
k:s r.
We now specialize the normal form formulation
in twodimensions for the Hopf bifurcation
below. Consider a system with a small perturbed
parameter IL and an equilibrium point with eigen
values ±iwo, wo > O.
where
i = f(x, IL), x E R2, IL E Rl
(13)
where we suppose, after shifting of origin and
canonical transformation,
Perform the Lie bracket operation to each basis
element on H 2,
thus,
LJ (Y01) =
[:0 ;0][Y01]  [2~1 ~][:o
(15)
L J (Y10Y2) = w0 (y~Y1YYl2)'
LJ (y0~) _ Wo
;0][~:J = Wo [2~r2J
and
( Y1Y2) ,LJ ( 0)
2 2 = Wo
Y2 Y 1
(Yi) ,
2Y1Y2
Furthermore, the Taylor expansion of Eq. (13)
gives
where
A ==' Dxf(O, 0) + D2x,J(0, O)IL
= [f31L
(aIL + wo)
(aIL + wo)],
f31L
F
2(Y)
= [F21(Y)] = [a2o a11 a02]{ Yl }. (16)
Fn(Y) b20b110b2 Y1Y2
We would like to find a coordinate transforma
tion (7) so that the nonlinear terms of order k + 1
in the new system of function Y vanish rather
than of order k in the original system of function
x. Inserting Y = (Yt. Y2Y E R2 and hk(y) = (hk1(y),
hk2(y)V E Hk into (13), we have the Lie bracket,
y~
SecondOrder Normal Form
If the smallest order of nonlinear terms appearing
in (13) is two, we try to find a transformation h2
of the form
x = Y + hz(y).
Because the det[Aj] = 8wo> 0, which does not
vanish for Wo > 0, we obtain the following null
complementary space for the homogeneous
equations Aj g = {O},
Therefore, we can eliminate the quadratic term
completely from Eq. (10) by means of the sec
ondorder coordinate transformation (19) and get
the secondorder normal form
{AY + Ah2(y) + ~ ~o Dr;'~!(y) [h2(y )]m}
M
Y = ~ (l)N[Dyh2(y)]N
N=O
M
= Ay + ~ F~l) (y)
n=2
where
F(1)(y) = [
n
F",{(y)]
F,,2.J.(y)
=
[~ aU) Yh~]
.
~ bU)Yh~ , (24)
in which, the number in the superscript parenthe
ses refers to the index of coordinate transforma
tion. From Eq. (22) we know that the comple
mentary space of secondorder normal form, G2,
is null so that we have the following equation
F~l)(y) = F2(y) + Jh2(y)  Dyh2(Y) . J . Y = O.
(25)
Then, we find coefficients Cij, dij in Eq. (19) by
solving the linear algebraic equation, Eq. (25),
C20
c",
CO2
d20
d"
d02
1
3wO
We further perform the nearly identity transfor
mation (7) to third order to obtain the more accu
rate nonlinear description of Eq. (14) in the
neighborhood of the original singular point.
where y represents the new coordinate and x is
its old one.
(26)
(27)
(29)
(30)
(22)
(23)
where
AJ = Wo
0 1
3
0
0
1
0
0
H3 = span {(~), (Y~2), (Y~~), (~), (~i)'
(yfY2)' (y~y~), (~~)}.
Similar to Eq. (20), the matrix representation of
LAH3) is given by
LJ (H3 = H3 . A}
Finally, through a proper coordinate transforma
tion, we determine the thirdorder normal form
of the original system,
where ai and bi are to be determined.
We now want to develop a systematic proce
dure to evaluate the coefficients of normal forms
and the normal transformations of third order.
We see from Eq. (11)
i =
M
L (1)N[Dh3(Y )]N
N=O
M
= Ay + F~I)(y) + L F~2)(y)
where
(31)
(32)
(33)
(34)
(35)
(37)
(38)
(39)
(40)
We observe that the homogenous equations
A 3Jg = {O}, g E R8 have two zero eigenvectors for
the matrix A} corresponding to its two zero ei
genvalues,
eT = {(1, 0, 1,0,0, 1,0, IV
(0, 1,0, 1, 1,0, 1, OV}.
Thus, A3 has a complementary space spanned by
the monomial basis (29),
The resultant linear equation is
A3· g = 'Y}
where
g = {C30, C2)' C\2, C03, d30 , d2), dlb d03Y
'Y} = {am  a), aW + b), alY  a), abY + b), b~~
 b), bW  a), blY  b), abY  alY
n = 3, .. " M, i = 2  j.
We know from Eq. (24) that the thirdorder nor
mal form is equal to its complement space so that
F~2)(y) = F~I)(y) + J. h3(Y)  Dh3(Y) . J . Y
= G3(y). (36)
where T is the canonical matrix consisting of the
eigenvectors of B2 + I. The second equation in
(39) becomes
in which
or
B=
d=
6=
1
°
2
0
°
2
0
1
a=
0
3
0
0
d30
d21
d l2
d03
b~~  bl
bW  al
bIY  bl
b&Y  al
c=
0
0
0
3
a~~  al
aW + bl
a\Y  al
abY + bl
C30
C21
CI2
C03
Dividing Eq. (38) into two parts, we obtain
{
WO(B2 + I) . c = B . ii.  6
wO(B2 + I) . a= B . 6 + ii.'
Note that there are two more unknowns existing
in Eq. (39) than the number of equations. How
ever, the coefficients of the normal form a), bl
can be evaluated independently by performing an
orthogonal linear transformation on one of the
equations in Eq. (39),
T=
0
3
LHS = TI (B2 + l) T . a= Wo
RHS = 11 (B . b + a) =
0
6
Sal + alY + 3am + 3b&Y + bW
By comparing the terms in the last two rows of
Eq. (41), we obtain
L
Y = A Y + 2: W2i+I(Y) + O(jyj2L+3), (44)
i=l
{
al = ~ (aW + 3abY + b~V + 3b&Y) .
bl = ~1 (aW + 3a&Y  bW  3b~U)
Substituting Eq. (42) into (39) and solving Eq.
(39), the thirdorder homogeneous polynomial
(2S) is then determined by
(42)
where
1
SWo
1
4wo
o
o
3a\Y  b&Y + 3a~U + bW
 3abY + b\Y + 3aW + b~U
alY  bbY + a~U  bW
3abY + 5b~U + aW  bIY
alY + 5bbY + 3a~U + bW
abY + b~U + aW + bIY
(43)
Higher Order Normal Form
The higher order normal form in a rectangular
coordinate system is derived accordingly and can
be written as follows,
A
[(a:: wo) (a~: wo)],
W2i+ l (y) = (YT + y~)i [::
~~i][::J
Y =
and L is a given degree of the normal form equa
tion. The validity of such simplification is guar
anteed by the implicit function theorem, as for
each /t near /to, there will be an equilibrium p(/t)
near p(/to) that varies smoothly with /t. The nor
mal form of Eq. (44) in polar coordinates can
then be written as the form
YI = r cos (), Y2 = r sin ()
L
; = r(f3/t + 2: air21) + hot.
L
{} = Wo + a/t + 2: bir2i + hot
i=l
(45)
A general formula for each ai, bi is not pres
ently available but we can derive their expres
sions up to any desired higher order using the
algorithm mentioned above recursively, such
that,
{
{
a3 = 1~8 (35a~d + 5a~~ + 3a~~ + 5am + 5b~) + 3b¥~ + 5b~{ + 35b~J),
b3 = ;;~ (5a~) + 3a~~ + 5a~~ + 35a~i  35b~~  5bf{]  3b~~  5b~~);
a4 = _2156_ (63a~~ + 7aW + 3a~~ + 3a~~ + 7a\~ + 7b~~) + 3b~7] + 3b~g + 7b~7 + 63bb'J),
b4 = ;~ (7a~V + 3a~7] + 3a~g + 7a&V + 63ab'J  63b~~  7bW  3b~~  3b~~  7bm);
and for the eleventh order,
a5 = 10~4 (231aJld + 21a?] + 7afJl + 5a~~
+ 7aW, + 21ai1l + 21bfI) + 7b~J + 5b~{
+ 7bf/ + 21b~J + 231bJ1)),
b5 = 1~;4 (21afI) + 7a~J + 5a~{ + 7arti
+ 21a~J + 231aJY 231bJlrl  21b?]
 7bfJl  5b~  7bCfJ  21M1l);
here the subscripts A, B refer to the indices 10
and 11, respectively.
The complexity of higher order normal forms
rapidly becomes apparent as pointed out by
Leung and Zhang (1994)
. Thus, the algorithm is
required to be implemented with the symbolic
manipulation in Mathematica
(Wolfram, 1991)
.
In the case of m > 5, however, it may also cause
overflow in the computers equipped with con
ventional memory if we apply the nearly identity
normal transformations directly in Eqs. (23) and
(34). A general explicit formula representing the
homogeneous polynomial terms is therefore very
useful to reduce the size of the problem.
If the nearly identical change of coordinate of
order k
(46)
is applied then the new homogeneous polynomi
als of order n can be obtained by
n = k, F:(y) = Fn(Y) + J. hk(y)
 Dhk(y) . J . y;
(48)
n > k, F:(y) = Fn(Y) + [DFnk+1(y) . hk(y)
 Dhk(y) . F nk+l(Y)]n"'2kl
where F:(y) denotes the resulting homogeneous
polynomials, km = min[kh k2], and
kl = floor [~ =~]. k2 = floor [IJ
The operator floor['] gives the greatest integer
less than or eual to the varible in the bracket. The
result confirms that only oddorder terms would
appear in the normal form equation, correspond
ing to the Poincare resonance
(Guckenheimer
and Holmes, 1983)
.
OSCILLATORS WITH ODD
NONLINEARITY
Here we consider the nonlinear differential equa
tions with flu) to be an odd function of u. To
illustrate this, two examples are presented. The
first one is the wellknown Duffing oscillator,
which arises from various physical and engineer
ing problems. The solution of this oscillator has
been thoroughly studied so that the accuracy of
the results can be compared with existing meth
ods. In the second example, we work with an
oscillator having a term of the fifth power.
(50)
(51)
(52)
Consider the Duffing oscillator
with the initial conditions
x + w5x = 8X3
x(O) = a, i(O) = 0,
where the parameter 8 > 0, Wo is the linear fre
quency. In classical perturbation methods, it is
usually assumed that 8 is small. In the present
case, however, 8 need not be small. With the
transformation x = XI. i = WOX2, we rewrite Eq.
(50) as
WoO][XI] + [8 0 ].
X2 Wo xi
Comparing Eq. (52) with the standard thirdorder
normal form {i} = {x} + F3(X) , we have b30 =
(8Iwo) and all the other coefficients are zero. If
we substitute these coefficients into Eq. (42), we
obtain al = 0, b l = (38/8wo).
The thirdorder nonlinear transformation of
coordinate {x} = [J]{y} + h3(y) can be deter
mined by Eq. (43) and the coefficients of h3(y)
are C30 = C03 = C21 = d30 = dl2 = 0,
8 58 8
CI2 = W8o2' d21 = W8o2' d03 = 42·
Wo
Here we use y to denote the new coordinates and
x the old ones. While under the fourth order non
linear transformation, all the coefficients of the
homogeneous polynomial are identically zero,
C40 = C31 = C22 = C\3 = C04 = d40 = d31
= d22 = d\3 = d04 = O.
According to the method mentioned in the last
section, such computation can be repeated up the
the desired orders. Finally, we found that all the
coefficients of the evenorder coordinate homo
geneous polynomials h2n(Y) are equal to zero and
all the evenorder terms are already in their sim
plest forms. Thus, only the oddorder homoge
neous polynomials need to be handled in our
computation.
The coefficients of normal form a2 and b2 can
be obtained by taking the change of variable
2182
{y} = {z} + hs(z); a2 = 0, b2 = 256w6.
Then, the fifthorder nonlinear coordinate
change hs(z) is given by
CSO = Cos = C23 = C41 = dso = d l4 = d32 = 0,
If further transformations beyond the ninth order
were performed, it is interesting to note that all
the following higher order coefficients of normal
forms ai, hi (i = 5,6, 7 . . .) as well as the coordi
nate transformation hk(y) (k = 10, 11, 12, . . .)
vanish. Therefore, the normal form of our cubic
nonlinear differential system has an ultimate de
gree of nine. The dynamic system can be repre
sented by an exact polynomial differential ampli
tude equation in a finite number of terms.
The asymptotic solution of this normalform
equation is
l
UI = r coso
r = 0
To get the steadystate periodic solution in the
original coordinate, we could trace back all the
transformations, for example,
{x} = {y} + h3(y)
= {z} + h5(z) + h3(z + h5(z»
= {w} + h7(W) + h5(W + h7(W» + h3({W}
= {u} + h9(U) + h7({U} + h9(U» + h5({U}
+ h7(W) + h5(U + h7(W»)
+ h9(U) + h7({U} + h9(U»)
+ h9(U) + h7({U} + h9(U»».
(55)
+ h3({U} + h9(U) + h7({U} + h9(U» + h5({U}
The steadystate periodic solution of the Duffing
equation is
x=
(
ea3 23e2a5 167e3a7 3431e3a9)
a  32w6 + 1024w3  16384wg + 524288w3 coso
( ea3 3e2a5 30ge3a7 7033e4a9 )
+ 32w6  128wg + 32768wg  1048576w3 cos 30
( e2a5 31e3a7 191e4a9 )
+ 1024w3 + 32768wg + 524288w6 cos 5(J.
_ (
(J
+ (65~~a:W6) cos 9(J
It is clear to see that the coefficients of the first
three orders of e in Eq. (57) are exactly the same
as those given in Eq. (55). Nevertheless, Eq. (55)
is different from Eq. (57) in the higher order
terms. The result derived by normal form theory
gives a finite, asymptotic power series, but using
the LP method the solution was represented by
an infinite one.
Oscillator with FifthPower Nonlinearity
Consider an oscillator with a fifthpower nonlin
earity as the second example,
i + w6X = ex5
with initial conditions (51). Using the normal
form method described previously, one obtains
the following results:
The terms of asymptotic normal form series ex
pansion of Eq. (59) is also finite with all higher
order normal forms equal to zero.
Ql = 0,
and the fifthorder nonlinear coordinate change is
given by
Then, taking seventhorder transformation of co
ordinates, we obtain a2 = 0, b2 = 0, with h7(W) =
{o}. After that, we directly take the ninthorder
change of variable which yields a3 = 0, b3 =
(215e2/3072wfi), and
The steadystate periodic solution of the oscilla
tor with fifth power nonlinearity is
and
OSCILLATORS WITH EVEN
NONLINEARITY
The normal form method can also be easily ex
tended to analyze an oscillator with even nonlin
earity. A quadratic nonlinearity equation is con
sidered in the following calculation.
Quadratic Nonlinear Oscillator
Consider a system having quadratic nonlinearity,
with initial conditions (51). According to the nor
mal form method, one needs to do a secondor
der homogeneous polynomial transformation of
coordinates to simplify the secondorder terms in
the original equations. The coefficients of this
change can be evaluated by Eq. (26), the results
are,
e 2e
C20 = W3o2' CII = 0, C02 = W3o2' d20 = 0,
2e
d ll = W3o2' d02 = 0.
The subsequent thirdorder variable change gives
(60)
al = 0,
+ ( 98e320a49w6) cos 90.
The fourthorder homogeneous polynomial of
transformation cannot be overlooked in this case
and finally its coefficients are found to be
C31 = CI3 = d40 = d22 = 0,
e3 8e3
C40 = 96' C22 = 96' C04 =
Wo Wo
Then normalform coefficients a2, b2 and the
fifthorder homogeneous polynomial of coordi
nate transformation can be evaluated as
a2 = 0,
The coefficients of the normal form a3, b3 are
Substituting the transformation into the original
equation, we have the steadystate periodic solu
tion of the oscillator with quadratic nonlinearity
__ ea2 _ 13e3a4 31ge5a6 ( _ e2a3
2w2 24w3 + 1728wAo + a 48w~
143e4a5 35005e6a7)
+ 20736w~  20736wb2 cosO
[YI] [eA 1][YI]
Y2 = 1 eA Y2
e [YIY~ + Y~]
0
, (63)
in which, the firstorder near identity transforma
tion is given by
where we choose XI = X, X2 = x. Note that the
higher order terms in e are neglected during the
computation. We take transformations up to or
der five. The resultant normal form of Eq. (63) in
polar coordinates is expressed as follows,
YI = r cosO, Y2 = r sinO
f = er(A  ~ r2)
The coefficients of the intermeidate coordinate
transformations are
3e
C30 = C03 = 0, C21 = CI2 = 8'
e 3e e
d30 = d03 = 4 ' d21 = 8' dl2 = 8 '
respectively. The asymptotic solution ofEq. (65)
could also be obtained if we have traced back all
the nearly identity coordinate transformations
x = ea(A + ~ a2) cos 0 + ..£. a 3cos 30
32 32
+ ( a + ge 3) . 3.
a sm 0 + e a sm 30
32 32
(65)
(66)
Y.
++~'I' Y,
Sink Center
Y,
A. <0
where the constant a denotes the amplitude de
fined by the initial condition. Further discussion
of the normalform amplitude equation of (65)
reveals the Hopf bifurcation of the Duffingvan
der Pol equation on a parametric plane r  A. As
for A < °and a source at (0, 0) surrounded by a
shown in Figure 1, there is a spiral sink at (0, 0)
limit cycle for A > o. The limit cycle evolves
continuously from the center at (0,0) for A = o.
The Hopf bifurcation is of importance in situa
tions where a flowinduced oscillator is subjected
to flutters or selfexciting movements. At such
circumstances, the orbits of the steadystate peri
odic solutions stay on the surface of the para
boloid rotated by A = (l/8) r2.
CONCLUSION
We have presented an arithmetic algorithm to
compute the higher order normal forms. By ap
plying the explicit formula proposed in this work,
we can achieve, in a standard manner, the de
sired higher order normal forms for nonlinear dif
ferential polynomial equations. We found that
the steadystate solution of the undamped Duff
ing equation can be represented by a finite cosine
Y.
Sink
y
•
r
Limit Cycle
Source
series with varying phase in finite polynomial
terms. To show the versatility of the algorithm,
we illustrated an example in which the order of
nonlinearity is not restricted to odd numbers.
The application to limit cycle bifurcation is also
demonstrated by the Duffingvan der Pol oscilla
tor.
The second author wishes to thank Dr. Q.C. Zhang for
his helpful discussion on the topic of normal form. The
research was supported by the Research Grant Coun
cil of Hong Kong.
Algorithm/or Higher Order Hop/ Normal Forms
Journal of
Engineering
Volume 2014
Volume 2014
Sensors
Hindawi Publishing Corporation
ht p:/ www.hindawi.com
Rotating
Volume 2014
International Journal of
Hindawi Publishing Corporation
ht p:/ www.hindawi.com
Active
and
Advances in
Civil Engineering
Hindawi Publishing Corporation
ht p:/ www.hindawi.com
Journal of
Robotics
Hindawi Publishing Corporation
ht p:/ www.hindawi.com
Advances
ctronics
Submit
your manuscr ipts
VLSI
Design
Hindawi Publishing Corporation
ht p:/ www.hindawi.com
Hindawi Publishing Corporation
Hindawi Publishing Corporation
Navigation and
Observation
Hindawi Publishing Corporation
ht p:/ www.hindawi.com
Modelling
Sim
ulation
&
Engineering
International Journal of
Engineering
Electrical
and
Computer
International Journal of
Aerospace
Engineering
Arnold , V. I. , 1983 , Geometrical Methods in the Theory of Ordinary Differential Equations , SpringerVerlag, New York.
Birkhoff , G. D. , 1927 , Dynamical Systems , Vol. 9 , AMS Collection Publications .
Blevins , R. D. , 1977 , FlowInduced Vibration , Van Nostrand Reinhold, New York.
Chow , S. N. , and Hale , J. K. , 1982 , Method of Bifurcation Theory, SpringerVerlag, New York.
Guckenheimer , J. , and Holmes , P. , 1983 , Nonlinear Oscillations , Dynamical Systems, and Bifurcations of Vector Fields , SpringerVerlag, New York.
looss , G., and Joseph , D. D., 1980 , Elementary Bifurcation and Stability Theory , SpringerVerlag, New York.
Leung , A. Y. T. , and Zhang , Q. C. , 1994 , "Normal Form Analysis of Hopf Bifurcation Exemplified by Duffing's Equation," Shock and Vibration , Vol. 1 , pp. 233  239 .
Poincare , H. , 1889 , Les Methods Nouvelles de la Mecanique Celeste , GauthierVillars, Paris.
Sethna , P. R. , and Sell , G. R. , 1978 , "Review of the Hopf Bifurcation and Its Applications," Journal of Applied Mechanics Vol. 45 , pp. 234  235 .
Takens , F. , 1973 , "Normal Form for Certain Sigularities of Vector Fields," Annals of the Institute of Fourier , Vol. 23 , pp. 163  195 .
Vakakis , A. F. , and Rand , R. H. , 1992 , "Normal Mode and Global Dynamics of a TwoDegreeofFreedom Nonlinear SystemI. Low Energies," International Journal NonLinear Mechanics , Vol. 27 , pp. 861  874 .
van der Beek , C. G. A. , 1989 , "Normal Form and Periodic Solutions in the Theory of Nonlinear OscillationExistence and Asymptotic Theory," International Journal NonLinear Mechanics , Vol. 24 , pp. 263  279 .
Wolfram , S. , 1991 , Mathematica: A System for Doing Mathematics by Computer , AddisonWesley, Redwood City, CA.
and Volume 2014