From single- to multiple-soliton solutions of the perturbed KdV equation

Physica D 237 (2008) 2987–3007

Contents lists available at ScienceDirect

Physica D journal homepage: www.elsevier.com/locate/physd

From single- to multiple-soliton solutions of the perturbed KdV equation Yair Zarmi ∗ Jacob Blaustein Institutes for Desert Research & Physics Department, Ben-Gurion University of the Negev, Midreshet Ben-Gurion, 84990, Israel

article

info

Article history: Received 5 June 2008 Accepted 10 July 2008 Available online 20 July 2008 Communicated by A.C. Newell PACS: 02.30.IK 02.30.Mv 05.45.-a Keywords: Perturbed KdV equation Solitons Anti-solitons Elastic and inelastic interactions

a b s t r a c t The solution of the perturbed KdV equation (PKDVE), when the zero-order approximation is a multiplesoliton wave, is constructed as a sum of two components: elastic and inelastic. The elastic component preserves the elastic nature of soliton collisions. Its perturbation series is identical in structure to the series-solution of the PKDVE when the zero-order approximation is a single soliton. The inelastic component exists only in the multiple-soliton case, and emerges from the first order and onwards. Depending on initial data or boundary conditions, it may contain, in every order, a plethora of inelastic processes. Examples are given of sign-exchange soliton–anti-soliton scattering, soliton–anti-soliton creation or annihilation, soliton decay or merging, and inelastic soliton deflection. The analysis has been carried out through third order in the expansion parameter, exploiting the freedom in the expansion to its fullest extent. Both elastic and inelastic components do not modify soliton parameters beyond their values in the zero-order approximation. When the PKDVE is not asymptotically integrable, the new expansion scheme transforms it into a system of two equations: The Normal Form for ordinary KdV solitons, and an auxiliary equation describing the contribution of obstacles to asymptotic integrability to the inelastic component. Through the orders studied, the solution of the latter is a conserved quantity, which contains the dispersive wave that has been observed in previous works. © 2008 Elsevier B.V. All rights reserved.

1. Introduction The generic form of the perturbed KdV equation (PKDVE) is [1–6]:

wt = 6ww1 + w3
+ ε 30α1 w2 w1 + 10α2 ww3 + 20α3 w1 w2 + α4 w5   140β1 w 3 w1 + 70β2 w 2 w3 + 280β3 ww1 w2 + 14β4 ww5 + 70β5 wx3 + ε2 +42β6 w1 w4 + 70β7 w2 w3 + β8 w7   4 630γ1 w w1 + 1260γ2 w (w1 )3 + 2520γ3 w 2 w1 w2 + 1302γ4 w1 (w2 )2 3  +ε +420γ5 w 3 w3 + 966γ6 (w1 )2 w3 + 1260γ7 ww2 w3 + 756γ8 ww1 w4 +252γ9 w3 w4 + 126γ10 w 2 w5 + 168γ11 w2 w5 + 72γ12 w1 w6 + 18γ13 ww7 + γ14 w9

|ε|  1, wp ≡ ∂xp w . + O ε4

(1.1)

The unspecified coefficients depend on the dynamical system for which Eq. (1.1) provides an approximation. For example, small-amplitude solutions of surface waves on a shallow-water-layer over a horizontal plane [7], of the ion acoustic wave equations in Plasma Physics [8,9], p and the Fermi–Pasta–Ulam problem [10] are approximately described by Eq. (1.1). (Throughout the paper, wp and up will represent ∂x w p and ∂x u, respectively.) Eq. (1.1) is integrable through O (ε ) [1–6]. Namely, if one expands w in powers of ε ,


w (t , x) = u (t , x) + ε u(1) (t , x) + O ε 2 ,

(1.2)

then the zero-order approximation, u, is determined by a Normal Form that is integrable, ut = 6uu1 + u3 + εµ1 30u2 u1 + 10uu3 + 20u1 u2 + u5 + O ε 2 ,




∗

Tel.: +972 8 659 6920. E-mail address: [email protected].

0167-2789/$ – see front matter © 2008 Elsevier B.V. All rights reserved. doi:10.1016/j.physd.2008.07.007




(1.3)

2988

Y. Zarmi / Physica D 237 (2008) 2987–3007

and the first-order correction, u(1) , has a closed-form expression as a differential polynomial in u:


u(1) = a1 u2 + a2 u1 q + a3 u2 q (t , x) = ∂x−1 u (t , x) .

(1.4)

The values of the constants µ1 and aj , j = 1, 2, 3, depend on the details of the perturbation scheme. Eq. (1.3) has the same soliton solutions as the KdV equation. The single-soliton solution is uSingle (t , x; k) = 2k2 / cosh [k (x + v t )]2 ,

(1.5)

and the two-soliton solution is given by [11]

( u (t , x) = 2∂ ln 1 + g1 (t , x) + g2 (t , x) + 2 x



k1 − k2

2

k1 + k2

) g 1 ( t , x) g 2 ( t , x)

(1.6)

(gi (t , x) = exp [2ki (x + vi t )]) . The perturbation in Eq. (1.3) merely updates the velocities, vi , in a known manner:

vi = 4k2i 1 + εµ1 4k2i + O ε2





.

(1.7)

As in the case of the KdV equation, solitons in the multiple-soliton solutions of Eq. (1.3) undergo elastic collisions [1–16]. Far from the soliton-collision region, the solution asymptotes into a sum of well separated pure single-soliton contributions: u ( t , x) →

X

uSingle (t , x − ηi ; ki ) .

(1.8)

i

Each soliton preserves its asymptotic characteristics before and after the collision. Soliton interactions affect only the phase shifts, ηi . Each ηi depends on the wave numbers of the other solitons. However, unless a2 6= 0, the first-order correction, given by Eq. (1.4), spoils the elastic nature of soliton collisions. The term u1 q contains a mixture of elastic and inelastic terms, due to the presence of the non-local entity q(t , x). This has been shown in the two-soliton case [17]. In the original approach [1–5], where µ1 and aj , j = 1, 2, 3, in Eq. (1.4) have the values 5 5 5 a1 = − α1 + α3 + α4 , 2 3 6

µ1 = α4 ,

a2 = −

10 3

α2 +

10 3

α4 ,

a3 = −5α1 +

5 3

α2 +

10 3

α4 ,

(1.9)

the u1 q-term exists because a2 does not vanish unless α2 6= α4 . The term u1 q can be eliminated from Eq. (1.4) if the coordinate transformation 5

τ = t + ε (α2 − α4 ) x, 3

ξ =x

(1.10)

is employed [18]. As a result of this transformation, Eq. (1.9) is replaced by

µ1 = α4 + 5 (α2 − α4 ) 5 5 20 5 a1 = − α1 + α3 + α4 + (α2 − α4 ) , 2 3 6 3

a2 = 0,

a3 = −5α1 +

5 3

α2 +

10 3

α4 +

40 3

(1.11)

(α2 − α4 ) .

In the new coordinates, (τ , ξ ), the velocity of each soliton in the zero-order approximation (solution of the Normal Form, Eq. (1.3)) is updated according to

vi = 4k2i 1 + ε (α4 + 5 (α2 − α4 )) 4k2i + O ε 2





,

(1.12)

instead of the usual [1–5]

vi = 4k2i 1 + εα4 4k2i + O ε 2





.

(1.13)

The u1 q-term has been eliminated from Eq. (1.4) in the new coordinates. However, in the multiple-soliton case, it is hidden in the zeroorder approximation, u (τ , ξ ), and reappears when the latter is expressed in terms of the original coordinates (in which the solitons are actually observed) and of observed wave numbers. Consider the single-soliton solution in the new coordinates. Substituting Eq. (1.10) in Eq. (1.5), one obtains: uSingle (τ , ξ ; k) = 2k2 / cosh [k (ξ + vτ )]2

  

5

= 2k / cosh k x 1 + ε (α2 − α4 ) 4k 2

2

3



+ 4k 1 + ε (α4 + 5 (α2 − α4 )) 4k 2

2





2 t

.

(1.14)

The observed wave number, K (the coefficient of x in Eq. (1.14)), contains a first-order correction:



K =k 1+ε

5 3

 (α2 − α4 ) 4k2 .

(1.15)

Re-writing Eq. (1.14) in terms of the observed wave number K , and the original coordinates, (t , x), to first-order in ε , the argument in Eq. (1.14) is reduced to the usual expression of Eq. (1.13): k (ξ + vτ ) = K x + 4K 2 1 + εα4 4K 2 + O ε 2




t .






(1.16)

Y. Zarmi / Physica D 237 (2008) 2987–3007

2989

In addition, when expressed in terms of the observed wave number, K , the amplitude of the solution in Eq. (1.14) contains an O (ε ) correction: 2k2 = 2K 2



1−ε

10 3

 (α2 − α4 ) K 2 .

(1.17)

Consider now the two-soliton solution, Eq. (1.6), in the (τ , ξ ) coordinates:

( u (τ , ξ ) = 2∂ξ ln 1 + g1 (τ , ξ ) + g2 (τ , ξ ) + 2



k1 − k2

2

k1 + k2

) g1 (τ , ξ ) g2 (τ , ξ ) .

(1.18)

Using Eq. (1.17) for k1 and k2 , one obtains



k1 − k2 k1 + k2

2

 =

K1 − K2 K1 + K2

2 

1−ε

80 3

(α2 − α4 ) K1 K2 + O ε 2




.

(1.19)

As a result of Eq. (1.19), Eq. (1.18) contains an inelastic first-order contribution. Re-writing Eq. (1.18) in terms of the physical coordinates, (t , x) and the observed wave numbers, K1 , K2 , and expanding the result through O (ε ), one obtains for K2 > K1 , u (τ , ξ )

→ uSingle (t , x − η1 ; K1 ) + uSingle (t , x − η2 ; K2 )

|t |→∞

−ε

10 3


(α2 − α4 ) ∂x 2K1 uSingle (t , x − η1 ; K2 ) − 2K2 uSingle (t , x − η2 ; K1 ) sgn(t ) + O ε 2 .

(1.20)

In Eq. (1.20), uSingle is given by Eq. (1.5), the argument is given by Eq. (1.16), and the phase shifts, ηi , have the standard KdV values [11–16]. The term in curly brackets in Eq. (1.20) is the inelastic part contained in the eliminated u1 q-term: ux (t , x; Ki ) of each soliton is multiplied by the asymptotic value of q(t , x; Kj ) of the OTHER soliton:




 qSingle t , x; Kj = ∂x−1 uSingle t , x; Kj = 2Kj tanh Kj x + vj t + ηj → 2Kj sgn(t ). |t |→∞

(1.21)

Thus, a pure inelastic first-order term is hidden in the solution obtained in [18], once it is expressed in terms of observable quantities. Last but not least, the validity of Eq. (1.20) is not limited to short times or distances, because it is obtained through an expansion of the amplitude of the solution, and not of the coordinates. In summary, one expansion scheme [1–5] generates first-order inelastic terms [17], whereas the other [18] avoids them. The two schemes generate inelastic terms beyond the first order and are equivalent and equally valid as perturbative analyses of the solution. Finally, the complication generated by the existence of first-order inelastic terms in [1–5,17] is replaced by a complicated computation of the phase shifts in [18], where they are modified beyond their KdV values. Motivated by the emergence of a pure inelastic term in the result of [18] (rather than the mixture of elastic and inelastic contributions contained in u1 q [17]), this paper offers a new approach to the construction of the solution Eq. (1.1) in the multiple-soliton case. The approach provides a systematic identification of pure inelastic contributions in every order of the expansion. The standard scheme [1–5], in which Eq. (1.9) holds, using the original physical coordinates, (t , x), and the observed wave numbers, will be adopted for the following reasons. To begin with, the first-order analysis is much simpler than that of higher orders, and yields the generic structure of inelastic terms that are generated in all orders. There are lessons to be learned, which may be obscured by the cumbersome nature of higher-order calculations (see Appendix I). Finally, this approach can be applied not only to the perturbed KdV and NLS equations, but also to equations that lack first-order asymptotic integrability (e.g., the perturbed Burgers [19,20] and mKdV [21] equations). In the multiple-soliton case, the solution of Eq. (1.1) is written as a sum of two components: elastic and inelastic. To construct the elastic component, one first obtains (Section 2) the series solution of Eq. (1.1) when the zero-order approximation is the single-soliton solution, uSingle , of the Normal Form (Eq. (1.3) or its higher-order updated version, Eq. (A.1)). Everywhere in that series, one then replaces uSingle by uMult , a multiple-soliton solution of the Normal Form. The result is identified as the elastic component (Section 3). This component does not generate changes in soliton parameters beyond the values they have in the zero-order term, uMult . The inelastic component is an O (ε ) effect. The first-order analysis generates a pure inelastic term in closed form. Its asymptotic shape coincides with the O (ε ) part Eq. (1.20). The latter describes the elastic collision of solitons and anti-solitons (negative amplitude solitons), which exchange signs upon collision (Sections 4, 5, 7 and 10). This is the generic type of inelastic terms, at least, through O (ε 3 ). Different initial data and/or boundary conditions generate other inelastic processes, such as soliton–anti-soliton creation or annihilation, soliton decay or merging, and inelastic soliton deflection (Sections 7–9). Invariably, the asymptotic profiles of these solitons and anti-solitons are proportional to those of the zero-order solitons of the Normal Form, Eq. (A.1). Hence, they do not generate changes in soliton parameters (wave numbers, velocities, phase shifts). Inelastic contributions are generated from parts in the perturbation that represent coupling between the KdV-solitons and non-KdV waves that are generated spontaneously by the perturbation (Section 6). Most of the second-order inelastic contributions to the solution can be written as differential polynomials in the zero-order approximation, u (Section 11). They represent two types of processes: the same soliton-anti soliton waves identified in the first-order analysis, or radiative tails that emanate from the origin and decay exponentially fast in all directions in the x–t plane. The obstacle to asymptotic integrability contained in the second-order perturbation [1–6] generates an inelastic contribution, which cannot be expressed as a differential polynomial, and is solved for numerically. This contribution contains the same soliton–anti-soliton wave encountered in all other terms, accompanied by a dispersive wave (Section 12). In the two-soliton case, this dispersive wave trails behind the slower of the two solitons. The dispersive wave has been observed numerically before [18,22,23]. This analysis identifies the specific term that generates it. The analysis uses extensively special polynomials. These are differential polynomials in u, the zero-order solution (of Eq. (A.1)), which vanish identically when u is a single-soliton solution. They are presented in Appendix F. By construction, the driving term in the dynamical

2990

Y. Zarmi / Physica D 237 (2008) 2987–3007

equation that determines the first-order inelastic contribution is a special polynomial. When u is a multiple-soliton solution, this driving term is localized around the soliton-collision region (a finite domain around the origin in the x–t plane) (Section 6). The solution it generates is also a special polynomial that represents a soliton–anti-soliton wave (Sections 5–7 and 10). In the second- (Sections 11 and 12) and third(Sections 13 and 14) order analyses, the driving terms that generate the inelastic contributions can be also constructed as localized special polynomials, allowing, again, for the description of the inelastic component of the solution in terms of soliton–anti-soliton waves. Finally, when the PKDVE is not integrable asymptotically, the expansion scheme presented here may be viewed as a transformation of Eq. (1.1) into a system of two equations: The Normal Form for the usual KdV solitons, and an auxiliary equation. The latter describes the contribution of obstacles to asymptotic integrability to the inelastic component. Its solution is a conserved quantity, which contains the dispersive wave that has been previously discovered [18,22,23] (Section 14). Numerical examples are provided for the two-soliton case, except for one three-soliton example (Section 10). The ideas presented in the following are applicable to solutions that involve more than two solitons, and an enormous richness of phenomena can be conceived of. The physical picture that emerges through O (ε 2 ), is unaltered when the analysis is extended to O (ε 3 ). 2. Solution of Eq. (1.1): Single-soliton case — no non-local terms In this section, the case when the zero-order approximation is uSingle , the single-soliton solution given by Eq. (1.5), is considered. The Normal Form, Eq. (A.1), and the terms that may appear in the expansion of the solution of Eq. (1.1) are presented in Appendices A–D. The velocity of the soliton is updated according to Eq. (A.5). Because of the simple structure of uSingle , coefficients of many of the terms enumerated in Appendix C remain free, to the extent that the solution can be obtained without invoking any terms that contain non-local quantities. This has been verified through third order. Hence, monomials that contain non-local terms are not included in the solution in the single-soliton case. The motivation is that the seriessolution of the single-soliton case will be used to construct the elastic component in the multiple-soliton case (Section 3), and non-local terms generate mixtures of elastic and inelastic contributions. Substituting Eq. (1.2) and the higher-order terms given in Appendices A–D in Eq. (1.1), and setting a2 , the coefficient of the non-local quantity in Eq. (1.4) to zero, the coefficients in the first-order correction to the solution are found to be (a superscript S has been appended to signify that these are the coefficients when u is a single-soliton solution): 5 10 5 5 (S ) a1 = − α1 + α2 + α3 − α4 , a(2S ) = 0, a(3S ) = −5α1 + 5α2 . 2 3 3 2 The series solution of Eq. (1.1) is Eq. (1.2), with a single-soliton zero-order approximation:

       
w (S ) uSingle = uSingle (t , x) + ε u(S )(1) uSingle + ε2 u(S )(2) uSingle + ε 3 u(S )(3) uSingle + O ε 4 .

(2.1)

(2.2)

For example, the first-order correction of Eq. (1.4) now has the following form: (S )

(S )

u(S )(1) uSingle = a1 uSingle,2 + a3



(S )



uSingle


2

,

(2.3)

(S )

with a1 and a3 given by Eq. (2.1). The coefficients obtained in the second- and third-order analyses are presented in Appendix E. There, as well, only local terms are included. Elimination of terms that contain non-local entities does not limit the choice of initial data or boundary conditions imposed on Eq. (1.1). The reason is that, in the single soliton case, these terms can be expressed as linear combinations of local terms. Consider, for example, the special polynomial R(3,1) (see Appendix F). It vanishes identically when u is a single-soliton solution. As a result, the u1 q-term is related to the local terms that appear in Eq. (1.4):

∂x R(3,1) [u] = ∂x {u1 + qu} = u2 + u2 + qu1

→ u→uSingle

0,


q = ∂x−1 u .

(2.4)

3. Multiple-soliton case: Two-component description and the elastic component When the zero-order approximation is a multiple-soliton solution, the full solution of Eq. (1.1) is written as a sum of an elastic and an inelastic component:

w (t , x) = wel [uMult ] + win [uMult ] ,

(3.1)

where uMult is the multiple-soliton zero-order approximation (solution of Eq. (A.1)). The elastic component, wel , is obtained from the series solution in the single-soliton case, Eq. (2.2). In that series, one replaces everywhere the single-soliton solution, uSingle , by uMult :


wel (t , x) = w(S ) [uMult ] = uMult (t , x) + ε u(S )(1) [uMult ] + ε 2 u(S )(2) [uMult ] + ε 3 u(S )(3) [uMult ] + O ε4 .

(3.2)

wel (t , x) preserves the elastic nature of soliton collisions because it contains only powers of uMult or of its spatial derivatives. As uMult breaks up asymptotically into a sum of well-separated pure single solitons, which undergo elastic collisions, so do all contributions in wel . Consequently, all the terms in wel have the same velocities, wave numbers and phase shifts of the zero-order solitons. 4. Multiple-soliton case: The first-order inelastic component Non-local terms in the full solution of Eq. (1.1) are accounted for by the inelastic component, win . In O (ε ), the inelastic term, (1) uin [u], is given by the difference between u(1) [u], the full first-order term of Eq. (1.4), and u(S )(1) [u], the elastic component, given by Eqs. (2.3) and (2.1) :

Y. Zarmi / Physica D 237 (2008) 2987–3007

2991

(1) uin [u] = u(1) [u] − u(S )(1) [u]

=−

10

 10 10 (α2 − α4 ) u2 + qu1 + u2 = − (α2 − α4 ) R(4,1) [u] = − (α2 − α4 ) ∂x R(3,1) [u] .

(4.1) 3 3 3 Here, and in the rest of the paper, u will represent a multiple-soliton solution. R(3,1) [u] is a special polynomial (see Appendix F). It vanishes if uSingle is substituted for u. Due to this property, it asymptotes into a purely inelastic term when u is an N-soliton solution: Each soliton is multiplied by wave numbers of other solitons. Using Eq. (D.10) for q(t , x), one finds: R

(3,1)

→ |t |→∞

N X

 Qi uSingle (t , x; ki ) ,

Qi =

 X

−2kj +

kj >ki

 k
i=1

j

X




i

 
 sgn(t ) . 2kj 

(4.2)

For N = 2, Eq. (4.1) becomes the inelastic term hidden in the formulation of [18] (see Eq. (1.20)): R(4,1) → ∂x 2k2 uSingle (t , x; k1 ) − 2k1 uSingle (t , x; k2 ) sgn(t ),




|t |→∞

(k2 > k1 ) .

(4.3)

5. Dynamical equation for the first-order inelastic component (1)

The first-order inelastic term, uin , has been found easily because Eq. (1.1) is integrable through O (ε ). It is still worthwhile to see how

(1)

uin is obtained as a solution of a dynamical equation, as this leads to the discovery of additional inelastic processes. Substituting Eqs. (1.2) and (1.3) in Eq. (1.1), the dynamical equation for the full first-order term, u(1) , is found to be:



∂t u(1) = 6∂x uu(1) + ∂x3 u(1) + 30α1 u2 u1 + 10α2 uu3 + 20α3 u1 u2 + α4 u5 . Writing u u

(1)

(1)

(5.1)

as a sum of the elastic and inelastic terms, both functionals of u:

= u(S )(1) [u] + u(in1) [u] ,

(5.2)

and using Eqs. (2.1) and (2.3) for u

(S )(1)

(1)

, Eq. (5.1) yields the following equation for uin :

∂t u(in1) = 6∂x uu(in1) + ∂x3 u(in1) + 10 (α2 − α4 ) R(7,1) [u]
R(7,1) [u] = ∂x R(6,1) [u] , R(6,1) = u3 + uu2 − u21 . 



(5.3)

R(7,1) and R(6,1) are special polynomials of scaling weights 7 and 6, respectively (see Appendix F). 6. Nature of driving term in Eq. (5.3) R(7,1) and R(6,1) are localized in the x–t plane when u is a multiple-soliton solution. This can be verified by direct substitution of any multiple-soliton solution, or demonstrated in two other ways. The first argument [24] is based on writing, R(7,1) in terms of symmetries (see Eqs. (A.2) and (A.3)): R(7,1) = G1 [u] S2 [u] − G2 [u] S1 [u] .

(6.1)

Away from soliton trajectories, R(7,1) falls off exponentially, because so do all entities in Eq. (6.1). However, thanks to the asymptotic behavior of Sn and Gn (Eq. (A.3)), R(7,1) vanishes asymptotically also along soliton trajectories. Hence, it is localized around the solitoncollision region. The second argument [25] is based on writing R(7,1) as R(7,1) [u] = u2 us ,

us = ∂x



G2



u

.

(6.2)

A detailed study shows that when u is an N-soliton solution, then us is composed of N interlaced non-KdV-solitons and anti-solitons. Each of these has a KdV-like profile. However, their amplitudes, wave numbers and velocities do not have KdV values. Each of these new signals propagates between two adjacent KdV solitons, and, hence, us overlaps with the KdV solution, u, only in the soliton collision region. As a result, the product in Eq. (6.2) is localized around that region. For example, when u is a two-soliton solution, given by Eq. (1.6), the asymptotic form of us is: A1

us −−−→

(cosh {K1 (x + V1 t )} + ξ1 )2 In Eq. (6.3), one has, for k2 > k1 , |t |→∞

+

A2

(cosh {K2 (x + V2 t )} + ξ2 )2

A1 = −sgn(t ) · 2K12 K2 ,

A2 = sgn(t ) · 2K1 K22

K1 = k1 + k2 ,

+ 3K22 ,

V1 =

K12

K2 = k2 − k1 ,

V2 = 3K12 + K22 ,

.

(6.3)

(6.4)

us asymptotes into a soliton and an anti-soliton that exchange signs upon an elastic collision. Their velocities and wave numbers do not obey the KdV dispersion relation. Hence, the driving term in Eq. (5.3) may be regarded as coupling between the multiple-KdV-soliton wave, u, and us , a non-KdV wave, generated spontaneously by the perturbation. These new solitons and the KdV solitons overlap only near the origin. As a result, R(7,1) is localized around the origin. Examples of u, R(6,1) , R(7,1) and us in the two-soliton case are shown in Figs. 1–4, respectively. The rate of exponential fall off of R(7,1) is discussed in Appendix I.

2992

Y. Zarmi / Physica D 237 (2008) 2987–3007

Fig. 1. Two-soliton solution (Eq. (1.6)), k1 = 0.1, k2 = 0.2.

Fig. 2. Local special polynomial, R(6,1) (Eq. (5.3)), for two-soliton case. Parameters as in Fig. 1.

Fig. 3. Local special polynomial, R(7,1) (Eq. (5.3)), for two-soliton case. Parameters as in Fig. 1.

7. Standard solution of Eq. (5.3) The fact that R(7,1) is a complete differential allows for simplification of the analysis. Defining (1)

uin = 10 (α2 − α4 ) ∂x ω,

(7.1)

one finds from Eq. (5.3) that the dynamical equation for ω is

∂t ω = 6u∂x ω + ∂x3 ω + R(6,1) [u] .

(7.2)

The homogeneous part of Eq. (7.2) is the equation adjoint to the homogeneous part of Eq. (5.3). The special polynomial, R(3,1) (see Appendix F) is a particular solution of Eq. (7.2): 1

1

3

3

ωSA = − R(3,1) [u] = − {u1 + qu} .

(7.3)

Y. Zarmi / Physica D 237 (2008) 2987–3007

2993

Fig. 4. Non-KdV wave, us , (Eq. (6.2)), k1 = 0.5, k2 = 0.75.

Fig. 5. First-order soliton–anti-soliton scattering wave, ωSA (Eq. (7.3)). Parameters as in Fig. 1.

Eq. (7.3) is in agreement with Eq. (4.1). In the two-soliton case, ωSA asymptotes into a linear combination of the zero-order single-soliton solutions of Eq. (A.1), the same as in Eq. (4.3):

ωSA −→ |{z} −

2 3

k2 uSingle [t , x − η1 ; k1 ] − k1 uSingle [t , x − η2 ; k2 ] sgn(t ),




(k2 > k1 ) .

(7.4)

|t |→∞

In Eq. (7.4), uSingle is given by Eq. (1.5) and the phase shifts - by:

ηi =

1 − (sgn(t ))i




2

 log

k2 + k1

  ki

k2 − k1

,

i = 1, 2.

(7.5)

The asymptotic form, Eq. (7.4) is linear in the zero-order single-soliton solutions. Hence, the addition of ωSA to the solution does not lead to modification of soliton velocity or phase shift. Moreover, Eq. (7.4) represents a purely inelastic term: Each soliton is multiplied by the wave number of the other soliton. A plot of ωSA is shown in Fig. 5. 8. Inelastic solutions of homogenous part of Eq. (7.2) Consider the homogeneous part of Eq. (7.2),

∂t ω = 6u∂x ω + ∂x3 ω.

(8.1)

Gn [u] of Eqs. (A.2) and (A.3) constitute an infinite hierarchy of solutions of Eq. (8.1) [12–16,26–28]. Employing appropriate linear combinations of Gn [u], one can construct solutions of Eq. (8.1) that correspond to inelastic processes, the nature of which is different from that of ωSA of Eq. (7.3). When u is an N-soliton solution. one can generate ‘‘amputated’’ solutions that contain any subset of these solitons, with arbitrary (positive or negative) amplitudes. This is a simple consequence of the asymptotic form of Gn [u] (Eq. (A.4)). For example, in the case of the two-soliton solution (Eq. (1.6)), the amputated solutions of Eq. (8.1) are: u˜ 1 = 4k22 G1 − G2 ,




u˜ 2 = G2 − 4k21 G1 .




(8.2)

In Eq. (8.2), each u˜ i , i = 1, 2, contains only soliton i with the phase shift given by Eq. (7.5), although G1 = u contains both solitons. Fig. 6 provides examples of u˜ 1 . Finally, the solution of Eq. (8.1), given by the linear combination u˜ = α G1 + β G2 ,

(8.3)

2994

Y. Zarmi / Physica D 237 (2008) 2987–3007

Fig. 6. ‘‘Amputated’’ solution of Eq. (8.1), u˜ 1 (Eq. (8.2)). Parameters as in Fig. 1.

Fig. 7. Soliton-anti-soliton annihilation (Eq. (9.1)). Parameters as in Fig. 1.

contains both solitons with arbitrary amplitudes. Their amplitudes will have opposite signs for

α · β < 0,

4k21 |β| < |α| < 4k22 |β| ,

(k2 > k1 ) .

(8.4)

In summary, the homogeneous part of the first-order dynamical equation (in fact, any order) can generate a plethora of inelastic processes when the zero-order approximation is an N-soliton solution, by appropriate adjustment of the initial data and/or boundary conditions. Finally, the waves of Eqs. (8.2) and (8.3) are solutions of Eq. (8.1) but not of the KdV equation. 9. Additional inelastic processes of physical interest: Two-soliton case Additional inelastic processes can be constructed as linear combinations of ωSA of Eq. (7.3) and solutions of the homogeneous equation, Eq. (8.1). The following are examples in the two-soliton case. The first example is the creation and annihilation of a soliton–anti-soliton pair:

ωCreation = ωSA + κ1 G1 [u] + κ2 G2 , ωAnnihilation = ωSA − κ1 G1 [u] − κ2 G2 !
2 k21 + k22 − k1 k2 1 κ1 = − , κ2 = . 3 (k2 − k1 ) 6 (k2 − k1 )

(9.1)

Their asymptotic structure is

ω −→ |{z} ±

4 3

k2 uSingle [t , x; k1 ] − k1 uSingle [t , x; k2 ] .




(9.2)

|t |→∞

The + and − signs in Eq. (9.2) correspond to annihilation and creation, respectively. Examples of these processes are shown in Figs. 7 and 8, respectively. The driving term in Eq. (7.2) acts like a ‘‘beacon’’; it ‘‘radiates’’ or ‘‘absorbs’’ the soliton–anti-soliton wave. Conversely, one has

ωSA =

1 2

(ωCreation + ωAnnihilation ) .

(9.3)

Fig. 9 shows a solution, in which soliton no. 1 is deflected into soliton no. 2, given by:

ω1→2 = −ωSA + λ1 G1 [u] + λ2 G2 ,

λ1 =

2 k21 + k22 + k1 k2 3 (k2 + k1 )


, λ2 = −

1 6 (k2 + k1 )

! .

(9.4)

Y. Zarmi / Physica D 237 (2008) 2987–3007

2995

Fig. 8. Soliton-anti-soliton creation (Eq. (9.1)). Parameters as in Fig. 1.

Fig. 9. Inelastic soliton deflection (Eq. (9.4)). Parameters as in Fig. 1.

Fig. 10. Soliton decay (Eq. (9.5)). Parameters as in Fig. 1.

Finally, up to overall multiplicative factors, the soliton decay and merging processes are given by:

ω1→1+2 = ωSA + ω1+2→2 = ωSA +



2 3



2 3



k1 − 4k22 β G1 [u] + β G2 [u] ,



k2 − 4k21 β G1 [u] + β G2 [u] ,

β<− β>−

1 6 (k1 + k2 ) 1 6 (k2 − k1 )

(9.5)

.

In Eq. (9.5), k2 > k1 . Examples are shown in Figs. 10 and 11, respectively. 10. One example of a three-soliton process In this section, one example of the three-soliton case is presented. The standard solution of Eq. (7.2) is Eq. (7.3), with u- a three-soliton solution, represented by the Hirota formula as [11]: u (t , x) = 2∂x2 ln f (t , x) ,

p Mij = δij + 2

k1 kj

ki + kj

f (t , x) = det M

g (t , x; i) g (t , x; j) ,

g (t , x; m) = exp km x + 4k2m t







, 1 ≤ i, j , m ≤ 3 .

(10.1)

2996

Y. Zarmi / Physica D 237 (2008) 2987–3007

Fig. 11. Soliton merging (Eq. (9.5)). Parameters as in Fig. 1.

Fig. 12. R(6,1) (Eq. (5.3)) in three-soliton case. k1 = 0.2, k2 = 0.15, k3 = 0.1.

To obtain a solution in which the ‘‘beacon’’ radiates waves, a particular solution of the homogeneous equation, Eq. (8.1) has been added to ωSA of Eq. (7.3):

ωCreation = ωSA + κ1 G1 [u] + κ2 G2 [u] + κ3 G3 [u] .

(10.2)

The coefficients κi , i = 1, 2, 3, are determined by the conditions required for zero-initial data:



2 2 κ1 + 4k21 κ2 + 4k21 κ3 = − (k2 + k3 ) 3

κ1 +

4k22




κ2 +

2 4k22




κ3 =



2 κ1 + 4k23 κ2 + 4k23 κ3 =

2 3 2

(k1 − k3 ) ,

(k1 < k2 < k3 )

(10.3)

(k1 + k2 ) . 3 Using Eq. (D.11), the asymptotic form of the solution, Eq. (10.2), is found to be 4 ωCreation (t , x) −→ |{z} 3 − (k2 + k3 ) uSingle [k1 ] + (k1 − k3 ) uSingle [k2 ] + (k1 + k2 ) uSingle [k3 ] .

(10.4)

t →+∞

As the solution is a special polynomial that vanishes identically in the single-soliton case, its asymptotic form represents a purely inelastic term: Each soliton is multiplied by wave numbers of the other solitons. Again, as Eq. (10.4) is a linear combination of the zero-order single-soliton solutions of the Normal Form, Eq. (A.1), no modification of soliton parameters is generated. R(6,1) , the driving term in Eq. (7.2) in the three-soliton case, is shown in Fig. 12 and the solution, Eq. (10.2) — in Fig. 13. ( 2)

11. Second-order inelastic term, uin : Differential-polynomial part Except for a dispersive wave that is generated by obstacles to asymptotic integrability, which will be studied in second order, the higherorder calculations yield the same physical picture that has been realized in the first-order analysis: Invariably, the asymptotic structure of inelastic terms corresponds to the same processes that have been identified in the first-order analysis. Sections 11 and 12 present the results of the second-order analysis. (2) The dynamical equation that determines the second-order inelastic component, uin , has the same structure as the first-order equation (Eq. (5.3)):

  ∂t u(in2) = 6∂x uu(in2) + ∂x3 u(in2) + R(2) [u] .

(11.1)

Y. Zarmi / Physica D 237 (2008) 2987–3007

2997

Fig. 13. One soliton and two anti-solitons creation (Eq. (10.2)) in three-soliton case. Parameters as in Fig. 12.

The driving term, R(2) [u], depends on the structure of the second-order perturbation in Eq. (1.1), of the first- and second-order elastic and (1) (2) (1) components (uel and uel , respectively), as well as of the first-order inelastic component, uin . (1)

(2)

(1)

uel is given by Eqs. (2.1) and (2.3), and uel by Eqs. (C.2), (E.1) and (E.3). For the first-order inelastic component, uin , the standard solution of Eq. (4.1) has been chosen. By construction, these choices guarantee that R(2) [u] vanishes identically if one substitutes for u the single(2) soliton solution, Eq. (1.5), as direct substitution confirms. This implies that one ought to construct uin [u] so that it also vanishes identically if the single-soliton solution is substituted for u. This task is possible thanks to the high level of freedom inherent in the perturbation expansion. One complication is expected owing to the existence of an obstacle to asymptotic integrability in the second-order analysis [2–6]. In (2) (2) general, R(2) contains a part, which drives a contribution to uin that cannot be written as a differential polynomial in u. Therefore, uin will be written as: (2)

(2)

uin [u] = u˜ in [u] + η(2) [u] .

(11.2)

(2)

(2)

In Eq. (11.2), u˜ in [u] is the differential polynomial part in uin [u], and η(2) [u] is the contribution that cannot be written as such a polynomial. (2)

(2) One can make uin [u] vanish identically in the single-soliton case by making u˜ in [u] and η(2) [u] vanish separately. The construction of (2)

u˜ in [u] is presented in this section, and that of η(2) [u] — in Section 12. (2)

(2)

The expression for u˜ in [u] is as in Eq. (C.2), with coefficients denoted by bin,k . To ensure that u˜ in [u] vanishes identically when u is a single-soliton solution, bin,k must obey the relations given in Eq. (G.1). The results reported in the following are based on bin,k that obey these relations. Eq. (11.1) provides additional constraints on bin,k . Still, ample freedom remains in their determination. The freedom has been exploited (2)

so as to reduce the number of non-local terms as much as possible, and to simplify the dynamical equation for η(2) [u]. In particular, u˜ in [u] becomes a sum of 18 special polynomials, each of which vanishes if uSingle is substituted for u. Each polynomial is multiplied by a coefficient that is quadratic in the first-order coefficients αi , or linear in the second-order coefficients, βj , of Eq. (1.1). These 18 terms represent two inelastic processes: Radiative tails and soliton–anti-soliton waves of the type encountered in the first-order analysis. The values chosen for bin,k are given in Appendix H. Radiative tails are generated by local special polynomials, i.e., polynomials that do not contain non-local entities of the types (2) enumerated in Appendix D. In the scaling weight of u˜ in [u], 6, there is only one local special polynomial, R(6,1) (see Eq. (F.5)). The radiative tail it generates is:



2 75α1 − 75α1 α2 +



−

140 3

100 9

β1 +

α22 −

140 3

50 3

β2 +

α1 α3 − 280 3

350 9

α2 α3 +

β3 − 98β4 −

200

140 3

9

α32 − 25α1 a4 +

β5 + 70β6 −

140 3

625 9

α2 α4 −

β7 + 28β8

 α42  (6,1) 9 . R

200

(11.3)

The remaining 17 terms are non-local polynomials. They all asymptote into soliton–anti-soliton waves. Hence, only two representative terms are shown here:

α3 α4

50 9

∂x3 R(3,1) [u] −→ |{z} α3 α4

100 9

 ∂x3 k1 uSingle [k2 ] − k2 uSingle [k1 ] sgn(t ),

|t |→∞

−

14 3

− + −

(β4 − β8 ) q (6uu1 + u3 ) − 7 (β1 + 2β2 − 6β3 + β4 + 4β5 − 6β6 + 4β7 ) q(3) u1 14 3 14 3 14 3

(10β1 − 10β2 − 20β3 + 21β4 + 10β5 − 15β6 + 10β7 − 6β8 ) u3 (5β1 − 20β2 + 10β3 + 7β4 − 10β5 + 15β6 − 10β7 + 3β8 ) u21 (12β1 − 6β2 − 32β3 + 27β4 + 18β5 − 27β6 + 18β7 − 10β8 ) uu2

(11.4)

2998

Y. Zarmi / Physica D 237 (2008) 2987–3007

Fig. 14. Second-order wave driven by obstacle to asymptotic integrability (Eq. (12.4)). Parameters as in Fig. 1.

7

(β1 + 2β2 − 6β3 + 3β4 + 4β5 − 6β6 + 4β7 − 2β8 ) u4 
 56 2 (β4 − β8 ) k22 k1 uSingle [k2 ] − k21 k2 uSingle [k1 ]
−→ |{z} 3 ∂x + (β1 + 2β2 − 6β3 + β4 + 4β5 − 6β6 + 4β7 ) k31 uSingle [k2 ] − k32 uSingle [k1 ] sgn(t ). −

3

(11.5)

|t |→∞

For the definition of q(3) , see Eq. (D.2). The asymptotic forms are for k2 > k1 . ( 2)

12. Second-order inelastic term, uin : Non-polynomial part

η(2) [u], the non-polynomial term in Eq. (11.2), is driven by the obstacle to asymptotic integrability, contained in driving term of Eq. (11.1) [24,25]. Because of the freedom in the choice of the coefficients bin,k , the functional form of the obstacle is not unique. With the choice of Eqs. (H.1)–(H.8), the obstacle to integrability becomes a gradient of a local special polynomial. The resulting equation for η(2) [u] is:
∂t η(2) = 6∂x uη(2) + ∂x3 η(2) + µ∂x R(8,1) + O (ε) . In Eq. (12.1), R 5

µ=

3

(8,1)

(12.1)

is a scaling-weight 8 local special polynomial (see Eq. (F.7)), and

3α1 α2 + 4α22 − 18α1 α3 + 60α2 α3 − 24α32 + 18α1 α4 − 67α2 α4 + 24α42




+ 7 (3β1 − 4β2 − 18β3 + 17β4 + 12β5 − 18β6 + 12β7 − 4β8 ) . Unlike the obstacle, the coefficient, µ, is unique, and has been found in previous works [4–6,25].

(12.2)

Focusing on the lowest-order part of Eq. (12.1), the analysis can be simplified by defining

η(2) [u] = µ∂x ω(2) [u] . The equation for ω(2) becomes

(12.3)

∂x ω(2) = 6u∂x ω(2) + ∂x3 ω(2) + R(8,1) .

(12.4) (2)

(2)

There are no differential polynomial solutions for either η [u] or ω [u]. (Thus, if µ 6= 0, then the driving terms in Eqs. (12.1) and (12.4) are, indeed, obstacles to asymptotic integrability [2–6]). Despite the fact that the driving term in Eq. (12.4) does not obey the conditions of the Fredholm Alternative Theorem, its solution is bounded. (See Appendix I for a discussion of this point.) Moreover, it vanishes for |x| → ∞. The reason is that the driving term in Eq. (12.4) is a local special polynomial (see Eq. (F.7)). It vanishes identically if the single-soliton solution is substituted for u. Consequently, it is localized around the origin in the x–t plane, and vanishes exponentially fast in all directions away from the origin in the multiple-soliton case. Moreover, it does not resonate with the homogeneous part of Eq. (12.4). Hence, it does not generate an unbounded contribution in ω(2) . This is borne out in an example: Eq. (12.4) was solved numerically, with zero-initial data for t → −∞. Fig. 14 shows the solution when u is a two-soliton solution with k1 = 0.1 and k2 = 0.2. It is comprised of a soliton and an anti-soliton, accompanied by a decaying dispersive wave. Within the accuracy of the numerical procedure, the soliton and anti-soliton have the same parameters (velocities, wave numbers and phase shifts) as the zero-order solitons in u. The velocities of the original KdV solitons (vi = 4k2i ), are equal to 0.04 and 0.16, respectively. In the range 100 ≤ t ≤ 1600, the velocities of the anti-soliton and the soliton were found to be 0.03991 ± 0.00008 and 0.16001 ± 0.00003, respectively. Up to overall amplitudes, determined from the numerical solution, the soliton and the anti-soliton are indistinguishable from profiles computed through Eqs. (1.5) and (7.5). Thus, like all other first- and second-order contributions, the second-order obstacle does not seem to lead to a modification in soliton parameters. The dispersive wave trails behind the slower of the two solitons. Away from that soliton, this wave falls off exponentially as exp(−α · x). α decreases with time. For example, at t = 800 and 1600 it has the values 0.02817 and 0.01469, respectively. To ascertain that this wave is not an artifact of the numerical scheme, Eq. (12.4) was solved numerically again, after reflection, x → −x. This would have allowed for

Y. Zarmi / Physica D 237 (2008) 2987–3007

2999

Fig. 15. Dispersive component in Fig. 14, t = 1600. Parameters as in Fig. 1.

the generation of a trailing wave behind the faster soliton. As this did not happen, the dispersive wave seems to be a genuine effect. Fig. 15 shows the dispersive wave for t = 1600 after the soliton–anti-soliton contribution has been subtracted. Such a dispersive wave had been found in previous numerical works [18,22,23]. The present work has identified its origin: the obstacle to asymptotic integrability, with a particularly simple form of the obstacle. Thanks to the fact that ω(2) is bounded, so is η(2) [u], as Eq. (12.3) implies. In fact, η(2) [u] must vanish for |x| → ∞. Consequently, in this order of the expansion, η(2) [u] is a conserved quantity: +∞

Z

d dt

η(2) [u] dx = 0 + O (ε) .

(12.5)

−∞

There will be higher order contributions to Eq. (12.1). Its O (ε ) part is determined in the analysis of the O (ε 3 ) contribution to the full solution of Eq. (1.1). This will be presented in Section 14. ( 3)

13. Third-order inelastic term, uin : Differential-polynomial part The results of the third-order analysis are similar to the second-order ones, with many more terms. Only the important features of the analysis are described in the following through examples of typical terms. The dynamical equation for the third-order inelastic (3) contribution, uin [u] is:

  ∂t u(in3) = 6∂x uu(in3) + ∂x3 u(in3) + R(3) [u] .

(13.1)

(3)

The structure of the driving term, R [u], is a consequence of the calculations in all-previous orders. It has scaling weight 11. Therefore, (3) uin [u] has scaling weight 8. By construction, because all previous orders were based on special polynomials, R(3) [u] vanishes identically if the single-soliton solution is substituted for u. This is verified by direct substitution. R(3) [u] also contains obstacles to asymptotic (3) integrability [2–5]. Therefore, as in Eq. (11.2), uin [u] is written as: (3)

(3)

uin [u] = u˜ in [u] + η(3) [u] .

(13.2)

(3)

(3)

The polynomial component, u˜ in [u], is a linear combination of special polynomials, and η [u] is a non-polynomial contribution, driven by the obstacles to asymptotic integrability contained in R(3) [u]. As in the second-order analysis, the two components in Eq. (13.2) are constructed so that each vanishes identically if the single-soliton solution is substituted for u. (3) In this section, results are presented for the polynomial component, u˜ in [u]. It has the structure given in Eq. (C.3), with the coefficients (3)

denoted by cin,k . To guarantee that u˜ in [u] vanishes identically when the single-soliton solution is substituted for u, cin,k have to obey the relations given in Eqs. (G.2)–(G.5). The results presented in the following, are based on cin,k that obey these relations. As in the second-order analysis, the additional constraints that Eq. (13.1) provides leave ample freedom in the determination of cin,k . The choice of coefficients (3)

was aimed at reducing the number of non-local entities in u˜ in [u], at constructing it as a sum of special polynomials (which each vanish in the single-soliton case) and at simplifying the expression for η(3) [u]. Of the 74 coefficients, 54 turn out to vanish and 20 are polynomials in the coefficients that appear in Eq. (1.1). There are three independent local special polynomials of scaling weight 8 (see Eq. (F.7)). They are all localized around the soliton(3) collision region in the x–t plane. Hence, the radiative tail contribution in u˜ in [u] is of the form: A8,1 R(8,1) [u] + A8,2 R(8,2) [u] + A8,3 R(8,3) [u] .

(13.3)

A8,k are known polynomials in the coefficients that appear in Eq. (1.1) (cubic in αi , bilinear in αi and βj , and linear in γk ). There are dozens of terms that are non-local polynomials, namely, special polynomials that contain non-local entities. Invariably, they all asymptote into soliton–anti-soliton contributions. Hence, only one term is shown:


−200 2q(5,1) − q(5,3) u1 + 186q(3) S2 [u] − 22qS3 [u] 229u4 + 448uu21 + 94u2 u2 − 215u22 + 217u1 u3 + 20uu4 5 (13.4) 

1792 20k31 − 11k1 k22
uSingle (k2 ) → γ12 k22 − k21 ∂x sgn(t ). 2 3 − 20k2 − 11k2 k1 uSingle (k1 ) |t →∞| 5 The asymptotic form is for k2 > k1 . For the definition of symmetries and non-local entities, see Appendices A and D, respectively. −

28

γ12



3000

Y. Zarmi / Physica D 237 (2008) 2987–3007

14. Equivalent system of two coupled equations Following the second-order analysis, one is tempted to solve for the non-polynomial component, η(3) [u], as a separate entity. Instead, it will be treated as a higher order effect on the dynamical equation for the second-order non-polynomial term, η(2) [u]. Namely, η(3) [u] determines the O (ε )-contribution to Eq. (12.1). The freedom in the choice of the coefficients cin ,k in the third-order polynomial component, (3) u˜ in [u], enables one to update Eq. (12.1) into: 6 uη(2) + ∂x3 η(2) + µR(8,1) + ε A η(2)  
10  (8,3)  + O ε2 . [u] ∂t η(2) = ∂x  q − α µ∂ R (α ) 4 2 x +ε  3
  +σ7,1 ∂x3 R(7,1) [u] + ∂x2 σ8,2 R(8,2) [u] + σ8,3 R(8,3) [u]












 

(14.1)

In Eq. (14.1), µ is the coefficient of the second-order obstacle to integrability (see Eq. (12.2)), the three coefficients σ7,1 , σ8,2 and σ8,3 are independent polynomials in the coefficients that appear in Eq. (1.1) (cubic in αi , bilinear in αi and βj , and linear in γk ) and A η(2) = 10 (α2 + 2α4 ) u2 − 20 (α2 − α4 ) qu1 − 5 (3α1 − 2α2 − 2a3 − α4 ) u2 η(2)








− 10 (α2 − 2α3 ) u1 ∂x η(2) + 10α2 u∂x2 η(2) + α4 ∂x4 η(2) .

(14.2)

The following points are worth noting regarding Eq. (14.1). First, as found in previous works [2–5], depending on the values of the αi , βj , and γk , there may be up to three third-order obstacles to asymptotic integrability (a non-zero value for each of σ7,1 , σ8,2 and σ8,3 amounts to one obstacle). Second, the O (ε )-part in Eq. (14.1) is linear in η(2) . Terms quadratic in η(2) will appear only in the next order. The most important feature of Eq. (14.1) is that Eq. (12.5) can be extended to O (ε ). The solutions of Eq. (14.1) are bounded and obey a conservation law. The reasons are the following. Like the zero-order terms in Eq. (14.1), the O (ε )-driving terms are all localized in the x–t plane, as they are all proportional to local special polynomials (see Appendix F). Moreover, through O (ε ), Eq. (14.1) is a continuity equation for η(2) . Thus, validity of Eq. (12.5) can be extended through O (ε ). Solutions for η(2) , which vanish for |x| → ∞, obey a conservation law: d dt

Z

+∞


η(2) (t , x) dx = 0 + O ε 2 .

(14.3)

−∞

This last result suggests a new interpretation for the perturbation scheme presented in this paper: That when the perturbed KdV equation is not asymptotically integrable, the expansion may be viewed as a transformation of Eq. (1.1) into a system of two coupled equations, the Normal Form (Eq. (A.1)) for the ordinary solitons, and Eq. (14.1) for the effect of the obstacles to asymptotic integrability. The KdV soliton solutions of the Normal Form obey the well-known infinite hierarchy of conservation laws [1–5,12–16], and the effect of the obstacles generates a new conserved quantity. This interpretation is motivated by the fact that the obstacles generate a new phenomenon that does not have the localized nature of the soliton solution. The examples presented in this paper show that, depending on the initial data and boundary conditions imposed on Eq. (1.1), the inelastic components may offer a rich variety of inelastic processes for N ≥ 2 solitons (in fact, in any order of the expansion). Finally, the results of this analysis suggest that, apart from a higher-order decaying dispersive wave, the effect of inelastic contributions does not modify soliton parameters beyond the values determined by the Normal Form. This qualitative conclusion is in agreement with the results of the numerical experiment performed in the case of the ion-acoustic wave equations in [29]. Acknowledgments Helpful discussions with G. Burde and L. Prigozhin, I. Rubinstein and B. Zaltzman are deeply acknowledged. Appendix A. Normal form expansion of the solution of perturbed KdV equation Through third order, the series-solution of Eq. (1.1) is given in Eq. (1.2). The zero-order term, u(t , x), obeys the Normal Form: ut = 6uu1 + u3 + εα4 S3 [u] + ε 2 β8 S4 [u] + ε 3 γ14 S5 [u] + O ε 4 .




(A.1)

Sn [u] are the symmetries of the KdV equation. They are known, and can be computed from the first member of the hierarchy of symmetries through a recursion relation [4,5,12–16,26–28]. The symmetries relevant for this paper are: S 1 = ux S2 = 6uu1 + u3 S3 = 30u2 u1 + 10uu3 + 20u1 u2 + u5 S4 = 140u3 u1 + 70u2 u3 + 280uu1 u2 + 14uu5 + 70u31 + 42u1 u4 + 70u2 u3 + u7

(A.2)

S5 = 630u4 u1 + 1260u (u1 )3 + 2520u2 u1 u2 + 1302u1 (u2 )2 + 420u3 u3 + 966 (u1 )2 u3 up ≡ ∂xp u .




The entities Gn = ∂x−1 Sn , will be required as well. The first three are G1 = u G2 = 3u2 + u2 3

G3 = 10u + 10uu2 +

(A.3) 5u21

+ u4 .

Y. Zarmi / Physica D 237 (2008) 2987–3007

3001

The following property of Sn and Gn will be exploited. Along each soliton, away from the soliton collision region [12–16,25], Sn →

N X

(4ki )n−1 ∂x uSingle (t , x; ki )

i=1

Gn →

N X

(A.4)

(4ki )

n −1

uSingle (t , x; ki ) .

i =1

The sums in Eq. (A.4) run over all solitons in an N-soliton solution. A soliton solution of the Normal Form, Eq. (A.1) will have precisely the same structure as the solution of the unperturbed KdV equation, with the velocity of each soliton updated according to




2
3 vi = 4k2i 1 + εα4 4k2i + ε 2 β8 4k2i + ε 3 γ14 4k2i + O ε 4 .

(A.5)

Appendix B. Scaling weight The perturbed KdV equation as well the dynamical equations derived in every order of the analysis can be characterized by a scaling weight. The latter is determined by assigning u, ∂t and ∂x the scaling weights +2, +3 and +1, respectively [10,26]. Corresponding to these, integration by x has scaling weight −1, and integration by t has weight −3. With these assignments, all terms in the KdV equation have scaling weight 5; the nth-order perturbation in Eq. (1.1) is constructed solely out of terms with scaling weight (2n + 5). The scaling weight of non-local terms is determined accordingly. For example, the non-local term q(t , x) = ∂x−1 u(t , x) has a scaling weight of +1. When u is an N-soliton solution, a differential monomial in u is a sum of terms that are proportional to products of powers of the wave PN m m m numbers, ki , 1 ≤ i ≤ N of the solitons: k1 1 k2 2 · · · kN N . All these terms have the same total power, Ω = i=1 mi . This is the scaling weight of the monomial. For example, when u is the single-soliton solution of Eq. (1.5), u has scaling weight 2, and is proportional to k2 . ∂x u has scaling weight 3, and is proportional to k3 . Appendix C. Structure of higher-order corrections in Eq. (1.2) In the higher-order corrections in Eq. (1.2), uk denotes ∂xk u, the non-local terms are denoted by q, with superscripts that specify their scaling weights (see Appendix D) and η(n) are terms that cannot be written as differential polynomials in u. The latter are driven by obstacles to asymptotic integrability that emerge from second order onwards [2–6,25]. Inclusion of many non-local terms is required in order to express the dynamical equations in terms of special polynomials (see Appendix F). The expressions used are: u(1) = a1 u2 + a2 u1 q + a3 u2 u

(2)


q (t , x) = ∂x−1 u (t , x) , 2

3

(3)

(C.1) 4

(3)

(4)

= b1 u4 + b2 u3 q + b3 u2 q + b4 u1 q + b5 u1 q + b6 uq + b7 uqq + b8 uq + b9 uu2 + b10 u21 + b11 uu1 q + b12 u2 q2 + b13 u3 + b14 q(6,1) + b15 q(6,2) + b16 q(6,3) + b17 qq(5,1) + b18 qq(5,2) + b19 qq(5,3) + b20 q2 q(4)
2 + b21 q3 q(3) + b22 q(3) + b23 q6 + η(2) ,

(C.2)

u(3) = c1 u6 + c2 u5 q + c3 u4 q2 + c4 u3 q3 + c5 u3 q(3) + c6 u2 q4 + c7 u2 qq(3) + c8 u2 q(4) + c9 u1 q5

+ c10 u1 q2 q(3) + c11 u1 qq(4) + c12 u1 q(5,1) + c13 u1 q(5,2) + c14 u1 q(5,3) + c15 uq6 + c16 uq3 q(3) + c17 uq2 q(4) + c18 uqq(5,1) + c19 uqq(5,2) + c20 uqq(5,3) + c21 uq(6,1) + c22 uq(6,2) + c23 uq(6,3)
2 + c24 u q(3) + c25 uu4 + c26 u1 u3 + c27 (u2 )2 + c28 uu3 q + c29 u1 u2 q + c30 uu2 q2 + c31 u21 q2 + c32 uu1 q3 + c33 uu1 q(3) + c34 u2 q4 + c35 u2 qq(3) + c36 u2 q(4) + c37 u2 u2 + c38 uu21 + c39 u2 u1 q + c40 u3 q2 + c41 u4 + c42 q(8,1) + c43 q(8,2) + c44 q(8,3) + c45 q(8,4) + c46 q(8,5) + c47 q(8,6) + c48 q(8,7) + c49 qq(7,1) + c50 qq(7,2) + c51 qq(7,3) + c52 qq(7,4) + c53 qq(7,5)
2 + c54 qq(7,6) + c55 qq(3) q(4) + c56 q2 q(6,1) + c57 q2 q(6,2) + c58 q2 q(6,3) + c59 q2 q(3) + c60 q3 q(5,1) + c61 q(3) q(5,1) + c62 q3 q(5,2) + c63 q(3) q(5,2) + c64 q3 q(5,3) + c65 q(3) q(5,3) + c66 q4 q(4)
2 + c67 q(4) + c68 q8 + c69 q5 q(4) + c70 qq(7,7) + c71 q(8,8) + c72 q(8,9) + c73 q(8,10) + c74 q(8,11) + η(3) (t , x) .

(C.3)

All bk , 14 ≤ k ≤ 23, and all ck , 42 ≤ k ≤ 74, have been set to zero. The freedom in the expansion allows this choice. Their inclusion turns out to be unnecessary and only complicates the analysis. Appendix D. Bounded non-local quantities The non-local quantities used in the analysis must be bounded throughout the x–t plane for soliton solutions. They are listed according to their scaling weights. Scaling weight 1 q = ∂x−1 u.

(D.1)

Scaling weight 3


q(3) = ∂x−1 u2 .

(D.2)

3002

Y. Zarmi / Physica D 237 (2008) 2987–3007

Scaling weight 4


q(4) = ∂x−1 u2 q .

(D.3)

Scaling weight 5


q(5,1) = ∂x−1 u3 ,


q(5,2) = ∂x−1 u2 q2 ,


q(5,3) = ∂x−1 u21 .

(D.4)

Scaling weight 6


q(6,3) = ∂x−1 u21 q .


q(6,2) = ∂x−1 u2 q3 ,


q(6,1) = ∂x−1 u3 q ,

(D.5)

Scaling weight 7


q(7,5) = ∂x−1 uu21 ,


q(7,4) = ∂x−1 u2 q(4)


q(7,3) = ∂x−1 u2 q4 ,


q(7,2) = ∂x−1 u3 q2 ,


q(7,1) = ∂x−1 u4 ,


q(7,6) = ∂x−1 u21 q2 ,


q(7,7) = ∂x−1 u22 .


q(8,2) = ∂x−1 u3 q3 ,


q(8,3) = ∂x−1 u3 q(3) ,

(D.6)

Scaling weight 8


q(8,1) = ∂x−1 u4 q ,


q(8,5) = ∂x−1 u2 q(5,2) , q

(8,8)

u21 q(3)

= ∂x

−1





q(8,6) = ∂x−1 u2 q(5,3) , (8,9)

,

= ∂x

−1

q

u21 q3





q(8,7) = ∂x−1 uu21 q

(8,10)

,


q(8,4) = ∂x−1 u2 q5

= ∂x

−1

q

u22 q




(D.7) (8,11)

,

q

= ∂x

−1

u31




.

Scaling weight 9


q(9,1) = ∂x−1 u5 ,


q(9,2) = ∂x−1 u4 q2 ,


q(9,5) = ∂x−1 u2 u21 , (9,8)

q

= ∂x

−1

u22 q2





q(9,3) = ∂x−1 u3 q4 ,


q(9,6) = ∂x−1 uu21 q2 , (9,9)

,

= ∂x

−1

q

u31 q





q(9,7) = ∂x−1 u21 q4 (9,10)

,


q(9,4) = ∂x−1 u2 q6

= ∂x

−1

q

u23




(D.8) (9,11)

,

q

= ∂x

−1

uu22




.

More scaling-weight 9 non-local terms exist, but are unnecessary for the analysis through O (ε 3 ). To isolate pure inelastic terms, one needs to specify the constants of integration in the integrals. The rule is that the integrals ought to have a well-defined symmetry under reflection of the argument. The guideline is the integral in the single soliton case, Eq. (1.5). For instance, the full definition of the first non-local quantity, q, which is anti-symmetric under reflection of its argument, is q (t , x) = −2k +

x

Z

u (t , x) dx = 2k tanh[k(x + v t )].

(D.9)

−∞

The generalization for the multiple-soliton case is q (t , x) =

N X

(−2ki ) +

Z

x

u (t , x) dx.

(D.10)

−∞

i =1

In Eq. (D.10), the sum runs over all the solitons in an N-soliton solution. Away from the soliton-collision region, this expression asymptotes into the sum of the single-soliton contributions: N X

q (t , x) →

(2ki tanh [ki (x + vi t + ηi )]) .

(D.11)

i =1

Appendix E. Second- & third-order coefficients in single-soliton case The coefficients in the second- and third-order corrections to the solution of Eq. (1.1) when u is a single-soliton solution are given here. The notation is as in Appendix C. The superscript S signifies that this is the single-soliton case. O(ε2 ) (S )

b1

=

75 8

−

α12 − 15α1 α2 +

25 6

α3 α4 +

15 8

65

b10 = 125α12 − 185α1 α2 +

(S )

196 3

β1 +

343 3

14

α42 −

(S )

−

α22 −

18

3 55

β2 +

b13 = −25α12 + 50α1 α2 −

3

6

β1 +

α22 −

3 196

50

35

3

α1 α3 + 49

6 170 3

β2 +

25

3

α1 α3 +

β3 − 133β4 −

α22 +

100 3

α2 α3 +

9 14

α1 α3 −

β3 − 40

3 49

3 50 3

35

18 28 3

α32 +

15 4

5

α1 α4 + α2 α4 3

7

7

7

6

2

6

β4 − β5 + β6 − β7

α32 + 15α1 α4 + 90α2 α4 − 30α42

(E.1)

β5 + 42β6 − 28β7 + 21β8 − b(9S )

α2 α3 − 25α1 α4 +

35 3

β1 −

70 3

β2 −

70 3

β3 + 35β4 + b(9S ) .

Y. Zarmi / Physica D 237 (2008) 2987–3007

3003

O(ε3 ) (S )

c1

= − − − + + + + − +

625

285

α13 +

16 395 36 305

4

α1 α32 +

α12 α2 −

40 9

α2 α32 +

209

55

α1 α22 +

12 25 54

α33 −

36 266

12 42

18

γ1 −

5 154 5

α32 α4 −

103

515

α23 +

27 75 2

24 21

16

2

α1 α4 − 173

215

α12 α3 −

α1 α2 α4 −

36

25

25

α1 α2 α3 +

18 589

α22 α4 + 309

α22 α3

18 121 12

α1 α3 α4

α1 α42 + α2 α42 + α3 α42 − α43 4 8 16   119 70 7 189 7 371 α1 − α2 − α3 + α4 β1 + − α1 + 21α2 + α3 + α4 β2 9 9 5  4 18 20   3  77 28 98 7 133 77 616 − α1 + α2 + α3 − α4 β3 + 28α1 + α2 + α3 − α4 β4 3 9 5  9 9 15  3   35 7 119 49 35 119 497 α1 − 7α2 + α3 − α4 β5 + − α1 + α2 − α3 + α4 β6 18 60 4 6 6 4    20 49 217 35 1169 7 273 α1 + α2 + α3 − α4 β7 + α1 − 28α2 + α4 β8

18

α2 α3 α4 +

325

42 5

γ11 −

(S )

c27 = −2375α13 +

18

84

γ2 +

186

5

60

217

γ3 −

γ12 +

5 12715

50

174

α12 α2 −

98

γ4 +

γ13 −

5

4450

3 560

5260

16

3

19 2

5

2 1127

γ5 +

50

147

γ6 −

5

10 63

γ7 +

5

γ8 −

21 5

357

γ9 −

10

γ10

γ14 ,

(E.2) 215

α1 α22 +

3

3560

α23 +

α12 α3 −

3

4040

100

α1 α2 α3 −

9 938

3

α22 α3

α1 α32 + α2 α32 − 100α12 α4 − 1051α1 α2 α4 − 838α22 α4 + α1 α3 α4 9 3 3 2 2 2 2 3 − 720   α2 α3 α4 + 528α3 α4 − 146α1 α4 + 2344α2 α4 −200α3 α4 + 888α4 3920 84 8666 3661 518 6363 16072 α2 − α3 + α4 β1 + − α1 + α2 − α3 + α4 β2 + 2436α1 − 9 9 5  3 3 9 5   5656 84 13314 + −1568α1 + 532α2 + α3 − α4 β3 + 1701α1 + 994α2 + 616α3 − α4 β4 9 5 5     1267 322 2107 8316 1610 α1 − α2 + α3 − α4 β5 + −707α1 + 280α2 − 1288α3 + α4 β6 + 3 9 15  5   3  917 2044 280 5544 6888 + α1 + α2 + α3 − α4 β7 + 112α1 − 1561α2 + 140α3 + α4 β8 3 3 3 5 5 2604 5208 6727 6076 34937 9114 3906 1512 2604 −

− − (S )

c38 =

5 8375

− + + + + + + (S )

4 6805

9 3104

5 50155 12

α1 α32 −

γ3 −

5 11592

γ11 −

5

3

α2 α32 −

400 9

10773

γ12 +

α12 α2 + 1975α1 α22 −

605

γ4 +

25

α33 +

5 1535

9 1775 2

α1 α3 α4 + 360α2 α3 α4 − 384α32 α4 +

3

2037 5 8631

γ1 +

4 1265 9 368 3

2037

γ10 −

α13 −

+

5 9408

γ5 +

5

γ6 −

25

5

γ7 +

5

γ8 −

5

γ9

(S ) (S ) γ13 − 504γ14 − c25 − c26 .

α23 −

4880

α12 α4 − 1827

α12 α3 +

3 669 2

13370 9

(E.3)

α1 α2 α3 −

400 3

α22 α3

α1 α2 α4 + 514α22 α4 5453

α1 α42 − α2 α42 + 200α3 α42 + 714α43  4  4081 15001 5495 4179 7553 3976 3542 − α1 + α2 + α3 − α4 β1 + α1 − 1456α2 − α3 − α4 β2 2 9 9 10  3 9 15    2408 7378 9331 3311 371 511 10059 1449α1 − α2 − α3 − α4 β3 + − α1 − α2 + α3 + α4 β4 3 9 15 2 3 3 10     1295 1106 224 896 6048 − α1 + α2 + α3 + α4 β5 + 686α1 − 385α2 + 924α3 − α4 β6 3 3 9 15  5    266 140 4032 5964 − α1 − 609α2 + α3 + α4 β7 + −266α1 + 1358α2 − 140α3 − α4 β8 3

5 3875

+

γ2 +

γ10 +

α13 −

+

c41 =

γ1 −

5 11067

γ2 −

5 7224

5 20915 12

α1 α32 −

3 4074

γ11 +

γ3 +

5 9576 5

5 3906

α12 α2 + 575α1 α22 +

65 3

α2 α32 −

200 9

α33 +

5

γ4 −

25 8064

γ12 −

5 20 9 25 2

4 

4753 5

γ5 −

28336 25

5

γ7 −

3213 5

γ8 +

(S ) (S ) (S ) γ13 + 462γ14 + 9c25 + 2c26 − c37 .

α23 − 530α12 α3 +

α12 α4 +

α1 α3 α4 + 180α2 α3 α4 − 192α32 α4 −

γ6 +

6867

249 4

1253 2

910 9

α1 α2 α3 + 100α22 α3

α1 α2 α4 + 257α22 α4

α1 α42 −

3889 4

α2 α42 + 100α3 α42 + 357α43

756 5

γ9 (E.4)

3004

Y. Zarmi / Physica D 237 (2008) 2987–3007



5999

6713

−



6 2191

9

3

+

1071 5 4788

γ1 +

966

γ2 −

5 3612

α2 +

77

3 2352

9

γ3 +

5 4788

α3 −





3619

112



8596

α4 β1 + α1 − 518α2 + α3 − α4 β2 9 15  3  1694 14653 1393 1057 13167 + α1 − 238α2 − α3 − α4 β3 + − α1 − α2 − α3 + α4 β4 9 15  2 3 3   10  3 3024 595 518 112 448 + − α1 + α2 + α3 + α4 β5 + 238α1 − 140α2 + 462α3 − α4 β6 3 3 9 15 5     448 2982 70 2016 + − α1 − 252α2 + α3 + α4 β7 + −28α1 + 679α2 − 70α3 − α4 β8

+

α1 +

1855

10  1897

5 1953 25

γ4 −

2534

4977

5

γ5 −

14168 25

γ6 +

(S )

3906 5 (S )

γ7 −

5 1134 5

γ8 +

378 5

γ9

(S )

γ10 − γ11 + γ12 − γ13 + 231γ14 − 3c25 + c26 + c37 . (E.5) 5 5 5 5 Some coefficients remain free in the second- and third-order analysis. The detailed calculation of the solution in the multiple-soliton case yields that they are always cancelled by identical terms that appear in the inelastic component. Hence, they have been set to zero: +

(S )

(S )

(S )

(S )

b9 = c25 = c26 = c37 = 0.

(E.6)

All remaining coefficients correspond to monomials that contain non-local entities, and have been set to zero. Appendix F. Special polynomials Special polynomials are differential polynomials in u, which vanish identically when u is a single-soliton solution (Eq. (1.5)). Special polynomials that contain non-local quantities (see Appendix D) will be called ‘‘non-local polynomials’’. Their asymptotic form is that of soliton–anti-soliton waves of the type discussed in Sections 4, 7, 10 and 11. Polynomials that do not contain non-local terms will be called ‘‘local polynomials’’. They are localized in the vicinity of the soliton collision region, and fall off exponentially in all directions in the x–t plane. The lowest scaling weight, for which a special polynomial exists, is Ω = 3. All special polynomials with 3 ≤ Ω ≤ 5 are non-local. For Ω ≥ 6, there are non-local as well as local polynomials. Scaling weight 3 R(3,1) [u] = ux + qu,

R(3,2) = 3qu + q3 − 3q(3) .

(F.1)

Scaling weight 4 1 1 R(4,1) = − q4 + qq(3) + qu1 , R(4,2) = q(4) + u2 . 3 2 Other independent weight-4 polynomials can be constructed from the weight-3 ones by R(4,3) = ∂x R(3,1) , Replacing R

(3,1)

R(4,4) = qR(3,1) .

by R

(3,2)

(F.2)

(F.3)

in Eq. (F.3) yields linear combinations of R(4,k) , 1 ≤ k ≤ 4.

Scaling weight 5 23 2 8 4 qu − qu2 R(5,1) = q(5,1) + q(5,2) + q(5,3) − q2 q(3) − 5 15 15 R(5,2) = q(5,1) − q(5,2) − q(5,3) − qu2 R(5,3) = q(5,2) − q(5,3)

(F.4)

4 1 R(5,4) = q(3) u − qu2 − qu2 . 3 3 Other weight-5 polynomials can be constructed from lower-weight, similar to Eq. (F.3). Scaling weight 6 R(6,1) = u3 + uu2 − u21 R(6,2) = q4 u − 7u3 − 11uu2 − u4 R(6,3) = q(6,1) + q(6,2) + q(6,3) + 2u3 + uu2 1 R(6,4) = q(6,1) + u3 3

(F.5)

R(6,5) = q(6,2) − q(6,3) R(6,6) = R(3,1)


2

.

Additional weight-6 six polynomials can be constructed from lower-weight ones. The number of independent special polynomials grows rapidly with scaling weight. Many of them emerge in the results of the secondand third-order calculations. In the following, only the local polynomials are shown, as they emerge as driving terms in the analysis through O (ε 3 ). Scaling weight 7 R(7,1) = ∂x R(6,1) .

(F.6)

Y. Zarmi / Physica D 237 (2008) 2987–3007

3005

Scaling weight 8 R(8,1) = uR(6,1) R(8,2) = u4 + 2uu21 + u1 u3 − u22 R

(8,3)

4

= 3u −

9uu21

− uu4 +

u22

(F.7)

.

Scaling weight 9 R(9,1) = u1 R(6,1) R(9,2) = ∂x uR(6,1)


(F.8)

R(9,3) = 30u3 u1 − 21u31 + u2 u3 − uu5 R(9,4) = 6u31 − u2 u3 + u1 u4 . ( 2)

( 3)

˜ in and u˜ in from special polynomials Appendix G. Construction of u (2)

(3)

The polynomial parts of the second- and third-order inelastic contributions, u˜ in and u˜ in , have the structure shown in Appendix C, with coefficients denoted by bin,k and cin,k , respectively. All bin,k , 14 ≤ k ≤ 23, and all cin,k , 42 ≤ k ≤ 74, have been set to zero. The freedom in the expansion allows this choice. Their inclusion turns out to be unnecessary and only complicates the analysis. (2) (3) For u˜ in and u˜ in to be composed out of special polynomials, their coefficients have to satisfy: 2 O(ε ) 1 4bin,1 − bin,8 − bin,9 + bin,13 = 0 2 11 5 1 1 1 5 37 − bin,1 + bin,2 − bin,3 + bin,4 − bin,5 − bin,6 + bin,7 + bin,8 + bin,9 3 3 3 3 3 2 3 2 2 2 + bin,10 − bin,11 + bin,12 − bin,13 = 0 (G.1) 3 3 3 3 5 3 17bin,1 − 4bin,2 + bin,3 + bin,5 − bin,7 − bin,8 − bin,9 − bin,10 + bin,11 − bin,12 + bin,13 = 0 4 2 2 20 7 4 1 1 1 5 1 1 1 1 bin,1 − bin,2 + bin,3 − bin,4 + bin,6 − bin,8 − bin,9 − bin,10 + bin,11 − bin,12 + bin,13 = 0. 9 9 9 3 3 12 18 9 9 9 6 O(ε3 ) 1 1 1 −180cin,1 + 34cin,2 − 4cin,3 − 4cin,5 + cin,7 + cin,8 − cin,11 + cin,12 + cin,17 − cin,18 2 2 2 1 1 + cin,22 + cin,23 − cin,24 + 11cin,25 + 4cin,26 + 2cin,27 − 4cin,28 − cin,29 + cin,30 + cin,33 2 2 − cin,35 − cin,37 − cin,38 + cin,39 − cin,40 = 0, (G.2) 2835

45 79 45 31 45 cin,3 + cin,4 + cin,5 − cin,6 − cin,7 − cin,8 + cin,9 8 8 8 8 8 8 15 5 3 3 45 15 5 + cin,10 + cin,11 − cin,12 + cin,13 + cin,14 − cin,15 − cin,16 − cin,17 + cin,18 8 4 8 8 8 8 4 3 3 11 − cin,19 − cin,20 − cin,22 − cin,23 + cin,24 − 22cin,25 − 8cin,26 − 4cin,27 + 8cin,28 8 8 8 + 2cin,29 − 2cin,30 − 2cin,33 + 2cin,35 + 2cin,37 + 2cin,38 − 2cin,39 + 2cin,40 = 0, 8

cin,1 −

951

8

cin,2 +

941

19

341

47

3

23

3

1

4

(G.3)

15

cin,2 − cin,3 + cin,4 − cin,5 − cin,6 + cin,7 + cin,8 + cin,9 12 12 4 4 4 4 3 4 3 1 1 5 5 15 3 1 1 − cin,10 − cin,11 − cin,12 − cin,13 − cin,14 − cin,15 + cin,16 + cin,17 + cin,18 4 3 6 12 12 4 4 3 6 5 5 1 1 1 28 14 10 14 + cin,19 + cin,20 + cin,22 + cin,23 + cin,24 + cin,25 + cin,26 + cin,27 − cin,28 12 12 3 3 4 3 3 3 3 8 2 2 2 2 8 − cin,29 + cin,30 + 2cin,31 − 2cin,32 + 2cin,34 − cin,37 − cin,38 + cin,39 − cin,40 = 0, 3 3 3 3 3 3 4 4 4 136cin,1 − 2cin,8 + cin,21 + cin,22 + cin,23 − 16cin,25 − 4cin,27 + 2cin,36 + 4cin,37 − 4cin,41 = 0 3 3 3 1 1 1 1 1 1 1 1 1 cin,1 − cin,2 + cin,3 − cin,4 − cin,5 + cin,6 + cin,7 − cin,9 − cin,10 8 8 8 8 24 8 24 8 24 1 1 1 1 1 1 1 1 1 − cin,12 − cin,13 − cin,14 + cin,15 + cin,16 + cin,18 + cin,19 + cin,20 + cin,24 = 0. 60 120 120 8 24 60 120 120 72

−

cin,1 +

499

(G.4) (G.5)

(G.6)

3006

Y. Zarmi / Physica D 237 (2008) 2987–3007

(2)

˜ in Appendix H. Coefficients used in u

(2)

Exploiting the freedom inherent in solving Eq. (11.1), the coefficients in u˜ in , the polynomial part of the second-order inelastic term, have been chosen so as to simplify the results, in particular, the structure of the obstacle to asymptotic integrability. The following bin,k are non-zero: 5

bin,1 =

24α1 α2 − 8α22 + 6α1 α3 − 30α2 α3 + 8α32 − 21α1 α4 + 24α2 α4 + 10α3 α4 − 13α42

9

− bin,2 = bin,3 =

25

(β1 + 2β2 − 6β3 + 3β4 + 4β5 − 6β6 + 4β7 − 2β8 ) ,

(α2 − a4 ) (3α1 − 2α3 + α4 ) −

9 50

bin,5 =

7 3




14 3

(H.1)

(β4 − β8 ) ,

(H.2)

(α2 − a4 )2 ,

9 5

(H.3)

9α1 α2 + 2α22 + 6α1 α3 − 20α2 α3 + 8α32 − 6α1 α4 − α2 α4 + 2α42

3




− 7 (β1 + 2β2 − 6β3 + β4 + 4β5 − 6β6 + 4β7 ) , bin,9 =

5

135α12 − 39α1 α2 − 12α22 − 6α1 α3 − 190α2 α3 + 72α32 − 129α1 α4 + 221α2 α4 + 40α3 α4 − 92α4

9

−

14

bin,10 = −

+ bin,11 = bin,13 =

(H.4)


2

3 25

9 14 3

100 3 25 9

−

(12β1 − 6β2 − 32β3 + 27β4 + 18β5 − 27β6 + 18β7 − 10β8 ) ,

(H.5)

27α12 − 66α1 α2 + 12α22 − 18α1 α3 + 32α2 α3 − 8α32 + 24α1 α4 − 6α3 α4 + 3α42

(5β1 − 20β2 + 10β3 + 7β4 − 10β5 + 15β6 − 10β7 + 3β8 ) ,

(H.6)

α1 (α2 − α4 ) − 28 (β4 − β8 ) ,

(H.7)

27α12 − 27α1 α2 + 4α22 − 6α1 α3 − 14α2 α3 + 8α32 − 9α1 α4 + 25α2 α4 − 8α42

14 3




(10β1 − 10β2 − 20β3 + 21β4 + 10β5 − 15β6 + 10β7 − 6β8 ) .


(H.8)

Appendix I. Lesson regarding Fredholm Alternative Theorem (1)

The fact that Eq. (5.3) is solved by uin , for which a closed-form expression exists (Eq. (4.1)), provides an important lesson with regard to the Fredholm Alternative Theorem. The driving term in that equation does not satisfy the conditions of the Theorem. Namely, it is not orthogonal to solutions of the equation adjoint to the homogeneous part of Eq. (5.3). It is ‘‘common wisdom’’ that in a situation like this, (1) the solution of the equation must contain an unbounded (secular) contribution. However, uin is clearly bounded throughout the x–t plane. In the analysis from second order onwards, the dynamical equations will also contain driving terms that are localized in the x–t plane. For some of these driving terms, there will be (bounded) closed-form solutions in terms of differential polynomials in u. For others (the ‘‘obstacles to asymptotic integrability’’), there will be no such closed-form solutions. Still, as the driving terms are localized in the plane, just like R(7,1) , the resulting solutions will be bounded. To see why this is possible, let us examine R(7,1) in detail. It is bounded and appreciable only in the soliton-collision region (a finite neighborhood of the origin) and decays exponentially fast away from that region (so does R(6,1) ). Hence, it does not resonate with the KdV solitons in u, and cannot generate unbounded contributions along their characteristic lines. For instance, in the two-soliton case, it decays along the line of soliton no. 1 as: R1 ∝ e

  −2 k2 x+4k22 t


k1 x + 4k2 t → C , 1


k2 x + 4k2 t → ∞. 2

(I.1)

The only potential pitfall is the possibility that R(7,1) may generate secular terms away from the soliton characteristic lines, namely, in regions in the x–t plane where u itself falls off exponentially. These are the triangular sectors {|4k2(i−1) t |  |x|  |4k2i t |} in the plane, bounded by the characteristic lines of adjacent solitons. In these regions, Eq. (5.3) is reduced to:

∂t u(in1) = ∂x3 u(in1) + 10 (α2 − α4 ) R(7,1) [u] .

(I.2) (7,1)

To resonate with the homogeneous part of Eq. (I.2), R must fall off as exp[−k|(x + k t )|] for some k. If it falls off at such a rate, then (1) uin (t , x) = (x + at ) exp[−k|(x + k2 t )|] solves Eq. (I.2). Unless a = k2 , this solution is unbounded in parts of the x–t plane. However, R(7,1) cannot generate such solutions. For example, in the two-soliton case, it falls off in the triangular sectors as: R1 ∝ e

    −2 k1 x+4k21 t −2 k2 x+4k22 t


k1 x + 4k2 t → ∞, 1

2


k2 x + 4k2 t → ∞. 2

(I.3)

Similar localized driving terms are encountered in the analysis of inelastic contributions in higher orders. For the arguments proffered here, they are expected to generate bounded inelastic terms.

Y. Zarmi / Physica D 237 (2008) 2987–3007

3007

References [1] [2] [3] [4] [5] [6] [7] [8] [9] [10] [11] [12] [13] [14] [15] [16] [17] [18] [19] [20] [21] [22] [23] [24] [25] [26] [27] [28] [29]

A.S. Fokas, Q.M. Liu, Phys. Rev. Lett. 77 (1996) 2347–2351. Y. Kodama, Phys. Lett. 112A (1985) 193–196. Y. Kodama, Physica 16D (1985) 14–26. Y. Kodama, Normal form and solitons, in: M.J. Ablowitz, et al. (Eds.), Topics in Soliton Theory and Exactly Solvable Nonlinear Equation, World Scientific, Singapore, 1987, pp. 319–340. Y. Hiraoka, Y. Kodama, Normal form and solitons, in: Lecture Notes, Euro Summer School 2001, The Isaac Newton Institute, Cambridge, 2002, p. 15–25. T.R. Marchant, Chaos Solitons Fractals 22 (2004) 261–270. D.J. Korteweg, G. De Vries, Philos. Mag. 39 (1895) 422–443. H. Washimi, T. Taniuti, Phys. Rev. Lett. 17 (1966) 996–998. V.I. Karpman, Non-Linear Waves in Dispersive Media, Pergamon, Oxford, 1975. N.J. Zabusky, M.D. Kruskal, Zarmi, Phys. Rev. Lett. 19 (1967) 1095–1098. R. Hirota, Phys. Rev. Lett. 27 (1971) 1192–1194. M.J. Ablowitz, D.J. Kaup, A.C. Newell, H. Segur, Stud. Appl. Math. 53 (1974) 249–315. V.E. Zakharov, S.V. Manakov, Sov. Phys. JETP 44 (1976) 106–112. M.J. Ablowitz, H. Segur, Solitons and the Inverse Scattering Transforms, SIAM, Philadelphia, 1981. S.P. Novikov, S.V. Manakov, L.P. Pitaevskii, V.E. Zakharov, Theory of Solitons, Consultant Bureau, New York, 1984. A.C. Newell, Solitons in Mathematics and Physics, SIAM, Philadelphia, PA, 1985. R.A. Kraenkel, Phys. Rev. E 57 (1998) 4775–4777. T.R. Marchant, N.F. Smyth, IMA J. Appl. Math. 56 (1996) 157–176. T. Fokas, L. Luo, Contemp. Math. 200 (1996) 85–98. R.A. Kraenkel, J.G. Pereira, E.C. de Rey Neto, Phys. Rev. E 58 (1998) 2526. T.R. Marchant, Chaos Solitons Fractals 32 (2007) 1328–1336. R.M. Mirie, C.H. Su, J. Fluid Mech. 115 (1982) 475–492. Q. Zou, C.H. Su, Phys. Fluids 29 (1986) 2113–2123. A. Veksler, Y. Zarmi, Physica D 217 (2006) 77–87. A. Veksler, Y. Zarmi, Nonlinearity 20 (2007) 523–536. C.S. Gardner, J.M. Greene, M.D. Kruskal, M. Miura, Comm. Pure Appl. Math. 27 (1974) 97. P.J. Olver, Applications of Lie Groups to Differential Equations, Springer-Verlag, New York, 1986. M.J. Ablowitz, P.A. Clarkson, Solitons, Nonlinear Evolution Equations and Inverse Scattering, Cambridge University Press, 1991. Y. Li, D.H. Sattinger, J. Math. Fluid Mech. 1 (1999) 117–130.