Moving discrete breathers in nonlinear lattice: Resonance and stability

Available online at www.sciencedirect.com

Wave Motion 45 (2007) 83–99 www.elsevier.com/locate/wavemoti

Moving discrete breathers in nonlinear lattice: Resonance and stability Kazuyuki Yoshimura a

a,*

, Yusuke Doi

b

NTT Communication Science Laboratories, NTT Corporation, 2-4, Hikaridai, Seika-cho, Soraku-gun, Kyoto 619-0237, Japan b Department of Adaptive Machine Systems, Graduate School of Engineering, Osaka University, 2-1, Yamadaoka, Suita, Osaka 565-0871, Japan Received 20 October 2006; received in revised form 1 February 2007; accepted 12 April 2007 Available online 25 April 2007

Abstract We have studied moving discrete breathers in the Fermi–Pasta–Ulam b lattice. The moving discrete breathers are not strictly localized but have spatially extended tails. It is shown that a resonance phenomenon in the tail amplitude occurs: the tail amplitude is enhanced for some particular values of the breather frequency and these values depend on the velocity. This resonance is caused by strong excitation of a particular set of the harmonic normal modes. We analytically derive a condition for this resonance. Moreover, we analyze the linear stability of moving discrete breathers and show that the instability appears when the breather frequency is close to the resonance points.  2007 Elsevier B.V. All rights reserved. Keywords: Fermi–Pasta–Ulam lattice; Moving discrete breather; Spatially extended tail; Resonance; Stability

1. Introduction Spatially localized excitations in nonlinear lattices have attracted great interest since the ground-breaking work of Takeno et al. [1]. The localized modes are called discrete breathers (DBs) or intrinsic localized modes (ILMs). DBs are time-periodic and spatially localized solutions of the equations of motion, which emerge due to nonlinearity and discreteness of the systems. DBs have the following characteristics: (i) they exist in nonintegrable lattices, (ii) they exist not only in one-dimensional lattices but also in high-dimensional lattices, and (iii) they exist in both of Hamiltonian and dissipative lattices. Because of these characteristics, existence of the DB solutions is quite general in nonlinear space-discrete dynamical systems. Stationary DBs, which do not propagate along the lattices, have been considered so far. Considerable numerical and theoretical progress has been achieved in understanding the nature of DB ([2] and references therein). Experimental evidence for the existence of DB has been also reported in various systems as diverse as a Josephson junction array [3], an optical waveguide array [4], and a micromechamical system [5]. *

Corresponding author. Tel.: +81 774 5136; fax: +81 774 93 5158. E-mail address: [email protected] (K. Yoshimura).

0165-2125/$ - see front matter  2007 Elsevier B.V. All rights reserved. doi:10.1016/j.wavemoti.2007.04.004

84

K. Yoshimura, Y. Doi / Wave Motion 45 (2007) 83–99

A generalization of the concept of DB is a moving DB that propagates along the lattice. It has been numerically shown that there exist moving localized modes, which propagate along the lattices and do not decay for a long time, in Fermi–Pasta–Ulam (FPU) [6] and nonlinear Klein–Gordon lattice models [7]. The FPU lattice has been extensively studied so far to understand the relaxation process from a nonequilibrium initial state to the equilibrium state since the work of Fermi et al. [8]. It has been also numerically found that moving localized modes emerge due to the modulational instability of the initially excited mode and they play an important role in the relaxation process [9–12]. As for the modulational instability, detailed studies can be found in [13– 16]. The moving localized modes found in these numerical calculations are not radiationless but emit small ripples. However, the above numerical results strongly suggest the existence of a radiationless moving localized mode, which we call a moving DB, although the mathematically rigorous existence proof has not yet been achieved. It is an interesting topic in the field of nonlinear lattice dynamics to search for the moving DBs and clarify their nature. Precise numerical calculations of the moving DBs have been carried out for nonlinear Klein–Gordon [17] and Salerno lattices [18] by using the Newton method. For the FPU lattice, no such calculation of the moving DBs has been done. These calculations have revealed that the profiles of the moving DBs are quite different from those of the stationary one: the moving DBs are not well localized but have spatially extended tails. In Ref. [18], a role of the spatially extended tail has been pointed out: the energy of the core of moving DB varies due to the interaction with its tail during the propagation by one lattice spacing; this interaction between the DB core and the tail is necessary for supporting the radiationless propagation. Two basic parameters of moving DB are the frequency and the velocity. One of the fundamental properties of moving DB is the dependence of its profile and stability on these two parameters. However, this fundamental property has not yet been clarified. There has been no extensive numerical calculation of the moving DB for various values of the frequency and velocity. Therefore, we carried out extensive precise numerical calculation of the moving DB in the FPU-b lattice, which has quadratic and quartic nearest-neighbor interaction potentials, in order to clarify how the profile and stability of the moving DB depends on its frequency and velocity. In the present paper, we show that a kind of resonance phenomenon occurs in the tail amplitude: the profile of the moving DB changes sensitively depending on the two parameters. Moreover, we show that the stability of the moving DB is significantly influenced by this resonance phenomenon. The present paper is organized as follows. In Section 2, we describe the FPU-b model and introduce the normal mode coordinates. In Section 3, we define the radiationless moving localized mode, i.e., the moving DB solution with rational velocity, and describe our numerical method to calculate it. In Section 4, we show numerical results and discuss the resonance phenomenon and the stability of moving DB. In Section 5, we develop a theory to describe the resonance phenomenon. A condition for the resonance is analytically derived. We make a comparison of this theoretical resonance condition with the numerical results and show their good agreement. Finally, conclusions are drawn in Section 6. 2. FPU-b model and normal mode We describe the FPU-b model and the normal modes. Our investigation is of the dynamical model described by the Hamiltonian  N N  X 1X 1 b 2 4 2 H¼ ðq  qn Þ þ ðqnþ1  qn Þ : p þ ð1Þ 2 n¼1 n n¼1 2 nþ1 4 This is referred to as the FPU-b model. This Hamiltonian describes a one-dimensional oscillator chain with nonlinear nearest-neighbor interaction and has been proved nonintegrable [19]. We employ the periodic boundary condition, i.e., qn+N = qn. The parameter b represents the nonlinear coupling strength. We set b = 1 in the numerical experiments presented later. We assume the system size N is even. The equations of motion derived from Hamiltonian (1) are h i d2 qn 3 3 ¼ q þ q  2q þ b ðq  q Þ  ðq  q Þ : nþ1 n1 n nþ1 n n n1 dt2

ð2Þ

K. Yoshimura, Y. Doi / Wave Motion 45 (2007) 83–99

85

We define the energy en of nth lattice site by the sum of the kinetic energy of the nth particle and the half of the interaction potential energies shared with its two nearest neighbors. The site energy en, n = 1, 2, . . ., N is given by i bh i 1 1h en ¼ p2n þ ðqnþ1  qn Þ2 þ ðqn  qn1 Þ2 þ ðqnþ1  qn Þ4 þ ðqn  qn1 Þ4 : ð3Þ 2 4 8 The orthonormal transformation q = (q1, . . ., qN) # Q = (Q(N/2-1), . . ., QN/2) defined by      N =2 n X ð1Þ 2p 2p mn  sin mn ; n ¼ 1; 2; . . . ; N ; Qm cos qn ¼ pffiffiffiffi N N N m¼ðN =21Þ

ð4Þ

gives the normal modes of the corresponding linear system [14]. Here, Qm is the amplitude of the mth normal mode. The natural frequency of the mth normal mode is given by pm xm ¼ 2 cos : ð5Þ N Note that the above definition of mode number is different from the usual one. In our definition, m = 0 mode corresponds to the zone boundary mode, which has the shortest wavelength and the highest natural frequency. _ Hamiltonian (1) is In terms of the normal mode coordinates Q and their conjugate momenta Pð¼ QÞ, rewritten as Nh 2

1 1 X b H ¼ P N =2 þ P þ x2m Q2m þ 2 2 m¼N h m 8N

Nh X

xi xj xk xl Qi Qj Qk Ql Dði; j; k; lÞ;

ð6Þ

i;j;k;l¼N h

where Nh = N/2  1 and D(i, j, k, l) represents the selection rule defining the interaction among the normal modes. It is given by Dði; j; k; lÞ ¼ Dði þ j þ k þ lÞ þ Dði þ j  k  lÞ þ Dði  j þ k  lÞ þ Dði  j  k þ lÞ; where the function D is defined by  m ð1Þ if r ¼ mN ðm 2 ZÞ; DðrÞ ¼ 0 otherwise:

ð7Þ

ð8Þ

In Hamiltonian (6), the variables for the N/2th mode, which represents the translational motion, are separated. Therefore, we can assume QN/2 = PN/2 = 0 and eliminate them in the equations of motion. The equation of motion for the mth normal mode is d2 Qm b þ x2m Qm þ 2 2N dt

Nh X

xm xi xj xk Qi Qj Qk Dðm; i; j; kÞ ¼ 0;

ð9Þ

i;j;k¼N h

where m = (N/2  1), . . ., N/2  1. Eq. (5) indicates that ±mth normal modes have the same natural frequency, i.e., xm = xm. There is no essential difference between these two normal modes. We define a modal energy Em by the sum of harmonic energies of these two modes. The modal energy Em, m = 0, 1, . . ., N/2  1 is given by

1

1 ð10Þ Em ¼ P 2m þ x2m Q2m þ P 2m þ x2m Q2m : 2 2 We will use Em to investigate the resonance phenomenon. 3. Moving DB and numerical method It is well known that there are two types of stationary DBs, which have different symmetries in their displacement patterns. One is called Sievers–Takeno mode or odd mode [1] and the other is called Page mode or even mode [20]. Figs. 1(a) and (b) show the profiles of the odd mode and the even mode, respectively. The DB center is located at a lattice site in the case of odd mode while it is located at a midpoint between

86

K. Yoshimura, Y. Doi / Wave Motion 45 (2007) 83–99

Fig. 1. Profile of stationary DB: (a) odd mode, and (b) even mode.

two lattice sites in the case of even mode. In both modes, each particle oscillates out of phase with its nearest neighbors and with the same frequency, which we denote by x. We call this oscillation the internal oscillation of DB. Figs. 1(a) and (b) illustrate the strongly localized profiles, which correspond to the case of large x. The localization of DB becomes weak as x decreases. It is well known that the particle amplitude decreases exponentially with increasing the distance from the DB center. It has been shown that the odd mode is linearly unstable while the even mode is linearly stable in the FPU-b lattice [21]. In this section, we define the moving DB with a rational velocity and describe an outline of our numerical method. The moving DB propagates from a lattice site to the next one, oscillating with the internal frequency x. The moving DB successively takes odd mode like, intermediate, and even mode like shapes during its propagation. We define the velocity V of the moving DB by the number of lattice spacings that it propagates during one internal oscillation period T  2p/x. Hereafter, we consider the case that the velocity V is a rational number: i.e., V¼

r s

ðsite=periodÞ;

ð11Þ

where r and s are integers. Eq. (11) means that the moving DB propagates by r lattice spacings during s internal oscillation periods sT. We call V simply the velocity although, precisely speaking, it should be called the periodwise velocity since V is defined for one internal oscillation period. At present, no numerical method is known for calculating the moving DB with an arbitrary velocity. However, it is possible to obtain a precise numerical solution for the moving DB with a rational velocity by using the Newton method when a good approximation for the moving DB can be obtained. Given x and V = r/s, the procedure for calculating the moving DB is as follows: Step 1 Calculate a precise numerical solution of the stationary DB with the given frequency x. Step 2 Construct an approximate solution for the moving DB with the given frequency x and the given velocity V by adding an appropriate perturbation to the stationary DB solution obtained in Step 1. Step 3 Use the solution constructed in Step 2 as an initial approximation and calculate a precise numerical solution of the moving DB with x and V by using the Newton method. Calculation of the stationary DB in Step 1 can be performed in several ways. For example, one can use a method proposed in [22], in which the stationary DB solution with a given x is obtained by continuation from the anti-continuous limit of a certain extended lattice model. An analytical approximate solution for the stationary DB has been obtained [23]. Therefore, it is also possible to use this solution as an initial approximation and numerically calculate a precise stationary DB solution by using the Newton method. It turns out that the domain of convergence of the Newton method in Step 3 is very small. Therefore, it is important to construct a precise approximate solution in Step 2. We have developed a simple method to determine the perturbation for constructing a good approximate solution. The detail of this method will be presented elsewhere. This method works successfully at least in the FPU-b lattice and our numerical results shown in the next section have been obtained by using this method. Chen et al. have proposed a method to determine a perturbation to obtain an approximate moving DB solution [7]. Their method numerically computes an eigenvector, which is called the pinning mode, of the Floquet matrix associated with the station-

K. Yoshimura, Y. Doi / Wave Motion 45 (2007) 83–99

87

ary DB solution. Compared with their method, our method is simple and useful because our method does not need numerical computation of the Floquet eigenvectors. We describe the details of calculation in Step 3. Let (q1(0), . . ., qN(0),p1(0), . . ., pN(0)) be a set of initial conditions. We denote it by (q(0), p(0)). If we integrate the equations of motion (2) over the period sT = 2p s/x with these initial conditions, then we obtain (q(sT),p(sT)). We define the time evolution map F : R2N ! R2N as follows: Fðq1 ð0Þ; . . . ; qN ð0Þ; p1 ð0Þ; . . . ; pN ð0ÞÞ ¼ ðq1 ðsT Þ; . . . ; qN ðsT Þ; p1 ðsT Þ; . . . ; pN ðsT ÞÞ:

ð12Þ

The initial conditions of the moving DB with given x and V = r/s have to satisfy the following condition r

Fðq1 ð0Þ; . . . ; qN ð0Þ; p1 ð0Þ; . . . ; pN ð0ÞÞ ¼ ð1Þ  ðq1r ð0Þ; . . . ; qNr ð0Þ; p1r ð0Þ; . . . ; pNr ð0ÞÞ;

ð13Þ

where mod N is taken in the indices. In the case of even r, Eq. (13) requires that the displacement and momentum pattern at t = sT has to coincide with that obtained by translating the initial displacement and momentum pattern by r lattice spacings. In the case of odd r, it requires that the displacement and momentum pattern at t = sT has to coincide with that obtained by reversing the phase of each particle after the translation by r lattice spacings. Eq. (13) gives a set of equations with respect to (q(0), p(0)). Therefore, if a good approximate solution to Eq. (13) is available, the numerical solution (q(0), p(0)) can be obtained precisely by using the Newton method. In this calculation, the map F and its Jacobian matrix DF are necessary. They can be obtained by numerically integrating the equations of motion (2) and the variational equations. Eq. (13) is degenerated because of the translational invariance of Eq. (2) and arbitrariness of the initial point along the trajectory of the moving DB. In order to remove this degeneracy, we carried out the calculation under the constraints QN/2 = PN/2 = 0 and Q1 = Q1. 4. Resonance and stability of moving DB 4.1. Examples of numerical solutions Examples of numerical moving DB solutions are shown in Fig. 2. The lattice size N and the velocity V are N = 60 and V = 1/3, respectively. The frequency x of the moving DB is x = 2.5 in Figs. 2(a)–(c) and x = 2.53 in Fig. 2(d). Fig. 2(a) shows the site energy en defined by Eq. (3) as a function of site number n and time t. The spatially localized excitation propagates smoothly with a constant velocity. This result demonstrates that a precise moving DB solution can be obtained by our numerical method. The profile of the same moving DB solution at t = 0 is shown in Fig. 2(b), where qn is plotted against n. This solution takes a profile similar to the odd mode at t = 0. However, the solution does not go to zero away from the DB center but has a small spatially extended tail. This extended tail can be clearly observed in Fig. 2(c), where the site energy en at t = 0 is plotted against n by solid line in a semi-logarithmic plot: en decreases exponentially with increasing the distance from the DB center in the region close to the center while en is almost constant with a nonzero value away from the center. It is well known that for a stationary DB en decreases exponentially with increasing the distance from the center as illustrated by dashed line in Fig. 2(c). Therefore, the existence of spatially extended tail is a remarkable feature of moving DB. Spatially extended tails of moving DBs have already been found also in nonlinear Klein–Gordon [17] and Salerno lattices [18]. Fig. 2(d) shows the profile of moving DB with x = 2.53. It can be found that the amplitude of the tail is much larger than that in the case of x = 2.5, although the difference in the frequency is small. This fact indicates that the tail amplitude shows a sensitive dependence on the parameters of moving DB. It is conceivable that a change in the profile of moving DB affects its stability. These issues have not been studied so far. In the present study, we focus on the parameter dependence of the tail amplitude and the stability. 4.2. Resonance in tail amplitude Let emax be the site energy at the DB center and etail be an average site energy of the tail part (see Fig. 2(c)). We use the ratio etail/emax to measure the tail amplitude. Fig. 3 shows etail/emax plotted against the inverse of velocity 1/V, where the other parameters are N = 40 and x = 2.4. It can be found that the tail amplitude

88

K. Yoshimura, Y. Doi / Wave Motion 45 (2007) 83–99

0.6

en

0.4

1 0.8 0.6 0.4 0.2 0

0

10

20

160 140 120 100 80 60 40 20

qn

0 -0.2

t

30

40

50

n

60

0.2

-0.4

0

0

10

20

30

40

50

60

0

10

20

30

40

50

60

n

10 0.6

e max

1

0.4

0.1

en

0.01

qn

0.001 0.0001

0 -0.2

e tail

1e-05 1e-06 0

0.2

-0.4 10

20

30

40

n

50

60

n

Fig. 2. Numerical solution for moving DB (N = 60, V = 1/3): (a) spatiotemporal plot of site energy, (b) initial displacements of the particles, and (c) semi-logarithmic plot of site energy at t = 0 (solid line) and illustration for site energy of stationary DB (dashed line); breather frequency is x = 2.5 in figures (a)–(c); (d) initial displacements of the particles for x = 2.53.

0.01 0.0001

etail /emax

1e-06 1e-08 1e-10 1e-12 1e-14 5

10

15

20

25

30

35

40

1/V Fig. 3. Tail amplitude vs. velocity (N = 40, x = 2.4).

shows a tendency to decrease as the velocity decreases. This tendency is natural because it is expected that the profile of moving DB approaches that of a stationary DB. An interesting feature is that etail/emax does not show a simple dependence but strongly depends on the velocity: a complicated fluctuation in etail/emax is observed.

K. Yoshimura, Y. Doi / Wave Motion 45 (2007) 83–99

89

Figs. 4(a)–(c) show etail/emax plotted as a function of x for V = 1/3, 1/5, and 1/7, respectively. The lattice size is N = 40. A strong dependence of the tail amplitude on x is clearly observed in each figure. In Fig. 4(a), the tail amplitude has a tendency to increase as x increases. The core of moving DB becomes more strongly localized with increasing x. Therefore, it can be stated that the tail amplitude tends to increase as localization of the moving DB becomes stronger. A striking feature is that the tail amplitude is sharply enhanced for some particular values of x: i.e., a resonance phenomenon in the tail amplitude occurs. It should be emphasized that this resonance is also clearly observed in both of Figs. 4(b) and (c). Thus, it may be concluded that the resonance phenomenon in the tail amplitude is one of the remarkable characteristics of moving DB, which occurs for any rational velocity. The resonance occurs for different values of x, depending on the velocity V. In Figs. 4(b) and (c), apart from the sharp peaks due to the resonance, the tail amplitude gradually increases on an average as x increases in small x region and then shows a sharp decrease. This sharp decrease occurs at x ’ 2.5 for V = 1/5 and x ’ 2.33 for V = 1/7. After the sharp decrease, the tail amplitude tends to increase with increasing x again, exhibiting some resonance peaks. We have calculated the modal energy Em defined by Eq. (10) to understand a mechanism for the resonance. Distributions of the modal energy are compared between a resonance point x = 2.415 (filled circle) and a nonresonance point x = 2.35 (open circle) in Fig. 5, where the other parameters are N = 40 and V = 1/3. The values of N and V are the same as in Fig. 4(a). The modal energy Em is plotted against the mode number m. In the inset, positions of the two chosen points are indicated on the resonance curve. In the case of the nonresonance point, Em decreases monotonically with increasing m. In contrast, Em for a particular mode number (m = 7) is significantly enhanced in the case of the resonance point: i.e., ±7th normal modes are strongly excited since E7 is the sum of harmonic energies of these two modes. This comparison indicates that the resonance in the tail

0.01

0.01

3 0.001

8

5

0.0001

e tail /e max

0.0001 1e-05 1e-06

4

2

1

6

1e-05 1e-06 1e-07

1e-07 1e-08 1e-08 1e-09 1e-09 2.1

2.2

2.3

2.4

ω

2.5

2.6

2.7

2.1

2.8

2.2

2.3

2.4

ω

2.5

2.6

0.01

1 0.001

4

3 2

0.0001

e tail /e max

e tail /e max

0.001

5

6

7

1e-05

5

1e-06 1e-07 1e-08 1e-09 2.1

2.2

2.3

2.4

2.5

ω

2.6

2.7

2.8

2.9

Fig. 4. Tail amplitude vs. frequency (N = 40): velocity is (a) V = 1/3, (b) V = 1/5, and (c) V = 1/7.

2.7

2.8

90

K. Yoshimura, Y. Doi / Wave Motion 45 (2007) 83–99 0.3 0.01

0.001

0.25 0.0001

1e-05

0.2

1e-06

Em

1e-07

1e-08

0.15

1e-09 2.1

2.2

2.3

2.4

2.5

2.6

2.7

2.8

0.1

0.05

0 0

2

4

6

8

10

12

14

16

18

20

m Fig. 5. Modal energy distribution for resonance and nonresonance points (N = 40, V = 1/3): resonance point x = 2.415 (filled circle) and nonresonance point x = 2.35 (open circle).

amplitude is caused by strong excitation of a particular set of the normal modes. Fig. 6 shows the modal energy distribution for the four different resonance points x = 2.255, 2.415, 2.56, and 2.689, where N = 40 and V = 1/3. Positions of the chosen points are indicated in the inset. The same feature of the modal energy

0.45 0.01

0.4

0.001

0.0001

0.35 1e-05

0.3

1e-06

1e-07

0.25

Em

1e-08

1e-09

0.2

2.1

2.2

2.3

2.4

2.5

2.6

2.7

2.8

0.15 0.1 0.05 0 0

2

4

6

8

10

12

14

16

18

20

m Fig. 6. Modal energy distribution for resonance points (N = 40, V = 1/3). The frequencies are x = 2.255 (filled circle), x = 2.415 (open circle), x = 2.56 (filled square), and x = 2.689 (open square).

K. Yoshimura, Y. Doi / Wave Motion 45 (2007) 83–99

91

distribution is clearly seen: i.e., a particular set of the normal modes is strongly excited. The modal energy is enhanced for m = 8, 7, 6, and 5 for x = 2.255, 2.415, 2.56, and 2.689, respectively. Theoretical analysis of the resonance phenomenon will be given in the next section. 4.3. Stability Linear stability of a moving DB can be examined by numerically performing the extended Floquet analysis, as it has been done in [17]. Results on the linear stability are shown in Figs. 4(a)–(c). The shaded bands indicate the intervals in x, in which the moving DB is linearly unstable. It should be noted that the instability appears for x close to the resonance points. This shows that the resonance of moving DB has a significant influence on its linear stability. In small x region, it is found that the unstable intervals in x are very narrow. That is, a weakly localized moving DB is stable for most values of x in this region. The unstable interval in x near each resonance point becomes wider as x increases up to the threshold, at which the sharp decrease occurs in etail/emax, in Figs. 4(b) and (c). In the region just above this threshold, the unstable intervals in x are narrow and thus the moving DB is stable for a large fraction of x values. It is not still clear whether this instability is destructive or not. We performed a preliminary numerical experiment on this issue. A small perturbation was added to an unstable moving DB at the initial and its time evolution was examined. We observed that the perturbed moving DB emits a small ripple and reduces its velocity, but it continues to propagate with a smaller velocity. Therefore, it seems that the instability is not destructive. A more complete study is necessary to clarify this issue. 5. Theoretical analysis of resonance We develop an approximate theory to derive a condition for the resonance. Then, we apply this resonance condition to the numerical results shown in the previous section. 5.1. Stationary DB solution We calculate an approximate solution for the stationary DB by using the method proposed in Ref. [24] for later use. It is useful to introduce an envelope function f(n, t) defined by f(n, t) = (1)n qn(t) to describe shortwavelength excitations in the lattices. The envelope function is assumed to be slowly varying with respect to n. Using the continuous approximation, we regard n as a continuous variable x and denote the envelope function by f(x, t). If we substitute qn(t) = (1)nf(x, t)|x=n into Eq. (2) and expand to an appropriate order, we can obtain the partial differential equation

ftt þ 4f þ fxx þ b 16f 3 þ 6f ðf 2 Þxx ¼ 0: ð14Þ We look for the stationary DB solution of Eq. (14) in the form f(x, t) = g(x)cos(xt). If we substitute this form into Eq. (14) and apply the rotating wave approximation (RWA), that is, cos3 (xt) is approximated by (3/ 4)cos(xt), then we have the first integral



2 4  x2 g2 þ 1 þ 9bg2 ðgx Þ þ 6bg4 ¼ C; ð15Þ where C is a constant. We consider Eq. (15) for x 2 (1, 1) although the periodic boundary condition is assumed for the original Hamiltonian (1). We look for the single stationary DB solution to Eq. (15), which satisfies the boundary conditions g ! 0 and gx ! 0 as x ! ±1. From these boundary conditions, we immediately find that C = 0. Let A be the maximal value of g. Since the derivative gx vanishes when g = A, we have a relation between x and A from Eq. (15) as x2 ¼ 4 þ 6bA2 :

ð16Þ

The frequency x is larger than the maximal natural frequency, i.e., x > 2. If we use Eq. (16) and C = 0 in Eq. (15) and perform the integration, we have

92

K. Yoshimura, Y. Doi / Wave Motion 45 (2007) 83–99

1 jx  x0 j ¼ pffiffiffiffiffiffi 6b

Z g

A

1 g

sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 1 þ 9bg2 dg; A2  g 2

ð17Þ

where x0 is the center of the stationary DB such that g(x0) = A. In the case of small amplitude (9bA2  1), we can compute the explicit form of g(x) by neglecting the term 9bg2 in the integrand and performing the integration in Eq. (17). Using the explicit form of g(x), we can obtain an approximate solution for the stationary DB as follows: hpffiffiffiffiffiffi i qn ðtÞ ¼ ð1Þn Asech 6bAn cosðxtÞ; ð18Þ pffiffiffiffiffiffi where we assumed the DB center is located at n = 0. In Eq. (18), sech½ 6bAn rapidly decreases as |n| increases. pffiffiffiffiffi ffi The width of this stationary DB is roughly given by 1= 6bA. The above stationary DB solution is converted into that in the normal mode coordinates via Eq. (4). If we substitute Eq. (18) into Eq. (4), consider n in the range N/2 + 1 6 n 6 N/2, and take the inverse transformation, we have     N =2 hpffiffiffiffiffiffi i  X A 2p 2p Qm ð0Þ ¼ pffiffiffiffi mn  sin mn ; ð19Þ sech 6bAn  cos N N N n¼ðN =21Þ where we set t = 0. The time dependence of Qm(t) is simply given by Qm(t) = Qm(0)cos(xt) from Eq. We pffiffiffiffiffi(18). ffi assume that the lattice system is much larger than the width of the stationary DB, i.e., N =2  1= 6b A. Since pffiffiffiffiffiffi sech½ 6bAn is very small for large |n|, the finite sum in Eq. (19) can be replaced by the infinite sum from n = 1 to n = 1. We approximate this infinite sum by the integral to have     Z 1 hpffiffiffiffiffiffi i  A 2p 2p p ffiffiffiffi Qm ð0Þ ’ mx  sin mx dx sech 6bAx  cos N N N 1  2  p pm ¼ pffiffiffiffiffiffiffiffiffi sech pffiffiffiffiffiffi : ð20Þ 6bN 6bAN From Eq. (20), we can find the following estimation for the mode amplitude: Qm ð0Þ ’ 2Kqjmj=N ;

ð21Þ

where K and q are defined by     p2 p2 ð22Þ q ¼ exp  pffiffiffiffiffiffi ¼ exp  pffiffiffiffiffiffiffiffiffiffiffiffiffiffi : 6bA x2  4 pffiffiffiffiffiffiffiffiffiffiffiffiffiffi pffiffiffiffiffiffi In the last equality, we used the relation 6bA ¼ x2  4, which is obtained from Eq. (16). Eq. (22) shows that 0 < q1 for x not much larger than 2 and q goes to zero in the limit x ! 2. p K ¼ pffiffiffiffiffiffiffiffiffi ; 6bN

5.2. Complex normal mode and translation map We introduce the complex normal modes Um 2 C, m = (N/2  1), . . ., N/2 1 defined by 1 i U m ¼ ðQm þ Qm Þ þ ðQm  Qm Þ: 2 2

ð23Þ

The variables satisfy the relation U m ¼ U m , where * stands for complex conjugate, since both Qm and Qm are real variables. The potential function U of Hamiltonian (1) is written in terms of Um as follows: U¼

Nh 1 X b xi xi U i U i þ 2 i¼N h 4N

Nh X

xi xj xk xl U i U j U k U l Dði þ j þ k þ lÞ:

i;j;k;l¼N h

The equations of motion for Um is given by

ð24Þ

K. Yoshimura, Y. Doi / Wave Motion 45 (2007) 83–99

d2 U m @U ¼ : 2 @U m dt

93

ð25Þ

Using Eq. (24), we have the equations of motion as follows: d2 U m b þ x2m U m þ N dt2

Nh X

xm xi xj xk U i U j U k Dðm  ði þ j þ kÞÞ ¼ 0;

ð26Þ

i;j;k¼N h

where Nh = N/2  1 and we used xm = xm. The use of the complex variables Um simplifies the selection rule for the interaction among the normal modes: the selection rule is simply described by the function D while it is described by D in Eq. (9), which has a more complicated form. Therefore, the complex normal mode coordinates defined by Eq. (23) are suitable for a theoretical analysis. Consider an arbitrary point (q, p) in the phase space of system (1), provided that QN/2 = PN/2 = 0 hold. Suppose that this point is expressed by fU m ; U_ m g in the complex normal mode coordinates, where m = (N/ 2  1), . . . , N/2  1. We define the translation map T : C2ðN1Þ ! C2ðN1Þ as follows: T : fU m ; U_ m g7!fU m exp½i2hm ; U_ m exp½i2hm g;

ð27Þ

where hm = pm/N. This map T corresponds to shifting the displacement and momentum pattern (q, p) by one lattice spacing and reversing the phase of each particle in the physical space. 5.3. Resonance condition Let us consider a moving DB solution with a rational velocity V = r/s, where r and s are integers. We assume the existence of moving DB solutions although it has not yet been rigorously proved. We use new variables wm instead of Um to describe this solution, which are defined by the transformation h mr i ð28Þ U m ðtÞ ¼ wm ðtÞ exp i xt ; Ns where wm(t) 2 C is a complex function of t. The relation wm ¼ wm holds since U m ¼ U m . Note that for t = sT, where T = 2p/x, the factor exp[i(mr/Ns) xt] in this equation coincides with the factor exp[i2hm r] that appears in the map Tr . Let (q(t),p(t)) be the moving DB solution with V = r/s represented in the physical space. Since V is a rational number, (q(t + sT),p(t + sT)) has to coincide with the displacement and momentum pattern obtained by reversing the phase of each particle r times after translating (q(t), p(t)) by r lattice spacings, where T = 2p/x: i.e., the condition fU m ðt þ sT Þ; U_ m ðt þ sT Þg ¼ Tr fU m ðtÞ; U_ m ðtÞg has to be satisfied. It follows from this condition that wm(t + sT) = wm(t): wm(t) is a periodic function with the period sT. It should be noted that transformation (28) separates the internal oscillation from the translational motion of the moving DB: wm(t) and exp[i(mrx/Ns)t] represent the internal oscillation and the translational motion, respectively. If we substitute Eq. (28) into Eq. (26), then we have  Nh mrx2  d2 wm 2mrx dwm b X 2 þ x  i  ¼  xm xi xj xk wi wj wk ei½fmðiþjþkÞgr=Nsxt Dðm w m m Ns dt Ns N i;j;k¼N h dt2  ði þ j þ kÞÞ:

ð29Þ wmð0Þ ; m

First, we consider Eq. (29) in the case of V = 0 for later use. Let ¼ ðN =2  1Þ; . . . ; N =2  1 be the set of solutions to Eq. (29) for V = 0 that represents the stationary DB. The equations of motion for wð0Þ m can be obtained by setting r = 0 in Eq. (29) and they read d2 wð0Þ b m þ x2m wð0Þ m ¼  N dt2

Nh X

ð0Þ

xm xi xj xk wið0Þ wð0Þ j wk Dðm  ði þ j þ kÞÞ:

ð30Þ

i;j;k¼N h

Let l be an integer such that l 2 {0, 1, . . ., 2N  1}. This set of equations has a solution of the form wð0Þ m ðtÞ ¼ /m ðtÞ exp½ilhm  for each l, where /m(t) 2 R is a real function of t. The solutions for even l correspond to odd modes while those for odd l correspond to even modes. The integer l specifies the position of the DB center in the physical space. The solution for l = 0 corresponds to the stationary DB with the odd sym-

94

K. Yoshimura, Y. Doi / Wave Motion 45 (2007) 83–99

metry centered at n = 0. We choose the solution for l = 0 as wð0Þ m . In the case of l = 0, Um(t) is a real function of t because U m ðtÞ ¼ wmð0Þ ðtÞ holds for V = 0 from Eq. (28) and wð0Þ m ðtÞ is real. Then, it follows from Eq. (23) that wmð0Þ ðtÞ ¼ Qm ðtÞ since Um(t) is real. Since wmð0Þ represents the stationary DB solution, wð0Þ m ðtÞ is periodic and can be expanded by Fourier series. This Fourier series includes the frequency components ±x, ±3x, ±5x, . . .. However, it is well known that only the fundamental frequency component is dominant in stationary DBs. We employ the approximate solution (18) to obtain wmð0Þ ðtÞ. Using Eq. (20) and taking account of the relation wð0Þ m ðtÞ ¼ Qm ðtÞð¼ Qm ð0Þ cosðxtÞÞ, we have an approximation for wð0Þ as follows: m X ð0Þ wm ðtÞ ¼ Am exp½ilxt; ð31Þ l¼ 1

where Am is defined by

 2  p pm Am ¼ pffiffiffiffiffiffiffiffiffi sech pffiffiffiffiffiffi ’ Kqjmj=N : 2 6bN 6bAN

ð32Þ

Next, we consider the case of moving DB and derive a condition for the resonance. The internal oscillation of a moving DB is similar to that of the stationary DB with the same x if the translational motion is separated. Therefore, it is useful to decompose wm into two parts as follows: wm ðtÞ ¼ wð0Þ m ðtÞ þ um ðtÞ;

ð33Þ

where wmð0Þ represents the stationary DB solution and um is the deviation of wm from wmð0Þ . It is expected that um is small unless x is too close to a resonance frequency. The deviation um(t) is also a periodic function with the period sT since both wm(t) and wð0Þ m ðtÞ have the period sT. Therefore, um(t) can be expanded by the Fourier series um ðtÞ ¼

1 X

am;n exp ½inðx=sÞt;

ð34Þ

n¼1

where am,n is the Fourier coefficient. If we substitute Eq. (33) into Eq. (29) and use Eq. (30), then we can obtain the following equation up to the first order of um:  mrx2  d2 um 2mrx dum 2mrx dwmð0Þ mrx2 ð0Þ 2 þ x þ  i  wm um ¼ i m 2 Ns dt Ns Ns Ns dt dt Nh



 b X ð0Þ  þ xm xi xj xk wið0Þ wð0Þ 1  eiðr=sÞxt DN þ 1  eiðr=sÞxt DN j wk N i;j;k¼N h 

3b N

Nh X

ð0Þ i½fmðiþjþkÞgr=Nsxt xm xi xj xk wð0Þ Dðm  ði þ j þ kÞÞ þ Oðu2 Þ: i wj uk e

ð35Þ

i;j;k¼N h

where D+N and DN are the functions of the four indices defined by  1 if m  ði þ j þ kÞ ¼ N ; D N ¼ 0 otherwise:

ð36Þ

ð0Þ We evaluate the order of the coefficient of each uk on the right hand side of Eq. (35). The order of wð0Þ i wj is evaluated as 1/N from Eqs. (31) and (32). If k is fixed, the number of terms included in the sum with respect to i and j is of the order of N. Thus, the coefficient of uk is of the order of 1/N because of the factor 3b/N. This indicates that for each k the term including uk on the right hand side is much smaller than the terms depending on um on the left hand side, whose coefficients are of the order of unity, for large N. In addition, it is expected that the whole contribution from the sum taken over all k remains small since for each frequency the Fourier components contained in the terms for different k oscillate with different phases and they may cancel with each other. Therefore, the terms of O(uk) on the right hand side may be neglected in the lowest order approxima-

K. Yoshimura, Y. Doi / Wave Motion 45 (2007) 83–99

95

tion. If we neglect these terms and substitute Eqs. (31) and (34) into Eq. (35), we arrive at the following equation 1   mr n mr2  X X mr  2l x2 Am eilxt x2m   x2 am;n einxt=s ¼ s Ns Ns Ns n¼1 l¼ 1 ( ) Nh X b X þ xm xi xj xk Ai Aj Ak ðDN þ DN Þ eiðl1 þl2 þl3 Þxt N l ;l ;l ¼ 1 i;j;k¼N h 1 2 3 ( ) Nh X b X  xm xi xj xk Ai Aj Ak DN eiðl1 þl2 þl3 þr=sÞxt N l ;l ;l ¼ 1 i;j;k¼N h 1 2 3 ( ) Nh X b X  xm xi xj xk Ai Aj Ak DN eiðl1 þl2 þl3 r=sÞxt : ð37Þ N l ;l ;l ¼ 1 i;j;k¼N h 1 2 3

The Fourier coefficients am,n are determined by equating the same Fourier components between the left and right hand sides. Eq. (37) shows that |am,n| becomes large and the resonance occurs when the coefficient of am, n becomes close to zero. Eq. (37) is a reasonable approximation unless the coefficient of am, n is too close to zero. Thus, the resonance condition is given by n mr   ð38Þ   x ¼ xm : s Ns The resonance occurs if n, m, and x approximately satisfy Eq. (38) for given N and V = r/s. We note that am, n remains finite due to the terms neglected in deriving Eq. (37) even if condition (38) holds exactly, although approximate equation (37) indicates divergence of am, n. Because of the RWA, only the Fourier components with the frequencies ±x, ±3x, ±(3 ± r/s)x, and ±(1 ± r/s)x appear in the right hand side of Eq. (37). Therefore, as for n, it is enough to consider the Fourier components for n = ±s, ±3s, ±(3s ± r), and ±(s ± r) in the left hand side. Hereafter, we assume that V = r/s is in the range 0 6 V 6 1/3, which correspond to our numerical experiments. It can be checked that condition (38) is not satisfied for n = ±3s and ±(3s ± r) by taking into account the facts m/N < 1/2 and xm < 2 < x. Solving Eq. (38) with respect to x, we have x = xm/|n/s  mr/Ns|. In the case of n = ±s, we have xm/|n/s  mr/Ns| < 2.012 for V 2 [0, 1/3]. Therefore, condition (38) is not satisfied for x used in our numerical experiments in the case of n = ±s. Eq. (38) is rewritten as x = xm/|1 + (1 m/N)r/s| when n = ±(s + r). Using m/N < 1/2 and xm < 2, we have xm/|1 + (1 m/N)r/s| < 2/(1 + r/2s). This implies that condition (38) does not hold also for n = ±(s + r) because 2/(1 + r/2s) < 12/7 ’ 1.714 for V 2 [0, 1/3] and x > 2. In summary, it has been shown that there does not exist a set (n, m, x) satisfying Eq. (38) for n = ±s, ±3s, ±(3s ± r), and ±(s + r). The above arguments show that condition (38) can hold only in the two cases n = ±(sr). We consider the first case of n = sr, for which Eq. (38) reads n  m  ro 1 1þ ð39Þ x ¼ xm : N s ~ m stands for ~ m Þ, where x A number of discrete pairs (m, x) satisfy this condition. We denote these pairs by ðm; x the resonance frequency for m. Some of the pairs have positive m and the others have negative m. There is the ~ m satisfies maximum of m in the pairs with positive m and we denote it by mp. In the case of positive m, x ~ mp < x ~ mp 1 < x ~ mp 2 <   . Similarly, we have x ~ mn < x ~ mn þ1 < x ~ mn þ2 <    in the case of negative m, where x we denote the minimum of m by mn. The values of mp and mn depend on V = r/s and N. In Eq. (37), the Fourier component for n = s  r has the frequency (1  r/s)x. If we equate the Fourier components with this frequency in the both sides, we have   Nh n  m  r o2 2 3b X x2m  1  1 þ x am;sr ¼  xm xi xj xk Ai Aj Ak DN : ð40Þ N s N i;j;k¼N h We estimate the order of the right hand side (rhs). Using Eq. (32), we can obtain the order estimation

96

K. Yoshimura, Y. Doi / Wave Motion 45 (2007) 83–99

xm xi xj xk Ai Aj Ak ’ K 3 qðjijþjjjþjkjÞ=N

ð41Þ

for each term in the sum. The sum in Eq. (40) includes only the terms satisfying i + j + k = N + m due to the function DN. It follows that |i| + |j| + |k| P |i + j + k| = |N + m|, where the equality holds if and only if i + j + k = N + m and i, j, k P 0. Therefore, the terms for i, j, k P 0 are of the largest order since q < 1. These terms are evaluated as xm xi xj xk Ai Aj Ak ’ K 3 qjNþmj=N ;

ð42Þ

2 The number of terms satisfying i + j + k = N + m and i, j, k P 0 is of pthe ffiffiffiffi order of N . If we evaluate the sum in Eq. (40) by the largest order terms and take account of K ¼ Oð1= N Þ, we have pffiffiffiffi ð43Þ ðrhsÞ ¼ OðqjNþmj=N = N Þ:

~ mþ Þ and ðm ; x ~ m Þ satisfying Eq. (39), where m+ > 0 and m < 0. Suppose that Let us consider two pairs ðmþ ; x ~ m . The values of q are almost equal to each other between these two cases from this assumption. Eq. ~ mþ ’ x x (43) indicates that the order of the right hand side is larger for m than m+. Consider Eq. (40) in the two cases: ~ mþ j ¼ d þ  1 and the other case is that m = m and jx  x ~ m j ¼ d   1. one case is that m = m+ and jx  x The order of am, sr is determined mainly by the right hand side of Eq. (40) when d+ ’ d. Therefore, |am, s  r| ~ m Þ leads to a is larger in the latter case than the former case when d+ ’ d. This implies that the pair ðm ; x ~ mþ Þ when n = s  r. It can be concluded that the pairs ðm; x ~ m Þ with negative m stronger resonance than ðmþ ; x correspond to the stronger resonances in the case of n = s  r. In the second case of n = (s  r), Eq. (38) reads n  m ro ð44Þ 1 1 x ¼ xm : N s ~ m Þ, which satisfy Eq. (44), as in the first case. Both of the pairs with There are the same number of pairs ðm; x positive m and negative m exist. The Fourier component for n = (s  r) has the frequency (1  r/s)x in Eq. (37). Equating the Fourier components with this frequency in the both sides, we have   Nh n  m ro2 2 3b X x2m  1  1  x am;ðsrÞ ¼  xm xi xj xk Ai Aj Ak DN : ð45Þ N s N i;j;k¼N h A similar order estimation of the right hand side results in pffiffiffiffi ðrhsÞ ¼ OðqjNmj=N = N Þ:

ð46Þ

~ m Þ with positive m corBased on a consideration similar to the above one, we can conclude that the pairs ðm; x respond to the stronger resonances in the case of n = (s  r). If we combine Eq. (39) with m < 0 and Eq. (44) with m > 0, we can obtain     jmj ð47Þ 1 1 V x ¼ xm ; N where we replaced r/s with V. Eq. (47) gives the condition for the fundamental resonances. The condition (47) means that if there exists a pair of m and xm that approximately satisfy Eq. (47) for a given x, then the corresponding modal energy Em becomes large. Eq. (47) indicates that the resonance occurs when the linear combination of the internal frequency x and the frequency Vx, which emerges due to the motion of DB, coincides with the natural frequency xm of a certain normal mode. We apply the condition (47) to the four fundamental resonances for V = 1/3. Each of these fundamental resonances is indicated in Fig. 4(a) by showing its resonant normal mode number m, which is observed in Fig. 6. Solving Eq. (47) with respect to x, we have x = xm/{1  (1  |m|/N)V}. Given m, the right hand side is determined. The right hand side gives the resonance frequency of DB such that the ±mth normal modes are strongly excited. We evaluate this theoretical value for the resonance frequency of DB by using the parameter values for the relevant four fundamental resonances. The resonance frequencies predicted from Eq. (47) and those found numerically are compared in Table 1. The theoretical values are roughly in agreement with the numerical one. This agreement shows that the resonance condition (47) can describe the fundamental reso-

K. Yoshimura, Y. Doi / Wave Motion 45 (2007) 83–99

97

Table 1 Comparison of theoretical and numerical values of breather frequency x for resonance Numerical

Eq. (47)

Eq. (50)

2.255 2.415 2.560 2.689

2.206 2.352 2.487 2.609

2.264 2.432 2.586 2.727

(m = 8) (m = 7) (m = 6) (m = 5)

nances. The condition (47) also applies to the other two cases of V = 1/5 and 1/7: the resonance peaks described by Eq. (47) are indicated by showing the resonant normal mode number m in Figs. 4(b) and (c). In the derivation of Eq. (47), the terms of O(uk) on the right hand side of Eq. (35) were neglected. This approximation may become better as N increases because the coefficients of these neglected terms are of the order of 1/N. Therefore, it is expected that the condition (47) describes the resonance frequencies more precisely for large N. There is a discrepancy between the theoretical and numerical resonance frequencies in Table 1. One of the reasons for this discrepancy may be a change in the natural frequency of the resonant normal mode due to the nonlinearity of lattice potential. To show this, we derive a resonance condition with a correction term and compare it with the numerical results. Consider the sum of the first order terms in uk on the right hand side in Eq. (35), which we denote by S1. If we substitute Eqs. (31) and (34) into this sum, we can obtain ! Nh   6b 2 X 2 2 ð48Þ S 1 ¼  xm x i Ai am;sr eið1r=sÞxt þ am;ðsrÞ eið1r=sÞxt þ ðother termsÞ: N i¼N h P h Noting the fact that J ¼ Ni¼N x2i A2i is the half of the lattice harmonic energy, we can compute this term by h using Eq. (18) as follows: Z 1 N =2 hpffiffiffiffiffiffi i X 1 1 pffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 2 J¼ 6bAx dx ¼ ðqnþ1  qn Þ ’ A sech2 ð49Þ x2  4; 4 n¼ðN =21Þ 3b 1 pffiffiffiffiffiffiffiffiffiffiffiffiffiffi pffiffiffiffiffiffi where we used 6bA ¼ x2  4. If we include the first two terms of S1 into the left hand side of Eq. (37), we obtain the improved resonance condition  1=2     jmj 2 pffiffiffiffiffiffiffiffiffiffiffiffiffiffi 1 1 : ð50Þ V x ¼ xm 1 þ x2  4 N N Note that the correction term vanishes in the limit N ! 1. If we solve Eq. (50) with respect to x for a given m, we can determine the resonance frequency of DB such that the ±mth normal modes are strongly excited. The resonance frequencies x predicted from Eq. (50) are shown in Table 1. These theoretical values are in good agreement with the numerical one, validating the above explanation to the discrepancy. We give an explanation to the sharp decrease of an average part of etail/emax, which is observed at x ’ 2.5 in Fig. 4(b) and x ’ 2.33 in Fig. 4(c), based on Eq. (47). In Eq. (47), xm satisfies 0 6 xm 6 2 because of the dispersion relation (5) and the inequality 1  V 6 1  (1  |m|/N)V also holds. Therefore, the condition (47) can be satisfied only when x 6 xc = 2/(1  V). For any x smaller than the upper bound xc, there exist a few normal modes which at least roughly satisfy Eq. (47) and have the enhanced Fourier coefficients am, n. The ratio etail/emax is rather large due to the contribution from these normal modes. In contrast, if x becomes slightly larger than xc, no normal mode can satisfy Eq. (47). In addition, the condition {1  (1 + |m|/N)V}x = xm for the weak resonance also cannot be satisfied in this region in x. Thus, the coefficient of am, n deviates from a vicinity of zero in Eq. (37) for all the nonvanishing Fourier components under the RWA and thus the modulus of any am, n is no longer large. This results in a considerable decrease of the tail amplitude. The upper bound xc is obtained as xc = 2.5 for V = 1/5 and xc = 2.33 for V = 1/7. These theoretical values are in good agreement with the numerical one. It may be stated that a sharp decrease of an average part of etail/emax occurs when x = xc.

98

K. Yoshimura, Y. Doi / Wave Motion 45 (2007) 83–99

5.4. Remarks We make some remarks on a relation between a symmetry of the potential function and the resonance of moving DB. Let us define the map Tk : CN1 ! CN1 by Tk : U m 7!U m exp½imk;

ð51Þ

where k is a real parameter and m = (N/2  1), . . ., N/2  1. This map is a continuous version of the map T restricted to the configuration space, where T is given by Eq. (27). In the real normal mode coordinates, Tk is expressed as follows: Tk : ðQm ; Qm Þ7!ðQm ; Qm Þ  RðmkÞ;

ð52Þ

where R(mk) is the 2 · 2 matrix   cosðmkÞ sinðmkÞ RðmkÞ ¼ :  sinðmkÞ cosðmkÞ

ð53Þ

The potential function U given by Eq. (24) can be divided into the two parts Us(U) and Ua(U), which are defined by Us ðUÞ ¼

Nh 1 X b xi xi U i U i þ 2 i¼N h 4N

Nh X

xi xj xk xl U i U j U k U l ;

ð54Þ

i;j;k;l¼N h iþjþkþl¼0

and Ua ðUÞ ¼ 

b 4N

Nh X

xi xj xk xl U i U j U k U l ;

ð55Þ

i;j;k;l¼N h iþjþkþl¼ N

where U = (UN/2+1, . . ., UN/21). The former part Us is the symmetric part with respect to the transformation Tk , i.e., Us ðTk UÞ ¼ Us ðUÞ for any k, and the latter part Ua is the asymmetric part. Eqs. (40) and (45) indicate that the resonance is caused by the terms including DN or DN. These terms come from the asymmetric part Ua. This fact indicates that the potential asymmetry with respect to Tk is the origin of the resonance of moving DB. It is expected that nonlinear lattice models in general have such asymmetry because a lattice has an additional integral if its potential does not include the asymmetric part. Therefore, the resonance of moving DB may be observed in various nonlinear lattice models such as the nonlinear Klein– Gordon lattice, although we have revealed the phenomenon only in the FPU-b lattice. The present theory, which is a perturbation approach after removing the translational motion of moving DB by using the complex normal mode coordinates, may be useful to derive the resonance conditions in various lattice models. It is mathematically possible to construct a lattice model, the potential of which is symmetric with respect to Tk , although the lattice potential may be physically unnatural. Indeed, if we remove the Ua from the potential of the FPU-b lattice, we can obtain a symmetric lattice model. It is conjectured that the resonance and the resultant instability do not occur in symmetric lattice models: there may exist moving DBs, which have no apparent tails and propagate stably. 6. Conclusions We studied the nature of radiationless moving DBs with rational velocities in the FPU-b lattice. The Newton method allows to find precise numerical solutions for these moving DBs. A moving DB is not strictly localized but has a spatially extended tail. The tail amplitude sensitively depends on the frequency and velocity of the moving DB. We found that a resonance phenomenon in the tail amplitude occurs: the tail amplitude is sharply enhanced for some particular values of the internal frequency x of moving DB and these values depend on the velocity. It was shown that the resonance phenomenon is caused by strong excitation of a particular set of the normal modes. We investigated the linear stability of moving DB by numerically performing the extended Floquet analysis. It was found that the instability appears for x close to the resonance points.

K. Yoshimura, Y. Doi / Wave Motion 45 (2007) 83–99

99

The resonance has a significant influence on the stability of the moving DB. Moreover, we developed a theory to describe the resonance phenomenon. An analytical expression for the resonance condition has been derived. This analytical resonance condition is in agreement with the numerical results. The present theory well describes the resonance phenomenon of moving DB. Acknowledgment One of the authors (K.Y.) thank the members of NTT Communication Science Laboratories for their continual encouragement. References [1] S. Takeno, K. Kisoda, A.J. Sievers, Intrinsic localized vibrational modes in anharmonic crystals: stationary modes, Prog. Theor. Phys. Suppl. 94 (1988) 242–269; A.J. Sievers, S. Takeno, Intrinsic localized modes in anharmonic crystals, Phys. Rev. Lett. 61 (1988) 970–973. [2] S. Flach, C. Willis, Discrete breathers, Phys. Rep. 295 (1998) 181–264. [3] E. Trias, J.J. Mazo, T.P. Orlando, Discrete breathers in nonlinear lattices: experimental detection in a Josephson array, Phys. Rev. Lett. 84 (2000) 741–744; P. Binder, D. Abraimov, A.V. Ustinov, S. Flach, Y. Zolotaryuk, Observation of breathers in Josephson ladders, Phys. Rev. Lett. 84 (2000) 745–748. [4] H.S. Eisenberg, Y. Silberberg, R. Morandotti, A.R. Boyd, J.S. Aitchison, Discrete spatial optical solitons in waveguide arrays, Phys. Rev. Lett. 81 (1998) 3383–3386. [5] M. Sato, B.E. Hubbard, A.J. Sievers, B. Ilic, D.A. Czaplewski, H.G. Craighead, Observation of locked intrinsic localized vibrational modes in a micromechanical oscillator array, Phys. Rev. Lett. 90 (2003) 044102; M. Sato, B.E. Hubbard, A.J. Sievers, B. Ilic, H.G. Graighead, Optical manipulation of intrinsic localized vibrational energy in cantilever arrays, Europhys. Lett. 66 (2004) 318–323. [6] S. Takeno, K. Hori, A propagating self-localized mode in a one-dimensional lattice with quartic anharmonicity, J. Phys. Soc. Jpn. 59 (1990) 3037–3040; S. Takeno, K. Hori, Self-localized modes in a pure one-dimensional lattice with cubic and quartic lattice anharmonicity, J. Phys. Soc. Jpn. 60 (1991) 947–959. [7] D. Chen, S. Aubry, G.P. Tsironis, Breather mobility in discrete U4 lattices, Phys. Rev. Lett. 77 (1996) 4776–4779. [8] E. Fermi, J. Pasta, S. Ulam, in: E. Segre´ (Ed.), Collected Papers of E. Fermi, University of Chicago, Chicago, 1965. [9] V.M. Burlakov, S.A. Kiselev, V.I. Rupasov, Localized vibrations of homogeneous anharmonic chains, Phys. Lett. A 147 (1990) 130– 134. [10] T. Cretegny, T. Dauxois, S. Ruffo, A. Torcini, Localization and equipartition of energy in the b-FPU chain: chaotic breathers, Physica D 121 (1998) 109–126. [11] K. Ullmann, A.J. Lichtenberg, G. Corso, Energy equipartition starting from high-frequency modes in the Fermi–Pasta–Ulam b oscillator chain, Phys. Rev. E 61 (2000) 2471–2477. [12] Y.A. Kosevich, S. Lepri, Modulational instability and energy localization in anharmonic lattices at finite energy density, Phys. Rev. B 61 (2000) 299–307. [13] N. Budinsky, T. Bountis, Stability of nonlinear modes and chaotic properties of 1D Fermi–Pasta–Ulam lattices, Physica D 8 (1983) 445–452. [14] P. Poggi, S. Ruffo, Exact solutions in the FPU oscillator chain, Physica D 103 (1997) 251–272. [15] K. Yoshimura, Parametric resonance energy exchange and induction phenomenon in a one-dimensional nonlinear oscillator chain, Phys. Rev. E 62 (2000) 6447–6461. [16] K. Yoshimura, Modulational instability of zone boundary mode in nonlinear lattices: rigorous results, Phys. Rev. E 70 (2004) 016611. [17] S. Aubry, T. Cretegny, Mobility and reactivity of discrete breathers, Physica D 119 (1998) 34–46. [18] J. Go´mez-Garden˜es, F. Falo, L.M. Florı´a, Mobile localization in nonlinear Schro¨dinger lattices, Phys. Lett. A 332 (2004) 213–219. [19] K. Yoshimura, K. Umeno, Nonintegrability of nonhomogeneous nonlinear lattices, J. Math. Phys. 45 (2004) 4628–4639. [20] J.B. Page, Asymptotic solutions for localized vibrational modes in strongly anharmonic periodic systems, Phys. Rev. B 41 (1990) 7835–7838. [21] K.W. Sandusky, J.B. Page, K.E. Schmidt, Stability and motion of intrinsic localized modes in nonlinear periodic lattices, Phys. Rev. B 46 (1992) 6161–6168. [22] J.L. Marin, S. Aubry, Breathers in nonlinear lattices: numerical calculation from the anticontinuous limit, Nonlinearity 9 (1996) 1501–1528. [23] K. Yoshimura, Analytical study of discrete breathers in Fermi–Pasta–Ulam lattice, unpublished. [24] Y.A. Kosevich, Nonlinear envelope-function equation and strongly localized vibrational modes in anharmonic lattices, Phys. Rev. B 47 (1993) 3138–3152.