Multiple limit cycles of some strongly nonlinear Liénard–Van der Pol oscillator

Multiple limit cycles of some strongly nonlinear Liénard–Van der Pol oscillator

Applied Mathematics and Computation 270 (2015) 620–630 Contents lists available at ScienceDirect Applied Mathematics and Computation journal homepag...

877KB Sizes 1 Downloads 42 Views

Applied Mathematics and Computation 270 (2015) 620–630

Contents lists available at ScienceDirect

Applied Mathematics and Computation journal homepage: www.elsevier.com/locate/amc

Multiple limit cycles of some strongly nonlinear Liénard–Van der Pol oscillator✩ Xianbo Sun∗ Department of Applied Mathematics, Guangxi University of Finance and Economics, Nanning, 530003 Guangxi, PR China

a r t i c l e

i n f o

Keywords: Melnikov function Limit cycle Cuspidal loop Homoclinic loop

a b s t r a c t In this paper, we investigate limit cycles of some Liénard–Van der Pol oscillator; the system has a double homoclinic loop and two cuspidal loops if the damping effect has vanished. By the related Melnikov function theory and bifurcation theories, the limit cycles near the singular circle and center with their distribution are found. The number of limit cycles obtained also reveals that some recent results on the lower bounds of the maximal number of limit cycles bifurcated from this kind of systems can be improved. © 2015 Elsevier Inc. All rights reserved.

1. Introduction and main results The following second order differential equation, generalized from the classical Van der Pol equation and classical Liénard system by the famous French physicist Alfred-Marie Liénard,

x¨ + f (x)x˙ + g(x) = 0,

(1.1)

is usually called Liénard–Van der Pol equation. It is the classical Liénard equation if g(x) = x and classical Van der Pol equation if g(x) = μx, f (x) = ρ(1 − x2 ). System (1.1) can be applied in many natural fields such as chemical reactions, radio and vacuum tube technology, predator-prey systems, growth of a single species, vibration analysis and so on, see [1–3] and references therein. When describing the motion of a nonlinear oscillator, x represents displacement function, −g(x) represents the reacting force, − f (x) represents the damping effects and limit cycles of (1.1) exhibit the phenomenon of nonlinear oscillations. In many cases the damping effect is very less, system (1.1) is always written in the form:

x¨ + ε f (x)x˙ + g(x) = 0, which can be transformed into the equivalent form:

x˙ = y, y˙ = −g(x) − ε f (x)y,

(1.2)

where |ε | is sufficiently small. There are several research topics regard to Liénard–Van der Pol system (1.2), here only relatively new works are listed, such as the existence and uniqueness of limit cycle [4,5], the center problem [6], algebra limit cycles and hyperbolic limit cycles [7,8], the period function problem [9], the Hopf bifurcation [10–14] and alien limit cycles [15], see these papers and their references for details. In the field of pure mathematics, system (1.2) is still very interesting because the number of limit cycles of system (1.2) is related to the second part of Hilbert’s 16th problem and its weaker versions (weak Hilbert’s 16th ✩ ∗

This work was supported by the National Natural Science Foundation of China (NSFC) (11461001) and Research project of GXUFE (2014SYS10). Tel no.:+86 15977781786. E-mail address: [email protected]

http://dx.doi.org/10.1016/j.amc.2015.08.049 0096-3003/© 2015 Elsevier Inc. All rights reserved.

X. Sun / Applied Mathematics and Computation 270 (2015) 620–630

Case 1: α < 0

Case 4: α = 1

Case 2: α = 0

Case 5: 1 < α < 4

621

Case 3: 0 < α < 1

Case 6: α = 4

Case 7: α > 4

Fig. 1. Seven cases for system (1.3)ε=0 with two cusps.

problem [16], Smale’s 13th problem [17]), see [16–26] and their references therein. For instance, in [26] the authors gave the best result that a cubic system that can have 13 limit cycles. This cubic system was obtained by “Treble-making” a system (1.2) of degree 3 which had been proved it can has 5 limit cycles. In general, both of the classical Van der Pol system and the classical Liénard system always have a unique limit cycle and the distribution can be determined by the parameter. However, many Liénard–Van der Pol equations usually have more than one limit cycle, it is believed that the number of limit cycles depends on the degree of f(x) and g(x). In order to state conveniently, let call system (1.2) of type (m, n) if g(x) and f(x) are polynomials of degree m and n respectively, H(m, n) denote the maximal number of its limit cycles. There are some results on the lower bound of H(m, n) for some concrete m, see the introduction part of a very new paper [19], here we only list the results for case m = 7 for briefness. For m = 7 and 1 ≤ n ≤ 5, Yu and Han [27] proved H(7, 1) ≥ 3, H(7, 2) ≥ 5, H(7, 3) ≥ 6, H(7, 4) ≥ 8 and H(7, 5) ≥ 9 by investigating the small amplitude limit cycles. Romanovski and Han [28] proved H (7, n) ≥ 3n 2 − 9 for n ≥ 7. Very recently, Xiong and Han [19] gave some estimation for more general m and n, particularly, it was proved that for any integer n ≥ m, m ≥ 7,

H (m, n) ≥

1 2





1 (m + 2) ln (m + 2) ln (m + 2) −2 n+ − ln (2) 3 ln (2)



1 ln (3) 5 + 6 3 ln (2)

 ln (3) ln (m + 2) 5 + − . ( m + 2) + ln (2) 2 ln (2)

The aim of this paper is not only to present a new result on multiple limit cycles in a Liénard–Van der Pol equation which has a more complex portrait, but also to try to get more limit cycles for a given this kind equation, that is the equation of type (7, n) with n ≥ 6:

x˙ = y, y˙ = −x(x2 − 1)(x2 − α)2 + ε f (x)y

(1.3)  where 0 < ε  1, f (x) = ni=0 ai x2i , ai are real bounded parameters. The portraits of system (1.3)ε=0 have 7 cases according the value α , see Fig. 1. As far as our experience is concerned, a system with more annuli would have more limit cycles, because in this case the Melnikov functions on different annuli have the same coefficients by our method. Obviously, cases 3 and 5 have 5 annuli, other cases have 3 annuli. In this paper we study case 5, case 3 can be investigated as the same routine. Let 1 < α < 4 that will ensure √ the local maximum of the Hamiltonian function of (1.3)ε=0 at the nilpotent singular point ( α , 0) is larger than that at (1, 0) and smaller than 0 and two cuspidal homoclinic loops exist for ε = 0 as described case 5. The methods used in our paper include the computation of related definite integrals, the parameter α should be fixed. In the following we fix α = 3 for convenience. Our main result is given as follows. Theorem A. Taking n = 2 n = 6, 8, 10, 12 and 14, system (1.3) can have 11, 13, 17, 17 and 20 limit cycles respectively, this implies that H(7, 6) ≥ 11, H(7, 8) ≥ 13, H(7, 10) ≥ 17, H(7, 12) ≥ 17, H(7, 14) ≥ 20. Remark. (1) The distribution of these limit cycles are given in four propositions in Section 4. As we can see, most of limit cycles obtained n+1 , that is because it plays no role for each case are near the singular circle and center. In the perturbation there is no x2 in the related Melnikov function and a symmetric perturbed system ensure us getting more limit cycles easily. (2) As far as we know, there are few results on limit cycles of this kind oscillator system that has a double homoclinic loop and two cuspidal loops. (3) Substituting m = 7 and n = 6, 8, 10, 12, 14 into the estimation of H(m, n) of [19] respectively, it shows that the lower bounds obtained in this paper are larger for the case n = 6, 8, 10, 12, 14. Therefore this work not only is a subsequent work of [27], but also reveals that the lower bounds of H(7, 8), H(7, 10), H(7, 12) and H(7, 14) given in [19] can be enlarged. (4) Our results are based on some relative new Melnikov function theory and bifurcation theory, especially the combination of more than one Melnikov functions is considered, however, used few in the field of oscillation. The method can be still applied to some other generalization form of oscillator equations such as Rayleigh–Liénard oscillator and Duffing–Van der Pol oscillator with multiple perturbation.

622

X. Sun / Applied Mathematics and Computation 270 (2015) 620–630

The rest of our paper is organized as follows: in Section 2, we give some preliminaries and the bifurcation theorem which will be used to study the number of zeros of Melnikov function. In Section 3, we investigate the expansions of Melniikov functions of system (1.3) and compute their first coefficients. In Section 4, we use those coefficients and bifurcation theory to find limit cycles for each case, except to case n = 12. The analysis of other four cases are similar but have their own computation skills, therefore we list one by one for readability. 2. Preliminary lemmas and the bifurcation theorem In this section we introduce some known results used in this paper and the bifurcation theorem. As we know the second part of the well-known Hilbert’s 16th problem concerns the number and configuration of limit cycles for the polynomial vector fields in the plane. The problem is very difficult and remains to be unsolved. A weaker version of this problem proposed by Arnold is to study the zeros of Abelian integrals obtained by integrating polynomial 1-forms over ovals of polynomial Hamiltonian, this is the weak Hilbert’s 16th problem. To study the second part of the Hilbert’s 16th problem and the weak Hilbert’s 16th problem, most mathematicians concern the following perturbed Hamiltonian system which is called a near-Hamiltonian system:

x˙ = Hy + ε p(x, y, δ), y˙ = −Hx + ε q(x, y, δ),

(2.1)

where p, q and H are Cω functions, ε is a small parameter, δ ∈ D ⊂ RN with D compact. The unperturbed system of (2.1) is

x˙ = Hy ,

y˙ = −Hx .

(2.2)

A general assumption is that the above system has a family of periodic orbits  h defined by {H (x, y) = h, h ∈ J}, which forms an annulus of (2.2), where J is an open interval. The so called Abelian integral, or the first order Melnikov function of system (2.1) is the following integral:

M(h, δ) =



h

qdx − pdy,

h ∈ J,

(2.3)

one zero of which corresponds to a limit cycle of system (2.1) under some conditions. Considering the near-Hamiltonian system (2.1) and its unperturbed system (2.2), we make the following assumptions on system (2.2): (A1) System (2.2) is symmetric with respect to y-axis,  (A2) System (2.2) has a double homoclinic loop  = 1 2 consisting of two homoclinic loops 1 =  |x≤0 and 2 =  |x≥0 with hyperbolic saddle O(0, 0), where  is defined by the equation H (x, y) = 0. (A3) Inside the homoclinic loops  1 and  2 , system (2.2) has two cuspidal homoclinic loops 1∗ and 2∗ defined by H (x, y) = hs passing through the nilpotent cusps S1 (xs1 , ys1 ) and S2 (xs2 , ys2 ), respectively. (A4) System (2.2) has two centers 1c (xc1 , yc1 ) and 2c (xc2 , yc2 ) with xs1 < xc1 , xs2 > xc2 and H (1c ) = H (2c ) = hc < hs . The portrait of system (2.2) is just the same as that of system (1.3)ε=0 of case 5 in Fig. 1. Then, there exist three families of periodic orbits h1 , h2 and h3 given by

hi : H (x, y) = h, H ∈ (hc , hs )



(hs , 0), i = 1, 2,

 : H (x, y) = h, H ∈ (0, +∞), 3 h

(2.4)

hi surrounds the cuspidal homoclinic loop i∗ and the center ic (xci , yci ), i = 1, 2. see Fig. 1 (case 5). h3 surrounds the double homoclinic loop  . Correspondingly, we have three Melnikov functions

Mi (h, δ) = M3 (h, δ) =



hi



h3

qdx − pdy, h ∈ (hc , 0),

i = 1, 2,

qdx − pdy, h ∈ (0, +∞).

(2.5)

By symmetry, we have M1 (h, δ) = M2 (h, δ). For the expansions of the Melnikov functions near the centers and the homoclinic loops there exists the following results. From [29], we have Lemma 2.1 ([29]). Consider the above Melnikov functions and suppose

H (x, y) =

λ 2

and

p(x, y, δ) =

(y2 − x2 ) +



h i j xi y j ,

λ = 0

(2.6)

i+ j≥3

i+ j≥0

ai j xi y j , q(x, y, δ) =

i+ j≥0

bi j xi y j .

(2.7)

X. Sun / Applied Mathematics and Computation 270 (2015) 620–630

623

for (x, y) near the saddle O(0, 0). Then,

M1 (h, δ) = c0 + c1 h ln |h| + c2 h + c3 h2 ln |h| + O(h2 ), 0 < −h  1, M3 (h, δ) = 2c0 + 2c1 h ln |h| + 2c2 h + 2c3 h2 ln |h| + O(h2 ), 0 < h  1, where

(2.8)

 1 1 c0 = M(0, δ) = qdx − pdy, c1 = − (a¯ 10 + b¯ 01 ) = − ( px + qy )(0, 0, 0), |λ| |λ| 1  c2 = ( px + qy )dt if (a¯ 10 + b¯ 01 ) = 0, 1

1 1 ¯ 30 − b¯ 21 + a¯ 12 + 3b¯ 03 ) − [(2b¯ 02 + a¯ 11 )(3h¯ 03 − h¯ 21 ) {( − 3a 2|λ|λ λ ¯ 1. + (2a¯ 20 + b¯ 11 )(3h¯ 30 − h¯ 12 )]} + bc

c3 = −

Lemma 2.2 ([29]). Consider the above expansions given (2.8) for the Melnikov functions, if there exists δ 0 , such that c j (δ0 ) = 0, j = 0, 1, . . . , k − 1 and ck (δ0 ) = 0 and

rank

∂(c0 , c1 , . . . , ck−1 ) = k, ∂δ

then system (2.1) can have [ 5k 2 ] limit cycles near  for some (ε , δ ) near (0, δ 0 ). From [30], for the expansions and bifurcation theory we have Lemma 2.3 ([30,31]). Consider the above Melnikov function M1 and suppose

H (x, y) =

y2 + h30 x3 + hi j x(x − xs1 )i y(y − ys1 ) j , h30 < 0, 2

(2.9)

i+ j≥3

and

p(x, y, δ) =





ai j x(x − xs1 )i y(y − ys1 ) j , q(x, y, δ) =

i+ j≥0

bi j x(x − xs1 )i y(y − ys1 ) j .

(2.10)

i+ j≥0

for (x, y) near S1 (xs1 , ys1 ). Then, there exist some constants B00 > 0, B∗00 < 0, B10 > 0, B∗10 < 0, b and b∗ such that

M1 (h, δ) = c0 + B00 c1 |h| 6 + (c2 + bc1 )h + B10 c3 |h| 6 − 5

7

1 11 B00 c4 |h| 6 + O(h2 ) 11

(2.11)

for 0 ≤ −(h − hs )  1, and

1 ∗ 11 B c4 h 6 + O(h2 ) 11 00

M1 (h, δ) = c0 + B∗00 c1 h 6 + (c2 + b∗ c1 )h + B∗10 c3 h 6 + 5

7

for 0 ≤ (h − hs )  1, where

c0 (δ) = c2 (δ) =

 

1∗

qdx − pdy,

√ −1 c1 (δ) = 2 2h303 (a10 + b01 ),

( px + qy − a10 − b01 )dt, √ −5

1 c3 (δ) = 2 2h303 h30 (2a20 + b11 − h12 (a10 + b01 )) + (h212 − 2h40 )(a10 + b01 ) , 1∗

3

c4 (δ) = 9μ

α μ ( μ − 20μ1 μ2 μ3 + 4μ μ4 )α00 + (4μ21 μ3 − 10μ1 μ22 )α10 + 4μ μ α μα 1 2 μ1 = 3 h30 , μ2 = − h−303 ( − 2h40 + h221 ), μ3 μ4

α00 α10

−1 1

−7 3 01 − 2 1 [ 20 2 2 3 1 2 20 − 1 30 ],

2 1

6 1 − 53 h [12h30 (h50 − h31 h21 + h12 h221 ) − 4h240 + 4h40 h221 − h421 )], = 36 30 1 − 83 h [5h621 − 40h340 + 30h40 h221 (2h40 − h221 ) + 144h30 h40 (h50 − h31 h21 + h12 h221 ), =− 648 30 + 72h30 h221 ( − h50 + h31 h21 − h12 h221 ) + 216h230 ( − h60 − h22 h21 ]2 + h41 h21 + h03 h321 ), + 108h230 h231 + 432h230 h12 h21 (h12 h21 − h31 )], √ = 2 2(a10 + b01 ), √ = 2 2[−h12 (a10 + b01 ) + 2a20 + b11 ],

(2.12)

624

X. Sun / Applied Mathematics and Computation 270 (2015) 620–630





3 α20 = 2 2 (a10 + b01 ) 3h03 h21 − h22 + h212 − 2h12 a20 − h12 b11 + 3a30 + b21 − a11 h21 − 2b02 h21 ), 2 

 √ 5 α30 = 2 2 (a10 + b01 ) 3h13 h21 + 3h03 h31 + 3h12 h22 − 15h12 h03 h21 − h312 − h32 2

+(3a11 + 6b02 )h12 h21 − (2(b12 + a21 )h21 − (a11 + 2b02 )h31 + 4a40 + b31 + (3b11 + 6a20 )h03 h21

α01





3 − (2a20 + b11 )h22 + 3a20 + b11 h212 − (3a30 + b21 h12 ) , 2 √ = 2 2[a12 + 2b03 − 2h03 a11 − 4h03 b02 + (a10 + b01 )(5h203 − 2h04 )].

Lemma 2.4 ([30]). Consider the above expansions given (2.11) and (2.12) for M1 , if there exists a δ 0 , such that c j (δ0 ) = 0, j = 0, 1, . . . , l − 1 and cl (δ 0 ) = 0 and

rank

∂(c0 , c1 , . . . , cl−1 ) = l, ∂δ

∗ then system (2.1) can have [ 3l+1 2 ] limit cycles near the cuspidal loop 1 for some (ε , δ ) near (0, δ 0 ).

For (x, y) near (xc1 , yc1 ), we may suppose

1 2  (x + y2 ) + hi j xi y j , 2 i+ j≥3   = ai j xi y j , q(x, y, δ) = b i j xi y j .

H (x, y) = p(x, y, δ)

i+ j≥0

(2.13)

i+ j≥0

About the expansion of M1 (h, δ ) near the elementary center (xc1 , yc1 ) we have (see [32])

M1 (h, δ) =



b j (δ)h j+1 , 0 < h  1.

(2.14)

j≥0

Under (2.13), the coefficients can be obtained by using the programs in [32]. Summarizing the above lemmas, we can get the following bifurcation theorem. Theorem 2.1. Let system (2.1) satisfy the assumptions (A1)–(A4) and the above expansions hold. Taking 1 , 2 , 3 , 4 and

5 positive and very small, h1 = hc + 1 , h2 = hs − 2 , h3 = hs + 3 , h4 = 0 − 4 , h5 = 0 + 5 and h6 > h5 positive, let N1 = 1−sgn(M1 (h1 ,δ0 )M1 (h2 ,δ0 )) , 2

1−sgn(M (h ,δ )M1 (h4 ,δ0 ))

1 3 0 N2 = 2 Suppose that there exists a δ 0 , such that

and N3 =

1−sgn(M3 (h5 ,δ0 )M1 (h6 ,δ0 )) . 2

c j (δ0 ) = 0, j = 0, 1, . . . , k − 1, and ck (δ0 ) = 0,

(2.15)

c j (δ0 ) = 0, j = 0, 1, . . . , l − 1, and cl (δ0 ) = 0,

(2.16)

b j (δ0 ) = 0, j = 0, 1, . . . , m − 1, and bm (δ0 ) = 0,

(2.17)

and

rank

∂(c0 , c1 , . . . , ck−1 , c0 , c1 , . . . , cl−1 , b0 , b1 , . . . , bm−1 , ) = k + l + m, ∂δ

then system (2.1) can have



5k 2





+2×

3l + 1 2

 + 2m + 2N1 + 2N2 + N3

∗ limit cycles for some (ε , δ ) near (0, δ 0 ), in detail, [ 3l+1 2 ] limit cycles near the cuspidal loop i and m limit cycles near the center ic (xci , yci ) for i = 1, 2, N1 limit cycle between the curve h1 = {(x, y)|H (x, y) = h1 } and the curve h2 = {(x, y)|H (x, y) = h2 } of each side of y-axis, N2 limit cycle between the curves h3 = {(x, y)|H (x, y) = h3 } and h4 = {(x, y)|H (x, y) = h4 } of each side of y−axis, another N3 limit cycle surrounds the double homoclinic loop  .

Theorem 2.1 can be proved similarly as the proofs of Lemmas 2.2 and 2.4 by implicit function theorem and symmetry, for briefness we omit its proof here.

X. Sun / Applied Mathematics and Computation 270 (2015) 620–630

625

Fig. 2. The portrait of system (1.3)ε=0 with α = 3.

3. Expansions of Melnikov functions of system In this section we investigate the expansions of the related Melnikov functions corresponding Liénard–Van der Pol system (1.3) of type (7, 14), then let some ai vanish, we discuss the type (7, 6), (7, 8), (7, 10) and (7, 14) respectively and prove the main results. System (1.3)ε=0 is a Hamiltonian system with Hamiltonian function 2 2 (x, y) = y − 9 x2 + 15 x4 − 7 x6 + 1 x8 ≡ y + A(x). H 2 2 4 6 8 2

(3.1)

(x, y) = h) of Hamiltonian function (3.1) are sketched in Fig. 2. It is easy to check that H (x, y) defines three The level sets (i.e. H (x, y) = h, h ∈ ( − 43 , − 9 )  ( − 9 , 0)} (i = 1, 2) and  3 = {(x, y)|H (x, y) = h, h ∈ (0, +∞)} which families of ovals hi = {(x, y)|H 24 8 8 h (x, y) = − 43 are closed clockwise orbits of system (1.3)ε=0 inside three periodic annuluses also denoted by hi on the plane. H 24  c c  defines two elementary centers 1 ( − 1, 0) and 2 (1, 0), H (x, y) = 0 defines the double homoclinic loop  = 1 2 , where 1 = {(x, y)|H (x, y) = 0, x < 0} and 2 = {(x, y)|H (x, y) = 0, x > 0}. 1∗ = {(x, y)|H (x, y) = − 89 , x < 0} and 2∗ = {(x, y)|H (x, y) = √ √ − 98 , x > 0} are two homoclinic loops passing through nilpotent saddle ( − 3, 0) and ( 3, 0) of system (1.3)ε=0 respectively. Correspondingly, there exist three Melnikov functions of system (1.3)

M1 (h, δ) = M2 (h, δ) = M3 (h, δ) =



n

h3 i=0



n

h1 i=0

43

ai x2i ydx, h ∈ −

24



,0 ,

ai x2i ydx, h ∈ (0, +∞).

(3.2)

By Lemmas 2.1 and 2.2, we have

M1 (h, δ) = c0 + c1 h ln |h| + c2 h + c3 h2 ln |h| + O(h2 ), 0 < −h  1,

 1 ∗ 9 5 7 11  1, B00 c4 h 6 + O(h2 ), 0 < h + M1 (h, δ) = c0 + B∗00 c1 h 6 + (c2 + b∗ c1 )h + B∗10 c3 h 6 + 11 8

 1 9 5 7 11 B00 c4 |h| 6 + O(h2 ), 0 < − h +  1, M1 (h, δ) = c0 + B00 c1 |h| 6 + (c2 + bc1 )h + B10 c3 |h| 6 − 11 8

 43 b j h j+1 , 0 < h +  1, M1 (h, δ) = 24 j≥0

M3 (h, δ) = 2c0 + 2c1 h ln |h| + 2c2 h + 2c3 h2 ln |h| + O(h2 ), 0 < h  1,

(3.3)

Applying Lemma 2.1, we have

 c0 =

7

1 i=0

ai x2i ydx = 2

i=0

√ 3

1054+162

where m0 =

7 

0 m0

ai x2i



0 − 2A(x)dx =

7

Ii ai ,

(3.4)

i=0

√ √ √ √ 3 43((1054+162 43)2/3 −26+28 1054+162 43) √ , Ii √ 3 3 1054+162 43

=2

0 m0

ai x2i



0 − 2A(x)dx.

Applying the formula for c1 , we get

c1 = −a0 .

(3.5)

7 0  a x2i  Let c1 = 0, we have c2 = Ii ai , where  Ii = 2 m √ i 0

i=1

Applying Lemma 2.2, we have

c0 (δ) =



7

1∗ i=0

2i

ai x ydx = 2

7  i=0



√ 3 3

√ − 3

the formulas of Ki are listed in Appendix A

0−2A(x)

 ai x

2i

dx.

9 − 2A(x)dx = Ki ai , 4 7



i=0

(3.6)

626

X. Sun / Applied Mathematics and Computation 270 (2015) 620–630

In order to get c1 (δ ), c3 (δ ) and c4 (δ ), we introduce the transformation u := x +

√ 3, v := y, then system (1.3) becomes

√ √ √ u˙ = v, v˙ = −u7 − 7 u6 3 + 56 u5 − 70 u4 3 + 120 u3 − 24 u2 3 + ε q(u, v)v

 (u, v) = with Hamiltonian function H



q(u, v) =

7



ai xi y

1 8 8u





3u7 +

28 3

(3.7)

√ √ u6 − 14 3u5 + 30 u4 − 8 3u3 + 12 v2 −

9 8

and,

|x=u−√3,y=v .

i=0

Applying the formula for c1 (δ ), c2 (δ ), c4 (δ ) and c5 (δ ), we have

√ 2 5 3 6 (a0 + 3a1 + 9a2 + 27a3 + 81a4 + 243a5 + 729a6 + 2187a7 ), c1 (δ) = − 3 √ 2√ 6 3(5a0 + 3a1 − 27a2 − 189a3 − 891a4 − 3645a5 − 13, 851a6 − 50, 301a7 ), c3 (δ) = 12 √ 2 5 3 6 (553a0 + 471a1 + 153a2 + 135a3 − 567a4 − 35, 721a5 − 395, 847a6 − 2, 893, 401a7 ). c4 (δ) = 1728

(3.8)

Let c1 = 0, e.i. a0 = −81a4 − 3a1 − 9a2 − 27a3 − 243a5 − 729a6 , we get, by Lemma 2.2,

c2 (δ) =



7

1∗ i=0

  − √3 7 7 3 ai x2i ai x dt = dx = 2 √ y − 3 2i

1∗ i=0

i=0



ai x2i − 94 − 2A(x)

dx =

7

Ji ai .

(3.9)

j=1

the formulas for Ji are listed in Appendix A. For the expansion of M1 (h, δ ) near the center 1c ( − 1, 0), we have

M1 (h, δ) =



bi hi+1 , f or 0 < h +

i=0

43 24



 1.

(3.10)

Developing the programs in [32], we can give the formulas of b j , j = 0, 1, 2, 3 as follows:

b0 (δ) = b1 (δ) = b2 (δ) =

b3 (δ) =

√ 2 π (a0 + a1 + a2 + a3 + a4 + a5 + a6 + a7 ), 2 √ 2 π (25a0 + 49a1 + 121a2 + 241a3 + 409a4 + 625a5 + 889a6 + 1201a7 ), 384 √ 2π (7105a0 + 16, 945a1 + 43, 009a2 + 101, 425a3 + 222, 145a4 + 448, 945a5 + 839, 425a6 221, 184 + 1, 465, 009a7 ), √ 5 2 π (1, 235, 045a0 + 3, 188, 045a1 + 8, 248, 709a2 + 20, 391, 437a3 + 48, 111, 077a4 254, 803, 968 + 108, 072, 461a5 + 230, 750, 597a6 + 468, 061, 517a7 ).

(3.11)

Last, taking h = 1000, we can obtain that

M3 (1000, δ) =

7

wi ai ,

(3.12)

i=0

 ai x2i y dx. where wi =  1000 4. Proof of main results Now we use the coefficients in (3.4), (3.6), (3.8), (3.9) and (3.11) to find limit cycles. Firstly, taking a5 = a6 = a7 = 0, and δ = (a0 , a1 , a2 ) then system (1.3) becomes of type (7, 6). Noting



det

∂(c0 , c1 , c2 ) ∂(a0 , a1 , a2 )

K0 √ 2 5  36 ×  1 =− 3 0 ≈ 6.450 = 0,

K1 3 J1

 √  64 2 5 9  = 3 6 (365s2 − 11, 458sr + 25, 893r2 ) 229, 635 J2 

K2 

X. Sun / Applied Mathematics and Computation 270 (2015) 620–630

627

∂(c ,c ,c )

which implies that rank ∂(a 0 ,a1 ,a2 ) = 3, here s and r are given in Appendix A. Solving the equation c0 = c1 = c2 = 0 gives 0

1

2





−106, 609, 270sr + 4, 436, 575s2 + 477, 275, 911r2 a4 1 4015s2 − 94, 438sr + 385, 063r2 1 a0 = − a − ≡ a∗0 , 3 11 365s2 − 11, 458sr + 25, 893r2 2145 365s2 − 11, 458sr + 25, 893r2 a1 =

1 115, 915s2 − 2, 514, 958sr + 6, 081, 523r2 a3 33 365s2 − 11, 458sr + 25, 893r2





−2, 596, 784, 560sr + 121, 603, 000s2 + 6, 471, 945, 592r2 a4 1 − ≡ a∗1 , 6435 365s2 − 11, 458sr + 25, 893r2 a2 = −

1 73, 435s2 − 1, 941, 182sr + 4, 462, 227r2 a3 33 365s2 − 11, 458sr + 25, 893r2





−1, 493, 648, 500sr + 60, 194, 450s2 + 3, 497, 816, 322r2 a4 1 − ≡ a∗2 . 6435 365s2 − 11, 458sr + 25, 893r2

(4.1)

Let δ = (a0 , a1 , a2 ) and δ ∗ = (a∗0 , a∗1 , a∗2 ), numeral computation gives:

b0 (δ ∗ ) ≈ 0.750a3 + 4.540a4 , c3 (δ ∗ ) ≈ 1.859a3 + 11.269a4 , e0 (δ ∗ ) ≈ −6.728a3 − 42.141a4 .

(4.2)

M3 (1000, δ ∗ ) ≈ 35011.409 + 6.152 × 105 a4 .

(4.3)

and

Now let a3 = 1, a4 = 0, and δ01 = (a∗0 , a∗1 , a∗2 )|a3 =1,a4 =0 , then system (1.3) becomes of type (7, 6). We take 1 , 2 , 3 , 4 and

5 positive and sufficiently small, and h1 = 0 + 1 , h2 = − 89 − 3 , h3 = − 89 + 3 , h4 = 0 − 4 , h5 = 0 + 5 , Then the first term determine the sign of the expansions of the Melnikov functions when h = hi , (i = 1, 2, . . . , 5). It is obvious that by (4.2) and (4.3),

M1 (h1 , δ01 ) = b0 (δ01 ) h1 +

 

43 24

M1 (h3 , δ01 ) = B∗10 c3 (δ01 )h3 +



 

+ h.o.t. > 0, M1 (h2 , δ01 ) = B10 c3 (δ01 )h2 +

7

7

96  + h.o.t > 0, 8

96  + h.o.t < 0, M1 (h4 , δ01 ) = c0 (δ01 ) + h.o.t < 0, 8

M3 (h5 , δ01 ) = 2c0 (δ01 ) + h.o.t < 0, M3 (1000, δ01 ) ≈ 35011.409 > 0, 1−sgn(M (h ,δ 1 )M (h ,δ 1 ))

1−sgn(M1 (h3 ,δ01 )M1 (h4 ,δ01 )) 2

1 1 0 1 2 0 which gives = 0, 2 Theorem 2.1 we have the following proposition.

= 0 and

1−sgn(M3 (h5 ,δ01 )M3 (1000,δ01 )) 2

= 1. Therefore, applying

´ der Pol system of type (7, 6), Proposition 4.1. Taking a4 = a5 = a6 = a7 = 0 and fixing a3 = 1, then system (1.3) is a Lienard–Van then there exists some (a0 , a1 , a2 ) near δ01 such that system (1.3) has 11 limit cycles, 5 limit cycles of which are near each of cuspidal loop 1∗ and 2∗ , another limit cycle is surrounding the before double homoclinic loop  . Now let a3 = 10, a4 = −1.61, and δ02 = (a∗0 , a∗1 , a∗2 )|a3 =10,a4 =−1.61 , then system (1.3) becomes of type (7, 8), numeral computation gives

b0 (δ02 ) = 0.239, c3 (δ02 ) = 0.558, c0 (δ02 ) = 0.144, M3 (1000, δ02 ) = −6.342 × 105 .

(4.4)

We take some other 1 , 2 , 3 , 4 and 5 positive and sufficiently small, and h1 = 0 + 1 , h2 = − 89 − 3 , h3 = − 98 + 3 , h4 = 0 − 4 , h5 = 0 + 5 , by (4.4), we have

7  

43 96 2 2  M1 (h1 , δ ) = b0 (δ ) h1 + + h.o.t. > 0, M1 (h2 , δ0 ) = B10 c3 (δ0 )h2 +  + h.o.t > 0, 24 8  76  9  M1 (h3 , δ02 ) = B∗10 c3 (δ02 )h3 +  + h.o.t < 0, M1 (h4 , δ02 ) = c0 (δ02 ) + h.o.t > 0, 2 0

2 0

8

M3 (h5 , δ02 ) = 2c0 (δ02 ) + h.o.t > 0, M3 (1000, δ02 ) ≈ −6.342 × 105 < 0, 1−sgn(M (h ,δ 2 )M (h ,δ 2 ))

1−sgn(M1 (h3 ,δ02 )M1 (h4 ,δ02 )) 2

1 1 0 1 2 0 which gives = 0, 2 Theorem 2.1 we have the following proposition.

= 1 and

1−sgn(M3 (h5 ,δ02 )M3 (1000,δ02 )) 2

= 1. Therefore, applying

Proposition 4.2. Taking a3 = 10, a4 = −1.6, a5 = a6 = a7 = 0, system (1.3) is of type (7, 8), then there exists some (a0 , a1 , a2 ) near

δ02 such that system (1.3) has 13 limit cycles, 5 limit cycles of which are near each of cuspidal loop 1∗ and 2∗ , 1 limit cycles between i∗ and  i for i = 1, 2. Another limit cycle is surrounding the before double homoclinic loop  .

628

X. Sun / Applied Mathematics and Computation 270 (2015) 620–630

Next, we fix a4 = 1, a6 = a7 = 0 and δ = (a0 , a1 , a2 , a3 , a5 ) then system (1.3) becomes of type (7, 10). Noting

  K0   1   3  0 √ √ ∂(c0 , c1 , c2 , c3 , b0 ) 5 1  det = − 2 × 36 × 2 × 36 ×  ∂(a0 , a1 , a2 , a3 , a5 )  5  12  √  2π 

K1

K2

K3

1

3

9

J1 1 4 √ 2π 2

J2 9 − 4 √ 2π 2

J3 63 − 4 √ 2π 2

    81    J5   1215  −  4  √  2π   K5

2 2 √   256 2π =− 33, 663, 388, 925s2 − 570, 458, 745, 710sr + 1, 139, 908, 283, 297r2 31, 031, 725, 725 ≈ −2491.684 = 0,

∂(c ,c ,c ,c ,b )

which implies that rank ∂(a 0,a1 ,a2 ,a3 ,a0 ) = 4 + 1, Solving the equation c0 = c1 = c2 = c3 = b0 = 0 gives 0

a0 = a1 = a2 = a3 = a5 =

1

2

3

5

−16, 388, 666rs + 31, 741, 073r2 + 918, 425s2 15 ≡ a#0 , −1, 664, 080, 670rs + 2, 818, 155, 179r2 + 186, 433, 475s2 −3, 667, 745, 270rs + 7, 154, 749, 847r2 + 226, 014, 575s2 − ≡ a#1 , −1, 664, 080, 670rs + 2, 818, 155, 179r2 + 186, 433, 475s2 2 −10, 943, 022, 650rs + 20, 530, 420, 673r2 + 839, 746, 025s2 ≡ a#2 , 3 −1, 664, 080, 670rs + 2, 818, 155, 179r2 + 186, 433, 475s2 2 14, 287, 270, 109r2 − 8, 027, 254, 130rs + 765, 923, 525s2 − ≡ a#3 , 3 −1, 664, 080, 670rs + 2, 818, 155, 179r2 + 186, 433, 475s2 −8, 857, 670rs + 14, 362, 943r2 + 1, 114, 775s2 −21 ≡ a#5 . −1, 664, 080, 670rs + 2, 818, 155, 179r2 + 186, 433, 475s2

Take δ03 = (a# , a# , a# , a# , a# ), then we have 0 1 2 3 5

b1 (δ03 ) ≈ 0.064 > 0, c4 (δ03 ) ≈ −6.857 < 0, c0 (δ03 ) ≈ −1.596 < 0

(4.5)

and

M3 (1000, δ03 ) ≈ 5.490 × 105 > 0 We take some other 1 , 2 , 3 , 4 and 5 positive and sufficiently small, and h1 = 0 + 1 , h2 = − 89 − 3 , h3 = − 98 + 3 , h4 = 0 − 4 , h5 = 0 + 5 , by (4.5), we have

M1 (h1 , δ03 ) = b1 (δ03 ) h1 + M1 (h3 , δ03 ) =

43 24

2

+ h.o.t. > 0, M1 (h2 , δ03 ) = −

 11





 11

1 9 6  B00 c4 (δ03 )h2 +  + h.o.t > 0, 11 8

1 ∗ 9 6  B00 c4 (δ03 )h3 +  + h.o.t > 0, M1 (h4 , δ03 ) = c0 (δ03 ) + h.o.t < 0, 11 8

M3 (h5 , δ03 ) = 2c0 (δ03 ) + h.o.t < 0, M3 (1000, δ03 ) ≈ 5.490 × 105 > 0, 1−sgn(M (h ,δ 3 )M (h ,δ03 ))

1 1 0 1 2 which gives 2 have the following proposition.

= 0,

1−sgn(M1 (h3 ,δ03 )M1 (h4 ,δ03 )) 2

= 1 and

1−sgn(M3 (h5 ,δ03 )M3 (1000,δ03 )) 2

= 1. Applying Theorem 2.1 we

Proposition 4.3. Taking a4 = 1, a6 = a7 = 0, system (1.3) is of type (7, 10), then there exists some (a0 , a1 , a2 , a3 , a5 ) near δ03 such that system (1.3) has 17 limit cycles, 6 limit cycles of which are near each of cuspidal loop 1∗ and 2∗ , 1 limit cycle is near each of the center, 1 limit cycle between i∗ and  i for i = 1, 2. Another limit cycle is surrounding the before double homoclinic loop  . Last, we fix a7 = 1 and δ = (a0 , a1 , a2 , a3 , a4 , a5 , a6 ) then system (1.3) becomes of type (7, 14). Noting

det

∂(c0 , . . . , c3 , c0 , c1 , b0 ) ≈ 1.916 × 105 = 0, ∂(a0 , a1 , . . . , a6 ) ∂(c ,c ,c ,c ,c ,c ,b )

which implies that rank ∂(a 0,a 1,a2 ,a3 ,a0 ,a1 ,a0 ) = 2 + 4 + 1. 0

1

2

3

4

5

6

Solving the equation c0 = c1 = c2 = c3 = c0 = c1 = b0 = 0 gives

a0 = a&0 = 0, a1 = a&1 ≈ 596.560, a2 = a&2 ≈ −1319.649,

X. Sun / Applied Mathematics and Computation 270 (2015) 620–630

629

a3 = a&3 ≈ 1093.059, a4 = a&4 ≈ −472.966, a5 = a&5 ≈ 118.636, a6 = a&6 ≈ −16.640, Take δ04 = (a& , a& , a& , a& , a& , a& , a& ), then we have 0 1 2 3 4 5 6

b1 (δ04 ) ≈ 1.082 > 0, c4 (δ04 ) ≈ −100.904 < 0, c2 (δ04 ) ≈ −16.097 < 0

(4.6)

and

M3 (1000, δ04 ) ≈ 1.184 × 108 > 0. We take some other 1 , 2 , 3 , 4 and 5 positive and sufficiently small, and h1 = 0 + 1 , h2 = − 89 − 3 , h3 = − 98 + 3 , h4 = 0 − 4 , h5 = 0 + 5 , by (4.6), we have

M1 (h1 , δ04 ) = b1 (δ04 ) h1 + M1 (h3 , δ04 ) =

43 24



2

+ h.o.t. > 0, M1 (h2 , δ04 ) = −

 11



 11

1 9 6  B00 c4 (δ04 )h2 +  + h.o.t > 0, 11 8

1 ∗ 9 6  B00 c4 (δ04 )h3 +  + h.o.t > 0, M1 (h4 , δ04 ) = c2 (δ04 )h4 + h.o.t > 0, 11 8

M3 (h5 , δ04 ) = 2c2 (δ04 )h5 + h.o.t < 0, M3 (1000, δ04 ) ≈ 5.490 × 105 > 0, 1−sgn(M (h ,δ 4 )M (h ,δ 4 ))

1−sgn(M (h ,δ 4 )M1 (h4 ,δ04 ))

1 1 0 1 2 0 1 3 0 which gives = 0, 2 2 By Theorem 2.1 we have the following proposition.

= 0 and

1−sgn(M3 (h5 ,δ04 )M3 (1000,δ04 )) 2

= 1.

Proposition 4.4. Fixing a7 = 1, system (1.3) is of type (7, 14), then there exists some (a0 , a1 , a2 , a3 , a4 , a5 , a6 ) near δ04 such that system (1.3) has 20 limit cycles, 6 limit cycles of which are near each of cuspidal loop 1∗ and 2∗ , 2 limit cycles near each of the homoclinic loop  1 and  2 , 1 limit cycle is near each of the centers, and another two limit cycles are surrounding the before double homoclinic loop  . Summarizing above Propositions 4.1–4.4, we have proved Theorem A. Remark 4.1. In this section, we have used numerous computations only for determining the signs of the relevant terms not for their approximate values, the signs of relative terms play an important role in the method. By the coefficients of the expansions of related Melnikov functions and some bifurcation theory, limit cycles of the Liénard–Van der Pol system (1.3) of type (7, 6), (7, 8), (7, 10) and (7, 14) are found, respectively. However we do not investigate the case of type (7, 12), for which case we cannot obtain more than 17 limit cycles. That is because the perturbation term x12 y is equivalent to some linear combination of x10 y and x8 y, which is found by computation. Therefore it plays no role in our method, and 17 limit cycles can be obtained by taking the case of type (7, 10) as the special case of (7, 12). Appendix A In this section, we list some expression of the coefficients Ki and Ji . Let s = coefficients Ki and Ji are as follows:

√ √ √ √ 3EllipticK ( 2 3 2 ), and r = 3EllipticE ( 2 3 2 ), the

304 544 6704 371, 584 16 64 80 29, 536 s+ r, K1 = − s+ r, K2 = − s+ r, K3 = − s+ r, 27 135 567 405 729 3645 168, 399 120, 285 27, 151, 408 164, 341, 088 164, 341, 088 2, 144, 752 12, 971, 648 1, 408, 816, 144 s+ r, K5 = − s+ r, K6 = − s+ r, − 6, 567, 561 4, 691, 115 19, 702, 683 14, 073, 345 1, 004, 836, 833 14, 073, 345 2, 258, 163, 462, 592 178, 220, 608, 864 s+ r. − 57, 275, 699, 481 40, 911, 213, 915 188 11, 336 41, 368 4 104 904 18, 404 s− r, J4 = − s− r, J5 = − s− r, − s, J2 = −4 s − 4r, J3 = − 3 9 9 27 135 189 135 3, 879, 212 431, 516, 948 207, 004 139, 706, 992 s− r, J7 = − s− r. − 729 3645 168, 399 120, 285

K0 = − K4 = K7 = J1 = J6 = References [1] [2] [3] [4] [5] [6] [7]

J. Nagumo, S. Arimoto, S. Yoshizawa, An active pulse transmission line simulating nerve axons, Proc. IRL 50 (1960) 2061–2070. T. Pham Dinh, J. Demongeot, P. Baconnier, G. Benchetrit, Simulation of a biological oscillator: the respiratory rhythm, J. Theor. Biol. 103 (1983) 113–132. M. Han, Bifurcation Theory of Limit Cycles, Science Press, Beijing, 2013. T. Carlettia, G. Villarib, A note on existence and uniqueness of limit cycles for Liénard systems, J. Math. Anal. Appl. 307 (2005) 763–773. D. Xiao, Z. Zhang, On the existence and uniqueness of limit cycles for generalized Liénard systems, J. Math. Anal. Appl. 343 (2008) 299–309. L. Cherkas, V. Romanovki, The center conditions for a Liénard system, Comput. Math. Appl. 52 (2006) 363–374. C. Liu, G. Chen, J. Yang, On the hyperelliptic limit cycles of Liénard systems, Nonlinearity 25 (2012) 1601.

630 [8] [9] [10] [11] [12] [13] [14] [15] [16] [17] [18] [19] [20] [21] [22] [23] [24] [25] [26] [27] [28] [29] [30] [31] [32]

X. Sun / Applied Mathematics and Computation 270 (2015) 620–630 X. Yu, X. Zhang, The hyperelliptic limit cycles of the Liénard systems, J. Math. Anal. Appl. 376 (2011) 535–539. M. Sabatini, On the period function of Liénard systems, J. Differ. Equat. 152 (1999) 467–487. Y. Tian, M. Han, Hopf bifurcation for two types of Liénard systems, J. Differ. Equ. 251 (2011) 834–859. M. Han, Y. Tian, P. Yu, Small-amplitude limit cycles of polynomial Liénard systems, Sci. China Math. 56 (2013) 1543–1556. P. Yu, Local and global bifurcations to limit cycles in a class of Liénard equation, Int. J. Bifurc. Chaos 17 (2007) 183–198. M. Caubergh, J. Francoise, Generalized Liénard equations, cyclicity and Hopf–Takens bifurcations, Qual. Theory Dyn. Syst. 5 (2004) 195–222. J. Yang, On the limit cycles of a kind of Liénard system with a nilpotent center under perturbations, J. Anal. Appl. Comput. 2 (2012) 325–339. B. Coll, F. Dumortier, R. Prohens, Alien limit cycles in Liénard equations, J. Differ. Equ. 254 (2013) 1582–1600. C. Li, Abelian integrals and limit cycles, Qual. Theory Dyn. Syst. 11 (2012) 111–128. S. Smale, Mathematical problems for the next century, Math. Intell. 20 (1998) 7–15. M. Caubergh, Hilbert’s sixteenth problem for polynomial Liénard equations, Qual. Theory Dyn. Syst. 11 (2012) 3–18. Y. Xiong, M. Han, New lower bounds for the Hilbert number of polynomial systems of Liénard type, J. Differ. Equ. 257 (2014) 2565–2590. J. Yang, M. Han, V. Romanovski, Limit cycle bifurcations of some Liénard systems, J. Math. Anal. Appl. 366 (2010) 242–255. M. Han, H. Yan, J. Yang, C. Lhotka, On the number of limit cycles of some Liénard systems, Can. Appl. Math. Q. 17 (2009) 61–83. A. Atabaigi, R. Zangeneh, Bifurcation of limit cycles in small perturbations of a class of hyper-elliptic Hamiltonian systems of degree 5 with a cusp, J. Appl. Anal. Comput. 1 (2011) 299–313. R. Kazemi, R. Zangeneh, Bifurcation of limit cycles in small perturbations of a hyper-elliptic Hamiltonian system with two nilpotent saddles, J. Anal. Appl. Comput. 2 (2012) 395–413. X. Sun, H. Xi, R. Zangeneh, R. Kazemi, Bifurcation of limit cycles in small perturbation of a class of Liénard systems, Int. J. Bifurc. Chaos 24 (2014) 1450004. L. Zhao, The perturbations of a class of hyper-elliptic Hamilton systems with a double homoclinic loop through a nilpotent saddle, Nonlinear Anal. 95 (2014) 374–387. C. Li, C. Liu, J. Yang, A cubic system with thirteen limit cycles, J. Differ. Equ. 246 (2009) 3609–3619. P. Yu, M. Han, Limit cycles in generalized Liénard systems, Chaos Solitons Fractals 30 (2006) 1048–1068. M. Han, V. Romanovski, On the number of limit cycles of polynomial Liénard systems, Nonlinear Anal.: Real World Appl. 14 (2013) 1655–1668. M. Han, J. Yang, A. Tarta, Y. Gao, Limit cycles near homoclinic and heteroclinic loops, J. Dyn. Differ. Equ. 20 (2008) 923–944. M. Han, H. Zang, J. Yang, limit cycle bifurcations by perturbing a suspidal loop in a Hamiltonian system, J. Differ. Equ. 246 (2009) 129–163. J. Yang, M. Han, Limit cycle bifurcations of some Liénard systems with a nilpotent cusp, Int. J. Bifurc. Chaos 20 (2010) 3829–3839. M. Han, J. Yang, P. Yu, Hopf bifurcations for near-Hamiltonian systems, Int. J. Bifurc. Chaos 19 (2009) 4117–4130.