Chaotic and homoclinic behavior for numerical discretizations of the nonlinear Schrödinger equation

PHYSICA

Physica D 57 (1992) 447-465 North-Holland

Chaotic and homoclinic behavior for numerical discretizations of the nonlinear Schr6dinger equation David W. McLaughlin Department of Mathematics and Program in Applied and Computational Mathematics, Princeton University, Princeton, NJ 08540, USA

Constance M. Schober Program in Applied Mathematics, Unil,ersity of Arizona, Tucson, A Z 85721, USA

Received 2 July 1990 Revised manuscript received 6 November 1991 Accepted 12 December 1991 Communicated by H. Flaschka

Certain conservative discretizations of the NLS can produce irregular behavior. We consider the diagonal discretization as a conservative perturbation of the integrahle discretization and study the homoclinic crossings in its nonlinear spectrum. We find that irregularity sets in when two homoclinic structures are present and, in this case, many and continual homoclinic crossings occur throughout the irregular time series. We indicate a Melnikov analysis to study the consequences of this homoclinic behavior.

1. Introduction T h e m o d u l a t i o n a l instability is an e x t r e m e l y i m p o r t a n t instability in n o n l i n e a r d i s p e r s i v e waves. F o r e x a m p l e , it is r e s p o n s i b l e for e n v e l o p e solitons a n d cavitons in p l a s m a s , for f i l a m e n t a t i o n o f laser b e a m s , a n d for t h e b r e a k up o f m o n o c h r o m a t i c w a v e t r a i n s in water. This instability p r o d u c e s c h a o t i c waves w h e n a conservative system is p e r t u r b e d by d a m p i n g a n d driving [1, 2]. R e c e n t l y it has b e e n shown to p r o d u c e striking effects in n u m e r i c a l a p p r o x i m a t i o n s o f c o n s e r v a tive waves [3]. A s a m p l e o f such n u m e r i c a l effects is r e p r o d u c e d in fig. 1 for o n e p a r t i c u l a r disc r e t i z a t i o n of the n o n l i n e a r S c h r 6 d i n g e r e q u a tion.

T h e m o s t basic m o d e l u s e d to study p h e n o m e n a r e l a t e d to the m o d u l a t i o n a l instability is the self-focusing n o n l i n e a r S c h r 6 d i n g e r e q u a t i o n (NLS):

i g , , = -~Ox. - 21~0/2~

(1)

with p e r i o d i c b o u n d a r y c o n d i t i o n s

6(x+C,t) = 6(x,t).

(2)

T h e x - i n d e p e n d e n t p l a n e wave solution ~9= a e 2ilap2t is k n o w n to be linearly u n s t a b l e . W h e n e x a m i n i n g p e r t u r b a t i o n s o f this solution of the

0167-2789/92/$05.00 © 1992-Elsevier Science Publishers B.V. All rights reserved

448

D.W. McLaughlin, C.M. Schober / Discretizations of the NLS (3) Understand the structure of level sets and orbits in a full phase space neighborhood of the critical level sets on which these stable and unstable manifolds reside. (4) Numerically detect all such instabilities precisely and continuously monitor them in time even in the presence of perturbations.

Fig. 1. Sample time series of one nonintegrable discretization of the nonlinear Schr6dinger equation, 0
# s ( x , t ) = [a + ec~(x,t)] e 2ilal2t,

(3)

subject to the initial condition 0 ( x , 0 ) = a + ca(x,0), it is found that for a ( x , t ) = e ' r ' e i ~ * 6 , the linear growth rate 0- is given by 0-2= --K2(K2--41a12). The long waves are unstable provided 0
(4)

Furthermore, the number of unstable modes is given by the largest integer N such that 0 < N <

lalg/v. Such instabilities for special plane-wave solutions could be studied for a wide class of partial differential equations, not just this integrable NLS equation. However, since NLS is integrable via soliton mathematics, one can understand much more about its instabilities and their consequences. In particular, for this integrable equation one can [10, 11, 18] (1) Identify all instabilities of a general (rather than just a plane wave) solution, and calculate the dimensions of its stable and unstable manifolds. (2) Follow (exponentiate) each local linear instability into a global homoelinic orbit.

More geometrically, the flow for the integrable system takes place on level sets of constants of motion. Typically, these level sets are tori, whose dimension equals one half the dimension of the phase space. (In the PDE setting, the phase space is a function space which is infinite dimensional, with infinite dimensional tori. However, each numerical approximation is finite ( 2 N ) dimensional. Of course, N increases with the inverse of the discretization scale h, N ~ h-~ and can be very large.) Denoting the phase space by S, the integrable flow resides on the level sets of N commuting constants

{Hi" S ~ ~IH, = h i , j = 1 ..... N } .

(5)

For almost all values of the real numbers h i, these level sets are N-tori T N. However, there can exist certain special "critical values" h/*, and associated "critical level sets" H*, with the following properties: (1) The invariant tori T '~1 in a critical level set no longer have maximal dimension; rather, their dimension is M which can take any integer value less than N. (2) At T M c H *, the constants {Hi: S ~ g{, j = 1,2 . . . . . N} no longer have linearly independent gradients. (3) The critical values he/ are not local maxima or minima; rather they are saddles for the constants Hi. (4) Locally, at least, the critical level set H* has the product structure T M x ~(N-M)

D.W. McLaughlin, C.M. Schober / Discretizations of the NLS

together with various " b o u n d a r y c o m p o n e n t s " of the form TMx~

(N J),

J=M+

1,...,N.

(5) The boundary component of smallest dimension is the invariant torus T ~t itself. All orbits in the critical level set approach the torus TM in both forward and backward time. In this note, we define the term "homoclinic crossing" as follows: Consider a trajectory u(t) in phase space which is evolving in time t according to a dynamical system which is a perturbation of an integrable system. When the constants of motion {H i} for the integrable system, evaluated on this perturbed trajectory, cross a critical value h*, we say that a homoclinic crossing has occurred. As described below, for soliton systems these homoclinic crossings can be precisely detected numerically by the spectral transform. Thus integrability via soliton mathematics provides a theoretical and numerical tool with which to analyze the presence of these instabilities as a step toward inferring their chaotic consequences. Previously, such studies were pursued for the damped driven sine-Gordon and NLS equations [6, 7]. Here we apply this program to the interesting numerical p h e n o m e n a in conservative discretizations of the NLS equation identified by Ablowitz and Herbst [3]. We emphasize that a thorough numerical study of chaotic behavior in conservative discretizations of NLS, for different discretizations with wide ranging temporal durations, p a r a m e t e r values, and initial conditions, is beyond the scope of this note. Such studies, as initiated by Ablowitz and Herbst, are certainly needed, particularly here in this delicate situation of chaotic behavior in high dimensional conservative systems. Our purpose here is merely to take sample numerical experiments, at p a r a m e t e r values suggested in ref. [3], and show that our (very different) diagnostic tools from integrable soliton mathematics can shed light on these particular experimental results. Thus, these integrable diag-

449

nostics should be included in subsequent, more extensive and detailed, numerical experiments. Specifically in this note, we study the natural diagonal discretization ( D D N L S ) iq,, = - q"+J - 2q~ + % _ 1 h2 - 21q. 12q,,

(6)

as a perturbation of the integrable discretization (IDNLS) [14] ion = -

q,, + j - 2 qn + q,, - I _ ]q.[Z(q.,+, + q ~ _ , ) . h2

(7) with h = L / N the mesh size and under periodic boundary conditions qj+N = qj. It is clear that eq. (6) can be viewed as an ef(h 2) perturbation of the integrable discretization (7), i.e. (6) can be written as ion + (q,,- l - 2 q , , +

qn+~)/h 2

+ [q,,12(q,,_l +q,,+t) = Iq~12(q~ , - 2 q ~ + q n + , ) ,

(8)

which is (7) with a perturbation of order h21ql2q**. The latter discretization (7) is integrable because it is equivalent to consistency of the following linear problems:

(9)

Vn+ , =L,,(z,qn(t))V,,,

(10)

I),= T~(z,q,(t))Vn, where

L.

,

V/i+ Iq.l"

(z

'

(11)

Tn=i ( l-z2-1(qnq*n-l-bqt*~qnl) z 1 - q n* z + q n - l / z " qn I --q*n/Z --l+z-Z+~(qrtqn t+q,Tqn t)

(12) (Here h has been set equal to 1.) Thus (7) can be integrated with the discrete inverse spectral

450

D. 14d McLaughlin, C.M. Schober /Discretizations of the NLS

transform associated with the Lax pair [ L , , T , ] . We developed this Lax pair, which is slightly different from an equivalent pair in the literature [14], specifically to be able to handle the spectral theory associated with periodic potentials [13] (see ref. [19]). Quite generally we are interested in understanding the role played by homoclinic orbits in integrable systems as causes, sources, and organizers of chaotic behavior in nearby systems. For soliton systems we have previously developed theoretical [10, 11] methods for the detection of homoclinic crossings and have implemented these methods numerically [1, 2] in damped-driven situations. The problem studied here differs because the nearby systems remain conservative; moreover, the perturbed systems themselves are of independent interest for numerical analysis. Specifically, certain difference schemes can produce irregular time series, a phenomenon which becomes important when different schemes are compared at the same fixed discretization size. We apply the numerical tools of refs. [1, 2] to these conservative difference schemes and show with our measurements that the irregular time series do indeed contain many and continual "homoclinic crossings". In addition, we find that irregular behavior in these conservative time series is associated with "multiple homoclinic structures" (M > 1) which are more complicated than single homoclinic orbits ( M = 1). When only a single homoclinic orbit is present, the time series tested behave regularly (although admittedly, we have only performed a limited number of experiments).

section we present the relevant features of the Floquet theory for the difference operator (9) which we use to analyze the numerical experiments. Define two solutions O(n, z), O(n, z) of the linear problem (9) by the initial conditions

Then the fundamental matrix for (9) is given by :

=~-~,_I(Z,t)[~(1,Z) ~(1,z)], n = 2 , . . . , N + 1,

(14)

where

J,,(z,t)

=

HLj(z,<(t))

(is)

j=l

3 N ( Z , t) is the transfer matrix which maps solu-

tions of (9) across one period N of the potentials q,(t). In general, any solution of (9) will have the

form

U,,+kN:K(n,z)(JN(Z,t))kc,

1 <_n<_N, (16)

where C is a constant vector. The behavior of the solution is then determined by the behavior of ( ~ v ) k. To ensure stability of the solution as it is mapped across the lattice, the eigenvalues of JN(Z, t) must lie within or on the unit disk. The spectrum of this operator is then defined to be o(,~N(z,t))

2. Spectral theory background = {z c C I J N ( z , t ) is real, - 2 _ ~ N ( z , t )

Discretization (7) can be integrated with inverse spectral theory, actually with Floquet theory since the coefficients of the linear difference operator in (9) are periodic of period N. The homoclinic crossings are detected by numerically implementing the inverse spectral theory. In this

_~ 2}, (17)

where AN(Z, t), the Floquet discriminant, is given by A N ( Z , f) = T r g - N ( Z , t ) .

(lS)

D. ldd McLaughlin, C M. Schober / Discretizations of the NLS

gDN(Z,t) is analytic in z in the finite complex plane \ {0} and

N-I A N ( Z , t ) = z -N z 2 N + •

k=l

) a, z2k+l

,

(19)

where a k = ak(q I. . . . . qN, q~ . . . . . q~v) is a polynomial in its arguments. The integrability of IDNLS (7) follows from the invariance of A N ( Z , t), which actually contains all N constants of motion for this integrable Hamiltonian system. A few additional spectral quantities will be needed: the "simple periodic spectrum" is o - " = { Z , I A N ( Z k , t ) = +2, A'N(Z,,t ):gO},

(20)

and the "double points" of the spectrum are

~r~ = { z , [ A x ( z , , t ) = +_2, A'u(Zk,t ) = 0 } ,

(21)

where the prime denotes the derivative with respect to z. The spectral plots presented later, display AN(A , t) instead of AN(Z, t) where z = e + ia. This is done to enable a convenient comparison with the spectrum of the NLS given by (1). All future references to the spectrum of the D N L S will be with respect to the complex variable A. Since o - ( 3 N) has the following symmetry A ~ r ( g f ) then A* ~ o'(=~-X), only the upper half A plane is displayed in the spectral plots.

Remark 1 (To interpret the spectral plots.) The entire real axis is always the continuous spectrum for ~Y-N for all N and no special marking is indicated. The discrete or periodic spectrum is marked as follows: [] (~ Q)

corresponds t o

AN(A , t) = + 2, A N ( A , / ) = --2, AN(A,/) = 0.

T h e r e is no notational distinction between whether an element of the discrete spectrum is a simple point or a double point, so care must be

451

exercised. Analytically one can show that simple points are always the end points of a spine of spectrum and double points must lie on a spine of spectrum. The large + and × correspond to a continuous spectrum where the discriminant A N is respectively positive or negative. The small + and × correspond to curves of real A N which are not spectrum and which are respectively positive or negative. The spectrum of the initial plane wave consists of the real axis and a single band of the spectrum along the imaginary axis in the upper half plane. Complex double points can reside on this band of spectrum. The significance of these complex double points is that there is a one-to-one correspondence between the number of complex double points and the number of exponential instabilities of the linearized IDNLS. Further, the tori in the level sets for even (about x = 0) potentials lose one dimension for each of these complex double points. In the case of two purely imaginary complex double points, these tori have dimension N - 2 , with stable (unstable) manifolds of the torus T N-2 in this level set having dimension N, two more than the dimension of the torus itself. Coordinates for these stable (unstable) manifolds can be constructed with Bficklund transformations which use the complex double points as parameters. This construction is given in continuous P D E cases in refs. [11, 18] and in the discrete ( M = 1) NLS case in ref. [19]. Under small, even perturbations q(x), purely imaginary double points open either into gaps in the spectrum (gap state) or form spines of spectrum (cross state). In this even, purely imaginary, situation, continuity forces spectral configurations which pass from a gap state to a cross state to pass through a double point and thereby cross the associated homoclinic structure. So by examining the changes in the spectrum under the perturbed flow we will be able to detect the homoclinic crossings. It is found upon examining an e-neighborhood of the homoclinic orbit, that the distinct (even) spatial structures that exist nearby correspond to

D.W. McLaughlin, C.M. Schober / Discretizations of the NLS

452

+

+

+

~

+

+

~

+

-04

-0.2

0

+

t

0.2

0.4

DI SC

(a)

(b)

1.0 0.5 /../"

0.0 "....

.

-0.5 -l.O -0.5

0.0

0.5 A

1.0

1.5

(c) Fig. 2. (a) The homoclinic solution of the I D N L S for 0 _< t _< 50 with initial condition q(x,O) = ~, + e0(l + / ) c o s px. Thc nonlinear spectrum with one purely imaginary double point. (b) The surface q(x, t). (c) The phase portrait.

distinct spectral configurations near to the complex double point. This association of spectral configurations with spatial structures has been discussed for the sine-Gordon [10, 11] and the NLS [12]. In the next three examples, for the I D N L S (7), we establish a similar correspondence using initial conditions that are prescribed in ref. [3]. Fig. 2a is the surface plot for the homoclinic orbit associated with one complex double point.

The initial condition is q,, = 0.5 + e0(1

x,, = - L / 2

+ i)

cos px,,,

+ (n - 1)L/N,

(22) n = 1,2 . . . . . N +

1,

(23) where e 0 = 10 -~, N = 32, p = 2 w / L , L = 2~/2w. The homoclinic orbit is characterized by a single m o d e which is centered in the middle and the

D.W. McLaughlin,C.M.Schober/ Discretizationsof the NLS asymptotic decay of the homoclinic orbit to the plane wave is clearly visible. Fig. 2b is the corresponding spectral plot which shows the one purely imaginary double point. Fig. 2c is the phase portrait where we see the right lobe of a separatrix. In the phase portraits we are considering the phase space is formed by the following two variables:

A ( t ) = Iq(0, t)]2 _ 0.25, At=

dA dt'

(24) (25)

where A(t) represents the departure of the norm of the amplitude from equilibrium. There exists a second homoclinic orbit corresponding to the double point shown in fig. 2b. It is obtained by a space translation and the spatial excitation is centered in the wings and its phase orbit is the left lobe of an asymmetric separatrix. In an e-neighborhood of this homoclinic orbit there are two distinct spatial structures which are specified by the initial conditions q,, = 0.5 + 0.1 cos px n

(26)

and qn = 0.5 + 0.1icos px,,

with its connection to the connectivity of the level sets, may be found in ref. [18].

3. Numerical diagnostics We begin by integrating discretizations (6) and (7) numerically using the R u n g e - K u t t a routine in the IMSL library with a relative error of 10 -6. This relative error was checked by successively halving the tolerance level of the integrator routine and by comparison with known analytic solutions. The discriminant AN(A,t) , and in particular, the periodic spectrum, is invariant in time when qn(t) evolves by IDNLS (7). In the Hamiltonian framework for periodic solutions the periodic spectrum provides the actions (one half of the degrees of freedom of the system) in an action-angle description and they are constants of motion. These N constants of motion provide an excellent check in the time integrator and when they remain constant, at least one half of the variables are being integrated accurately. For the DDNLS (6) there is much less information to check as there are only two constants of motion, the L 2 norm

(27)

N-1

I= where L, p, x n and N are the same as before. Their spectral and surface plots and phase orbit are shown in figs. 3 a - 3 c and figs. 4a-4c, respectively. Figs. 3a-3c show the correlation between a "gap" spectral configuration, the spatial excitation locked in the center, and the phase orbit trapped within the right lobe of the separatrix. Figs. 4a-4c show the correlation between a "cross" spectral configuration, the waveform traveling between center and wings and the phase orbit lying outside the separatrix. Moon [16] illustrates that for the d a m p e d - d r i v e n NLS, frequent homoclinic crossings correspond to a random flipping between these waveform types. A bifurcation description of the phenomena, together

453

Y'. Iqj[ 2

(28)

j - 0

and the Hamiltonian

N-1 ( Iqj+ Oj12 H=-iY'~ l--

Iq l

) .

(29)

j=o

We are considering periodic even solutions and we want this feature to be preserved in our numerical solutions. This spatial symmetry allows us to work on half the x-interval, O < x < L / 2 . By imposing this symmetry condition at x = L / 2 we prevent the loss of symmetry as time evolves. In our numerical experiments, time series of H and

D.W. McLaughlin, C.M. Schober /Discretizations of" the NLS"

454

0 4

0.2

0 •

0.2

0.4

Ol S ( "

(a)

(b)

20

1.5

1.0 .-....,.

./:v. ..,...'...... . .."i •"_•.. '.•...

...,~: •s

o

E--,

"...•

g-.':

~"

k -.5

t%

/

>: j . , •...•.

-1.0

"....,,. ,c? : -...>

2/ ..9 "

-1.5

-~o

,

.6

, ~ - ~ - ~

-.4

-.2

•

0

.2

.

.

.4

.

.

.5

.

.

.

.8

.

'

1.0

1•2

14

'

~

1.6

-

1,B

2.0

A (c) Fig. 3. (a) The solution of the IDNLS for 0 _< t _< 50 with initial condition q(x,O) = 2 + 0.1 cos px. The nonlinear spectrum with one imaginary gap. (b) The surface q(x, t) showing an excitation locked in the center• Also there exists a second case, which is not shown, with the excitation locked in the wings. (c) The phase portrait•

I remained constant to six significant digits. Once the time series is generated, we analyze it with both diagnostics from dynamical systems theory and from inverse spectral theory. The dynamical systems diagnostics presented here are the phase plane z(t) = (z(t), zt(t)) for z ( t ) = I q ( 0 , t)[ 2 0.25 and the temporal power spectra A~ where A = Iq(0, t)l. In our numerical experiments, we consider small values of the p a r a m e t e r 6 = h2/12, so we

are in a near-integrable situation. It seems reasonable then to apply the nonlinear spectral transform and we consider this diagnostic even more illuminating. At each time t, we take the n u m e r i c a l l y g e n e r a t e d solution {q,,(t)ln = 1 . . . . . N} and perform the spectral transform by numerically determining ~N(,~). For solutions generated by the I D N L S (7), the spectrum is constant. But for solutions generated by D D N L S (6), the perturbation causes the spectrum to

D.W. McLaughlin, C.NL Schober / Discretizations of the NLS

+

+

+ ,-

0,4

/

455

,

~,

-0.2

P

'~

,

'~

02

2

4

+ I

,~,

0.4

DI SC

(a)

(b)

1.5 ..•'

1.11

.:..'.

,:,

,',

..:2"

•5

.,,',

0

n:l

~y,e

<

-1.0

i.

,

-1.5

-2.0

1_141 .b

.

12L .

I fJ

i

2L

i

.

I .4

i

I .b

i

18 i .

J 1,0

i

I 1.2

i

i 1.4

I 1.6

i

I 1.8

i 2.B

A (c) Fig. 4. (a) The solution of the I D N L S for 0 _< t _< 50 with initial condition q(x, 0) = ~1 + 0.1i cos px. The nonlinear spectrum with one complex cross. (b) The surface q(x, t) with a standing wave whose maximum travels between the center and the wings. (c) The phase portrait.

change and so the spectrum must be measured successively in time. The spectrum ON(A ) p r o v i d e s a measurement of the components of the waveform at time t. Further, by observing when the spectral configuration passes through a double point we are able to detect the homoclinic crossings. We are thus able to determine (as a function of time) the content of the waveform

and the homoclinic crossings for the perturbed system. The spectral plots are crucial to the identification of the homoclinic crossings and these were obtained using Overman's spectral solver code as was done in refs. [1, 2]. We have used a slight adaptation for the discrete case but a general description is provided in ref. [6]. Without the

456

D.W. McLaughlin, C.M. Schober /Discretizations of the NLS

spectral solver machinery, it is only known that the double points (and instabilities) are present initially. The numerics presented here show that the crossings persist throughout the time series and can be correlated with the onset of irregular behavior. This detection of homoclinic crossings is precise because we are in a near-integrable scenario and we have the nonlinear spectral transform to provide an accurate measurement of the spectrum. For the DNLS, this is even simpler than for NLS as we are dealing with difference operators and AN(A, t) = qSI(N + l, A) + &2(N + 1, A) is determined by direct matrix multiplication.

4. Homoclinic crossings We consider the diagonal DNLS as a perturbation of the integrable DNLS. This allows us to apply the previous spectral theory even though the spectrum is no longer constant under the perturbed flow. Upon examining the spectrum of the perturbed flow we obtain the homoclinic crossings as a function of time for a variety of initial spectral configurations (variations on the one or two double point scenario) and mesh sizes (perturbation strengths) and correlate these findings with the physical features of the solutions in an effort to understand the effects of the underlying homoclinic structures on the system. We would like to make a few preliminary remarks. One can speak of the class of initial conditions that are within an E-neighborhood of the one double point regime or of the two double point regime. From the more extensive numerical experiments we have performed that are not mentioned explicitly here, we have found that initial conditions within the same class have the same qualitative behavior in terms of homoclinic crossings or breakdown of spatial structures. Secondly, we do not classify the one or two double point regime by their initial spectral configuration alone. For example, one may initially have a gap and a double point on the imaginary axis. However, if as time evolves the gap closes down to a

double point we are in the two double point regime.

4.1. The one complex double point regime We consider the homoclinic orbit which corresponds to one complex double point and onc unstable mode. The initial condition is given by % = 0.5 ÷ e( 1 + i) cos px,,,

(3o)

where L, e, p, x,,, N are the same as before. We have already seen the solution for this under the integrable discretization (7) (see figs. 2a-2c). The solution to the diagonal discretization (6) with initial condition (30) is well behaved as can be seen in figs. 5 a - 5 c which are the surface, phase portrait, and power spectrum for this case. Initially the nonlinear spectrum has one imaginary double point (the spectral plot is the same as fig. 3b). Under the diagonal discretization, as time evolves, this double point opens up into a small cross (see fig. 5d for a representative spectral plot) and then closes down to the double point. We examined the spectrum for 0 < t < 100 and later at a time slice 975 < t < 1000 and we observed the pattern (double point to cross to double point to cross) repeating itself. However we did not detect any homoclinic crossings. The spectrum is in a cross state for a very short time and the double point lasts, comparatively, for a much longer time. Hence we do not see the traveling between center and wings as you would if the cross persisted for a sufficiently long time. (It must be emphasized that the time interval between spectral plots was chosen sufficiently small to ensure that no crossings were missed.) In examining the double points along the real axis, we observed that no spreading to higher nonlinear modes in A-space occurs. When we decrease the size of the mesh to N = 16 (and increase the perturbation strength) and examine the solution for the initial condition (28), in the same time frame as before, it is still well behaved. Four homoclinic crossings occur for ()_< t < 100. The spectral configuration is mostly in the cross state

D.W. McLaughfin, C.M. Schober / Discretizations of the NLS

457

1.0 0.5 //

<

0.0

-,,..

-0.5 -1.0 -0.5

(a)

i

.

.

.

.

0.0 (b)

2.0t''

. . . . . .

'

I

'

.

.

.

.

.

.

.

.

I

'

'

'

'

'

'

'

'

'

I

'

'

. . . . . . .

I

'

'

'

.

.

~

.

1.0

0.5 A

1.5

. . . . . .

H ~

1.0

.~ 0 . 5

0.0 0

10

20

30

n

40

50

(c)

+

04

-0.2

+

0.2

0.4

DI SC

(d) Fig. 5. (a) The solution of the D D N L S for 0 _< t _< 100 and N = 32 with initial condition q(x,O) = ~1 ~_ e(1 + i)COS pX. The surface q(x, t) showing an excitation recurring. (b) The phase portrait of q(O, t) with the inner right lobe visible. (c) The power spectrum of q(0, t). (d) The nonlinear spectrum at t = 1000 showing one very small cross.

458

D.W. McLaughlin, C.M. Schober / Discretizations of the NLS

(a)

4

2 0

-2 -4

I

_ _

,

I

0

2 A (b)

+

~

"~

*

-~

+

*+

4.

4~

4

+ /

-0.4

0.2

O

0.2

O.4

DI ~C

(c) Fig. 6. (a) An e - n e i g h b o r h o o d of the homoclinic solution of the I D N L S with initial condition q(x,O) = all + 4i(e 1 sin (hi ei'St cos plx + e 3 sin (b3 e i~3 cos p3x)] for 0 _< t _< 100 and N = 32. The surface q(x, t). (b) The phase portrait of q((], t) showing a double separatrix. (c) The nonlinear spectrum with two purely imaginary double points.

459

D. H~ M c L a u g h l i n , C.M. S c h o b e r / D i s c r e t i z a t i o n s o f the N L S 7

b 5

4 ." 3

?

.

.

. "'

. ,"

."

... '.:i"7,." i;- . .i..: .' ;: : ". ~ ..~ ~.:

•~ , . . ~ ) ' . ~ k ~

~ 1

..

t;'~:"."'.'i:

~_'::.

,

..

.

..

..

"..

"

. •

.

.

.

~ . ~

~:~:i: : "'i .' -" . ' .

~".~,'~..

":i'~.:?~'"~ :-:

:

.i

-2 -3

-4

-...

-5

(a)

-f2 -~ ,3

.

t

I -.1

i

I ~

I

I ,1

12

I .

r

I

,3

,a

i

15 .

i

I .b

i

,

.7

(b) 2.0

.3

1.5

o

1.0

°IV '.0

,

o

, , , I , ,

10

20

30

40

50

(c) Fig. 7. (a) The solution of the D D N L S for 0 ~ t < 100 and N = 32 with initial condition q ( x , 0) = a[1 + 4i(E 1 sin (~)1 ei'~l COS p ] x + e 3 sin &3 ei63 cos P3X)]. The surface q ( x , t). (b) The chaotic motion in the phase space of q(0, t). (c) The power spectrum. (d) The nonlinear spectrum at t = 32 showing a cross-gap configuration. (e) The nonlinear spectrum at t = 42 showing a gap-cross configuration. (f) The nonlinear spectrum at t = 46.6 showing a c r o s s - c r o s s configuration. (g) The nonlinear spectrum at t = 49 showing a g a p - g a p configuration.

and when it does switch to a gap, the gap is very small. We do not consider these homoclinic crossings significant since the spectral configuration is in the gap state for such a short time that the waveform does not randomly flip between center and wings and in general does not produce irreg-

ular behavior. We do not present the surface, phase or spectrum, as they are qualitatively the same as figs. 5a-5d. The instability and homoclinic orbit associated with one complex double point seems an insufficient mechanism to trigger temporal irregularity in this time series.

D.W. McLaughlin, C.M. Schober /Discretizations of the NLS

460

,

0.4

--0,2

0

0.2

0.4

-04

0.2

,q

DI SO

~1 ,~ 02

-

~ 0.4

k Dic~"~

(e)

(d)

+

, $

~ 04

~

M -0.2

~

+~-~-~-,~-0

02

04

0.4

0.2

02

0.4

D~ %C

D'SC

} (g)

(f) Fig. 7. (Continued)

4.2. The two complex double point regime

where N = 3 2 ,

L=4~/2"rr,

2v/L, x,,=-p/2+(n-

T h e second and much m o r e interesting case is for a combination homoclinic orbit which corresponds to two complex double points and whose initial condition, within an e-neighborhood, is given by qn = a l l + 4i(e I sin db ei~S~cos plx, + e 3 sin &3 ei<~ cos p 3 x ~ ) ] ,

(31)

a=0.5,

ps=2pl=

1)p/N, e I = 10 - 4 , e3 = 10 -s, and the @s are given by pj = 2a sin d~j, j = 1,3. This initial condition is d e t e r m i n e d in ref. [3]. T h e solution to I D N L S (7) for initial condition (31) is p r e s e n t e d in fig. 6a. T h e homoclinic orbit is characterized by the two m o d e s cos p~x and cos p3x and that the homoclinic orbit exponentially decays to the plane wave. Its corresponding phase orbit which a p p e a r s as a double separatrix can be seen in fig. 6b and its nonlinear

D. W. McLaughlin, C.IYL Schober / Discretizations of the NLS

spectrum showing the two pure imaginary double points in fig. 6c. When (31) evolves under the D D N L S (6) we obtain considerably different results. Figs. 7 a - 7 c are the surface, phase portrait, and power spectrum for (31). One sees a breakdown in spatial structure, the abrupt temporal onset of irregular behavior and the corresponding irregular motion in phase space. Figs. 7 d - 7 g show the various spectral configurations that we found for (31). Moreover we can correlate the changes in the spectrum with the behavior of the surface. Initially there are two imaginary double points. For 0 < t _< 30, only two homoclinic crossings occur. However, for t > 30 there is a rapid generation of homoclinic crossings and their density increases markedly. This occurs simultaneously with the breakdown in spatial structure. The time line (fig. 8) shows that many homoclinic crossings were observed as there are 40 for mode 1, 38 for mode 2 (0_< t_< 100), and they occur (apparently) irregularly throughout the time series. Inspecting the time line of the homoclinic crossings we see that rarely do mode 1 and mode 2 simultaneously cross. However, by varying the initial conditions within the same class it is possible to obtain a case where they are simultaneously crossing. In the time range studied, the third to fifth nonlinear modes become excited and there is an oscillation in their excitation but we did not detect any spreading in A-space further than that. This excitation of the third to fifth nonlinear modes makes it difficult to observe a flipping between waveform patterns located spatially in the center versus ones located in the wings, as was seen for the d a m p e d - d r i v e n NLS [16]. Upon decreasing the perturbation, i.e. N = 42, the solution to D D N L S (6) with initial condition (31) exhibits this selection in the waveform pattern until, as in the 32-particle case, the higher nonlinear modes become active. This flipping can be seen in fig. 9a and the crossings of the homoclinic orbit in phase space in fig. 9b. Here, as expected, irregularity takes longer to set in. Nevertheless, we find 13 crossings for mode 1 and 4

Mode 1 Mode 2

Mode 1 Mode 2

I t-'25

[I I

t='50

Mode 1

I,,1od~2

I I

t='75

I ,I

3b

II

I

413

6'O

II

I I ]

II

II I I 8'5

[1

65

I]

[ I

II

90

I 4'5

I I I I II II II II I I I I I

5~5

80

II

35

I

I[

461

70

II

I] II

11 II

5'0

7'5

I [

I 9'5 I

160

Fig. 8. The time line of the homoclinic crossings of a combination homoclinic orbit ( N = 32),

for mode 2 (0 _< t < 100) as shown in the time line (fig. 9c). When we increase the mesh further to 64, the solution to D D N L S (6) for initial conditions (31) has no significant homoclinic crossings. Small gaps have opened in the spectrum and so the spatial excitations repeat themselves instead of dying down to the plane wave. See fig. 10 for the surface to D D N L S (6), for N = 64, with initial conditions (31). The variance in the regularity of the solution, as the mesh is varied, indicates that the homoclinic crossings depend on the perturbation parameter ~2h 2 = & If ~ is taken too small, ~ < 6critical , o n e obtains regular behavior, while if 6 is too large, 6 > 6 . , we cannot apply our integrable diagnostics. For the one double point regime, as long as we are in the near-integrable scenario, 6 < a . , we obtain regular behavior. But for two double points, there is a critical region, 6critical < < a . , where we are in the near integrable scenario and yet obtain irregular behavior. The many homoclinic crossings that occur irregularly throughout the chaotic time series for two double points, for ~critical < ~ ~ (~* and the corresponding breakdown in spatial structure is in striking contrast to the regular behavior and few homoclinic crossings that were observed for one double point for 6 < 8 . . Thus having detected many homoclinic crossings for two double points, it is

462

D. W. McLaughlin, C.M. Schober / Discretizations of the N L S

0.4 0.2

, . • :.:':"/::••.

.

•

-0.0

..:,.._.-.'.~::1~.. . .. ~ .~ f!

•

> 5 :i;:/['," •

-0.2. ======================...... i!!!:!::::2:,. :::::::::::::::::::::::::

!!!!!iiiiii~iiiiiiii ......

;IIIII

i

I

i....•

....•

'"

"

-0.4

II

[,,llii44

(a)

-0.4

-0.2

0.0

0.2

A (b)

Mode 1 Mode 2

I t--25

Mode 1 Mode 2

3'0

I

t-50

5'5

[

3'5

40

I

4'5

I

60

Mode 1

6'5

[

5'0

I

7'0

[ II

7'5

I "l

I

95

100

Mode 2 t=75

80

8.5

90

(c) Fig. 9. (a) The solution of the D D N L S for 0 _< t _< 100 and N = 42 with initial condition q(x,O) = a[l + 4i(% sin 4'] e i'~ cos p l x + e 3 sin 4'3 ei~3 cos P3X)]. The surface q(x, t) showing an irregular flipping between center and edges. (b) The motion in the phase space of q(0, t). (c) The time line of the homoclinic crossings of a combination homoclinic orbit•

plausible that this more complicated homoclinic behavior is the source of the irregular temporal evolution.

5. Melnikov analysis Certainly our numerical spectral transform establishes that the irregular time series, when two double points are present, contains many homo-

clinic crossings; thus it seems reasonable to pursue the consequences of the presence of such homoclinic structures. In low (2-3) dimensional O D E ' s the usual approach for establishing the existence of homoclinic chaos is through a "Melnikov analysis". In studies of the d a m p e d driven sine-Gordon and NLS equations [1, 2, 12] it is found that the underlying homoclinic structures do play a central role in the associated chaotic p h e n o m e n a and a Melnikov analysis of

D. W.. McLaughlin, C.M. Schober / Discretizations of the NLS

463

and

j=(o

Fig. 10. T h e s o l u t i o n o f the D D N L S N = 64 with initial c o n d i t i o n

q(x,

O) - a[1 + 4 i ( ~ 1sin ~bI e i~' cos

f o r 0 _ < t _ < 100 a n d

pjx

The system obtained by setting 6 = 0 is referred to as the unperturbed system. The unperturbed NLS is a completely integrable system and possesses an infinite number of conserved quantities and we make use of the following three conserved quantities: I = rjocltole d x ,

+ e 3 sin &3 ei#~)cos p 3 x ) ] . The surface

c5 =

q(x, t).

finite dimensional models of these equations has been initiated in ref. [15]. A similar analysis is of interest here were the original system is high dimensional while the perturbation is a conservative one. Here we merely mention our initial ideas, deferring the details to a subsequent publication [8, 9]. The analytical formulas are most simple in a continuous framework. Hence, we begin by first viewing the DDNLS as a perturbation in the continuum limit of the NLS. For sufficiently smooth data, DDNLS is equivalent to an autonomous Hamiltonian perturbation of the NLS: iO, = - & * * - 2}&} 2& - 64,**~x,

(32)

with 6 = ~ h 2. The system can be written in the form

io,(g,.)=j( ~Ha/~¢' 8Ha/St~ *)'

(33)

where the Hamiltonian H a is defined by

= -£L(r
(34)

No=

Z [2 I~,.I

(36)

+ 214,

-- ( ~ ' t ~ 1 2 ) 2 -- 6,0x,2,~12] dx,

(37)

-LL(]Oxl2-

(3s)

I~tl4)dy.

The perturbed system has only the two invariants I and H a. To justify a Melnikov analysis for understanding the type of chaotic behavior we have observed numerically, a description of the geometrical structure of the relevant phase space of the system is needed. To achieve this we reiterate the important association between the complex double points, unstable modes (hyperbolic directions in the infinite dimensional phase space) and homoclinic orbits to justify that the interesting dynamics occurs in a low dimensional phase space. Corresponding to each complex double point there is associated one stable and one unstable direction (and one homoclinic orbit). There exists then one active degree of freedom for each complex double point. Further for the band of spectrum there corresponds one excited mode (center-like) or one active degree of freedom. Real double points do not add to the dimension count as they correspond to closed degrees of freedom where dynamically the motion takes place on a much longer time scale (i.e. excitation

464

D. l~E McLaughlin, C.I1/L Schober / Discretizations of the NLS

of these modes is ~'(e) when the others are G(1)) [10, 11]. Appealing to the nonlinear spectral theory in this manner, it can be argued that for solutions of the NLS with N complex double points, the significant dynamics is well represented in a finite 2 ( N + 1)-dimensional phase space. For a concrete representation of the finite phase manifold see refs. [8, 9]. As a first nonrigorous approach then, we are limiting our considerations to this finite 2 ( N + 1)-dimensional phase space. We do not explicitly construct a finite dimensional model but we use it to picture the situation. All the calculations are carried out by means of infinite dimensional P D E objects and we make an analogy to the finite dimensional theory in order to give interpretation to our calculations. For the N = (1,2) complex double point regimes considered in the numerical experiments, the unperturbed system has a hyperbolic periodic orbit u(t) = a e i(2a:'+~') connected to itself by the homoclinic manifold U(x; a, ~) = W~[TM where is a ( N + 1)-vector of parameters. Explicit formulas for the unperturbed homoclinic manifolds (figs. 2b and 6a) can be constructed via Biicklund transformations from the inverse scattering theory and are given in refs. [8, 9]. Noting that the norm I is invariant for the perturbed system we can reduce out I and its conjugate variable 0. Thus the reduced systems we work in are 2N-dimensional. Under the perturbed flow a hyperbolic fixed point persists and for 6 sufficiently small W~~'~ are 6 close to W(['u, respectively. In order to detect the existence of homoclinic orbits in the perturbed system we compute a Melnikov function which measures perturbatively the distance between Wg" and W~su, up to first order in the perturbation p a r a m e t e r 6, along the direction normal to U(x; ~) within the energy manifold. For the one double point regime, the reduced system is two dimensional and the homoclinic manifold is one dimensional. The tangent space of Wo~'~ = span(J VH 0) and VH 0 is orthogonal to the tangent space of W(~TM. The Melnikov function which measures the distance between Wg~

and Wau at a point P0 on the unperturbed homoclinic orbit along the direction VH 0 is

M, = f

121,,(U(x;a,r,yo))d'c.

(39)

Because of energy conservation, the Melnikov function, M~ is identically zero. This is consistent with the possibility that no splitting occurs in this N = 1 case, and is also consistent with the numerical studies of the one double point regime: the time series is regular and there are no significant homoclinic crossings and no chaos. For the two double point regime, the reduced system is four-dimensional and U(x; a, T/) is twodimensional. The tangent space of W~[TM = span(J VH 0, J VCs). The space complementary to the tangent space of W~ TM is given by span(VH0,VCs). Since H~ is equal to a constant, VH 0 is a direction complementary to the energy surface H s (VHcj is parallel to VHa to ~4/(6)). Thus it is unnecessary to measure the splitting of the stable and unstable manifolds along VHc~ as it gives no information at ~ ( 6 ) [20]. It is only required to measure the splitting along the direction VCs and the Melnikov function which does so is sc

M?= f

(~5(U(x;a,7,a,y,,))d'c.

(40)

sc

We have determined that M z has a unique nondegenerate zero. The geometric interpretation of this zero, the dynamical consequences and the correlation with the results from the numerical study of nonlinear spectrum are given in refs. [8, 9].

6. Conclusion

In observing the behavior of time series for the DDNLS, it is found that one obtains smooth qualitative behavior for a short time before the solution becomes abruptly "chaotic". The compli-

D.W. McLaughlin, C.M. Schober /Discretizations of the NLS

cated homoclinic structure of either the integrable NLS equation or the integrable I D N L S equation can be a potential source of this irregular behavior. We have studied numerically, for several representative sample experiments the nonlinear spectrum and its cross-gap structure in the limited p a r a m e t e r regime of one or two complex double points. We definitely establish that these homoclinic crossings occur throughout the time series and can quantify them as a function of time. Our results show that for the DDNLS, in the p a r a m e t e r regime of one unstable mode, few homoclinic crossings occur where on the other hand, with two unstable modes we detect many continual homoclinic crossings. Moreover we correlate this generation of many homoclinic crossings with the temporal onset of irregular behavior and a breakdown in spatial structure. Finally, we have indicated a Melnikov analysis, that is underway, to investigate analytically such homoclinic behavior.

[5]

[6]

[7]

[8]

[9]

[10]

[11]

[12]

Acknowledgements [13]

We wish to thank E. Overman for the use of his spectral solver package and for his helpful correspondences,

a n d A . C a l i n i a n d N. E r e o l a n i

for m a n y useful discussions. T h e w o r k r e p o r t e d here

was

supported

#DMS8703397

[14]

by

the

NSF

under

[15]

grant

and the United States Air Force

under grant #AFOSRF49620-86-C0130.

[16] [17]

References [1] A.R. Bishop, M.G. Forest, D.W. McLaughlin arid E.A. Overman II, A quasi-periodic route to chaos in a near integrable PDE, Physica D 23 (1986) 293. [2] A.R. Bishop, M.G. Forest, D.W. McLaughlin and E.A, Overman II, A quasi-periodic route to chaos in a near integrable PDE: Homoclinic crossings, to appear Physics Len A (1988) 126, 335-340. [3] M.J. Ablowitz and B.M. Herbst, On homoclinic structure and numerically induced chaos for the nonlinear Schr6dinger equation, Phys. Rev. Len. (1990). [4] B.M. Herbst and M.J. Ablowitz, Numerically induced

[18]

[19]

[20]

[21]

465

chaos in the nonlinear Schr6dinger equation, Phys. Rev. Lett. 62 (1989) 2065-2068. M.J. Ablowitz and B.M. Herbst, On homoclinic boundaries in the nonlinear Schr6dinger equation, University of Colorado preprint (199(I). E.A. Overman II, A.R. Bishop and D.W. McLaughlin, coherence and chaos in the driven damped sine-gordon equation: measurement of the soliton spectrum, Physica D 19 (1986) 1-41. G. Terrones, D.W. McLaughlin, E.A. Overman II and A. Pearlstein, Stability and bifurcation of spatially coherent solutions of the damped driven nonlinear Schr6dinger equation, SIAM J. Appl. Math. 50 (1990) 791-818. Constance M. Schober, Numerical and analytical studies of the discrete nonlinear Schr6dinger equation, Ph.D. Thesis, University of Arizona (1991). A. Calini, N. Ercolani, D.W. McLaughlin and C.M. Schober, Melnikov analysis of numerical discretizations of the nonlinear Schr6dinger equation, in preparation, University of Arizona (1991). N. Ercolani, M.G. Forest and D,W. McLaughlin, Geometry of the modulational instability. Part I: Local Analysis; Part II: Global Analysis, University of Arizona preprint (1987). N. Ercolani, M.G. Forest and D.W. McLaughlin, Geometry of the modulational instability. 111. Homoclinic orbits for the periodic sine-Gordon equation, Physica D 43 (1990) 349-384. N. Ercolani, M.G. Forest and D.W. McLaughlin, The origin and saturation of modulation instabilities, Physica D 18 (1986) 472-474. A.C. Newell, in a private communication, suggested this change. M.J. Ablowitz and J.F. Ladik, A nonlinear difference scheme and inverse scattering, Stud. Appl. Math. 55 (1976) 213-229. N. Ercolani, M.G. Forest and D.W, McLaughlin, Notes on Melnikov integrals for models of the driven pendulum chain, preprint (1988). Hie Tae Moon, Homoclinic crossings and pattern selection, Phys. Rev. Lett. 64, no. 4 (1990). S. Wiggins, Global bifurcation and chaos (Springer, New York, 1988). N.M. Ercolani and D.W. McLaughlin, Toward a topological classification of integrable PDE's, in: The Geometry of Hamiltonian Systems, ed. T. Ratiu (Springer, New York, 1991) 111-129. C. Schober, A B~icklund transform for the discrete nonlinear Schr6dinger equation, preprint, University of Arizona (1991). L.M. Lerman and Ia.L. Umanski, On the existence of separatrix loops in 4-dimensional systems similar to integrable Hamiltonian systems, PMM USSR, 47 (1984) 335-340. E. Previato, Hyperelliptic quasi-periodic and soliton solutions to tile nonlinear Schr6dinger equation, Duke Math. J. 52 (1985) 329-377.