NUMERICAL SIMULATION OF A NONLINEAR WAVE EQUATION AND RECURRENCE OF INITIAL STATES by ALBERT GRANT BUCKLEY B.Sc, Univ. of Calgary, Calgary, Alta., 1966 M.Sc., Univ. of Alberta, Edmonton, Alta., 1968 ..A THESIS SUBMITTED IN PARTIAL FULFILMENT OF .THE REQUIREMENTS FOR THE DEGREE OF DOCTOR OF PHILOSOPHY in the Department of MATHEMATICS We accept this thesis as conforming to the required standard THE UNIVERSITY OF BRITISH COLUMBIA October, 1972 In presenting this thesis in partial fulfilment of the requirements for an advanced degree at the University of British Columbia, I agree that the Library shall make i t freely available for reference and study. I further agree that permission for extensive copying of this thesis for scholarly purposes may be granted by the Head of my Department or by his representatives. It is understood that copying or publication of this thesis for finaneial gain shall not be allowed without my written permission. Department of The University of British Columbia Vancouver 8, Canada Date Oct 3 , m * -Supervisor: Dr. G. D. Johnson ABSTRACT In 1955 Fermi, Pasta and Ulam (FPU) [7] observed an unusual recurrence to i n i t i a l state in numerical solutions of a nonlinear wave equation. Zabusky and Kruskal (ZK) [47] have subsequently found an explanation for this phenomenon based on special travelling wave solutions ("solitons") of the (nonlinear) Korteweg de Vries (KdV) equation. In this thesis we extend ZK's explanation to a similar nonlinear wave equation given by Johnson[14]. We investigate existence and unique-ness of solitons for a (nonlinear) generalization of the KdV equation. (Chapter II) and present computational results to illustrate ZK's soliton explanation of the recurrence, both for FPU's equation and Johnson's equation (Chapter III). In Chapter IV we give some results concerning the stability of the difference schemes used to obtain solutions to the nonlinear partial differential equations. i i i . TABLE OF CONTENTS CHAPTER I THE RECURRENCE PHENOMENON §1. Introduction- The FPU Problem p. 1 §2. The Obvious Causes of the Recurrence p. 3 §3. Zabusky and Kruskal's Approach p. 7 CHAPTER II EQUATIONS OF EVOLUTION §1. Derivation of the Second Equation of Evolution p.13 §2. Existence Results for the Solitons p.18 §3. Uniqueness of the Solitons p.23 §4. A Uniqueness Result for the Evolution Equation p.32 CHAPTER III NUMERICAL RESULTS §1. Numerical Simulation of the KdV Equation p.,3.8 §2. Solitons and the BFdV Equation p.50 §3. Computation of the Individual Solitons p.54 CHAPTER IV STABILITY §1. Frozen Coefficient Stability for the Wave Equations p.60. §2. Numerical Stability of Gaussian Elimination for the p. 71 Evolution Equations §3. Kreiss' Stability Theorem: The Energy Method p.79 CHAPTER V CONCLUSIONS p.89; FIGURES p. 9 2 BIBLIOGRAPHY p. ^ 4 ACKNOWLEDGEMENTS I would like to thank my^supervisor, Dr. G.D. Johnson, for originally suggesting the problem to me and for his guidance and assistance during the preparation of the thesis. My thanks also go to Dr. J.M. Varah for several helpful suggestions. Finally, I wish to acknowledge gratefully the financial support of the H.R. MacMillan Family and the National Research Council of Canada. I The Recurrence Phenomenon §1 Introduction: The FPU Problem We begin by reviex^ing the problem.that Fermi, Pasta and Ul'am (FPU) studied in 1955 [7 ]. They wished to compute the rate of therma-lization for a weakly nonlinear, clamped, vibrating string; that i s , the rate at which the energy would be distributed amongst the various Fourier modes, given some i n i t i a l form for the string. In particular, they considered the string model ytt = (1+eyx)yxx t > 0 > x-e I0*1] y(x,0) = f(x) x e [0,1] yt(x,0) = <}>(x) x e [0,1] y(0,t) = y(l,t) =0 t > 0 (2) V where e is a small parameter. When f(x) = sin T T X and <f> = 0 , the energy is in i t i a l l y in mode 1. Figure 1 was obtained from an approximate solution to (1) which was computed from this i n i t i a l data. (At each of a sequence of times, the energy of the approximate solution to (1) was decomposed (numerically) into its Fourier modes. Each mode was then plotted as a function of time to give Figure 1. This also repeats figures given in FPU [ 7 ] and Zabusky [47].) Clearly the energy is at first shared: i t moves into modes 2 through 5. However i t is then almost completely returned to mode 1. FPU found this behavior quite surprising; however they did not succeed in explaining i t . 2. In order to compute approximate solutions to (1), FPU discre-tized both x and t and solved the resulting difference equations. These approximations are computed on a finite set of points; this fact alone suggests possible explanations for the recurrence. We briefly examine these in §2. We may, however, consider the problem from a different point of view: we can ask i f the fault lies in the model i t s e l f . In particular, notice that the sound speed eY ) i -n (1) c a n become pure imaginary, leading to an illposed elliptic problem. Observing this, Johnson [14] suggested that (1) does not really describe, a vibrating string. He then derived the following model, which i t was hoped would provide a better description: y = [1 + e(l i _ ) ] y t * 0 ' x e > ( 3 ) tt ( 1 + y2}% xx with initial/boundary conditions (2). Here the sound speed is uniformly bounded between 1 and We s t i l l observe a recurrence for this model. The quantita-tive behavior changes but the qualitative behavior is the same: the system returns to its i n i t i a l state. See Figure 2 (which is interpreted in the same way as Figure 1). Thus recurrence does not occur simply because we have chosen a poor model for our string. In order to investigate the recurrence in (1) more closely, Zabusky and Kruskal (ZK) [26,27,44,46,47,49 ] proposed a method of isolating the effect of the nonlinearity from the essentially linear oscillations of the string (i.e. the periodic oscillations of the string (1) with e = 0 ). In this way they hoped that the cause of the recurrence would become more evident. It did, and they exhibited a mechanism to explain the recurrence. In §3 we outline this procedure. Can the recurrence observed in computed solutions to (3) be explained using ZK's method? One of our goals is to show (see Chapter III) that i t may. As indicated earlier, the approximate -solutions to these partial differential equations were computed using finite difference techniques. Thus the numerical stability of the difference equations is important; we examine this question in Chapter IV. §2 The Obvious Causes of the Recurrence We now look at FPU's discretization of (1). We use the notation hD+y(x,t) = y(x+h,t) - y(x,t) hD_y(x,t) = y(x,t) - y(x-h,t) 2hDQy(x,t) = y(x+h,t) - y(x-h,t) . Solutions of (1) are computed on a uniform grid = jh j = 0,1,...,J , (3.1a) t = nk n = 0,1,2,... . (3.1b) where x = 0 and x = 1 . The difference scheme used by FPU is obtained U «J by using second order difference quotients to replace derivatives in (1) at the grid points (x,t) : y(x,t+k) - 2y(x,t) + y(x.t-k) = [1 + eDoy(x,t)]D+D_y(x,t) . (4) We set y(x_,,t ) = y(x_,t ) = 0 for a l l . n . Then the solution to (4) O n J n is well defined once the starting values at t = 0 and t = k are speci-fied. (Because (4) involves values of y at three time levels, starting data is needed at two times. We call (4) a 3-level scheme.) We obtain this from (2), y( X j,0) = f(x ) j - 1,2,. • • » J-l and from Taylor's Theorem (by truncating): j = 1,2,. • • » J - l (We have assumed that <j> = 0 and that y has 3 continuous derivatives at each point (x.,0).) Consider the effects of this discretization. First, we use finite precision arithmetic to compute these approximations, so recurrence could have resulted from the systematic accumulation of roundoff error. However the energy of a conservative system - such as (1) (see Zabusky [47J) - is time invariant. We computed the total energy of the approximate solutions obtained to (.1) at frequent times t . Since there was no evidence of growth (or decay) of the computed energy, we conclude that roundoff plays a negligible role. Secondly we consider aliasing, which is a result of computing a finite Fourier series to represent the leadivig coefficients of the true infinite Fourier series of a function. The phenomenon is discussed in Johnson [14] (see p. 10), where he shows how the (sufficiently) high modes (of the true Fourier series") appear as contributions to the low modes of the computed series. Now Zabusky [45] shows breakdown of the solution to (1): at some time t , y becomes discontinuous. (See also Lax [30 ] for a simpler proof of breakdown.) For such non-smooth functions, It is well known that the higher modes of the Fourier series may not be small. So, at (and after) t,, , the computed energy decomposition may not be correct. In particular, the feed back from the higher modes could appear as the energy recurrence. A primitive test for the presence of aliasing is to sample at unequally spaced points (see Tukey [42], p. 399). So we recomputed the approximation (4), but using a non-uniform x-spacing. The recurrence behavior was not changed; see Figure 3. This suggests that aliasing is not causing the recurrence, but does not prove i t (since aliasing could possibly occur with unequal spacing). In what follows we disregard aliasing, for two reasons. First, there is the indication of Figure 3. Secondly, we will see later how problem (1) may be modified so as to improve the smoothness of its solutions but s t i l l preserve the recurrence behavior. For this modified problem ZK have found an explanation for the recurrence. Again this suggests that aliasing was not the cause of the recurrence for (1) (or for the modified form of (1)). We remark that the ZK explanation (see §3) sheds no light on the non-uniform recurrence of Figure 3. Thirdly, Ford [ 8 ] observed that FPU always used J equal to a power of 2 (perhaps to use a Fast Fourier Transform; see Cooley and Tukey [ 6 ]) and suggested that this accounted for the recurrence. Figures 4a and 4b give the energy decomposition for solutions to (1) with J = 30 and 35 ; a recurrence is observed in each case. Clearly then Ford's suggestion is not correct. Finally, we discuss a finite mechanical system of J - l masses th described by ZK [44]. With y^t) describing the position of the i mass, this system has equations of motion given by y.(t) = K(l + a(yi + 1 - y ^ K y ^ " 2y. + y.^) (5) for i = 1,2,...,J-1 , where a is a constant. Set YgCt) H 0 = ' But n o w -*-et be discretized only in the x-direction (see (3.1a)) and let u^(t) = y(x^,t).. Then the.system of ordinary differential equations h2u.(t) - [1 + f^(u.+1 - u ^ K u ^ - 2u + V l ) (6) is obtained. By comparing (5) and (6) i t is clear that the solution to (4) describes the motion of a finite mechanical system (with K = ~ • , h e = 2ah . The solution to (4) is thought of as the solution to (6) computed at discrete times t = 0,k,2k,... ) . From this point of view, the solution to (4) defines the trajectory of a point in phase space. By the classical "Wiederkehrsatz", or Recurrence Theorem, this point must.return to its origin after some sufficiently long time. (More precisely, given any neighborhood of the original point, there exists a time t at which the trajectory will again lie within that neighborhood.) This is known as a Poincare" recurrence; in our problem, this would certainly appear as the energy recurrence. However, as FPU point out, the observed recurrence time is much too short to be a Poincare" recurrence. Thus we see that i t is difficult to explain the recurrence as an immediate result of discretizing a continuous problem. In the next section we consider ZK's explanation. §3 Zabusky and Kruskal's Approach We now outline how ZK explained the recurrence observed in approximate solutions to (1). In the next two chapters we extend this explanation to include the recurrence observed for the model (3). To begin, ZK introduced the equation 8. yt t = ( 1 + e yx) yxx + Byxxxx ( 7 ) for t > 0 and x e [0,1] , with initial/boundary conditions (2) and 3 > 0 a constant. We examine their motivation for the fourth derivative term. We mentioned earlier that Zabusky [45] has shown that (1) fails to have a global solution (for a class of i n i t i a l conditions including the example of §1);. y becomes singular at t . Kruskal [26] observes that (1) describes a xx B perfectly flexible string; i t allows unlimited bending and therefore the singularity. The fourth derivative term resists bending and therefore prevents formation of the singularity. (ZK do not prove existence of smooth solutions to (7). They simply assume them and consider the fact that an explanation for the recurrence is obtained as justification for the assump-tion.) h2 2 In fact ZK assume g = = 0(e ) . For, recall the mechanical system of §2. Suppose there is a function y(x,t) defined for t 5- 0 and x e [0,1] with y(x^,t) = y^(t) . (x^ given in (3.1a)). Assume motion of the (i±l)St mass may be described in terms of the 1^ by Taylor's Theorem: y(*i± 1,t) - y(x.,t)± hyx(x.,t) + ... + j£ yx x x x(x.,t) + o(h5) . Substitute into (5) and recall K = -^ y , e = 2ah . Then h ytt(x.,t) - [l+ e{yx(x.,t) + 0( h2) } ] x [yx x(X i,t) + ^ yx x x x( x . , t ) +0 ( h4) ] 9. 2 Dropping the subscript i and ignoring terms 0(h ) or higher gives (1). 2 But retention of the 0(h) terms (recalling e = 0(h)) gives (7) with h2 h2 $ = — . We assume throughout the following that B = -£2^ 0 • ' We see then that (1) and (7) approximately describe the motion of the same mechanical system. For this reason ZK suggested that computed solutions to (7) would exhibit recurrence. Figure 5 illustrates that they do. ZK then explained this recurrence as follows. They sought to mathematically isolate the effect of the nonlinearity from the oscillations (i.e. the linear oscillations when e = B = 0 ) of the string. Chapter II, §1 describes this procedure for a more general form of (7); here we simply state that they obtained the Korteweg de Vries (KdV) equation u + uu + 62u =0 (8) t X X X X to describe the effect of the nonlinearity. They then solved (8) numerically. A recurrence appeared. As well, a mechanism to explain the recurrence was observed. (This mechanism was consistent with several observations made in solutions to (1); we discuss this point in §1 of Chapter III.) The following took place. From the i n i t i a l data given for (8) several individual wave pulses appeared. Each of these solitary wave pulses, for which Zabusky has coined the term "soliton", moved with constant velocity and amplitude and preserved its identity even through interactions with other solitons. Interaction of 10. these solitons (in a particular phase relationship) caused the recurrence. (This may be easily observed in the computations of Chapter III, §1). Finally ZK established that these wave pulses were in fact par-ticular solutions to (8). This conclusion was motivated by the fact that the soliton interactions were demonstrating an essentially linear behavior. They proposed solutions u(x,t) = s(x-ct) to (8), where c is a constant and lim s(x) = s_ , s_ a constant |x|-*=° CO CO lim s'(x) = lim s"(x) = 0 |x|-x» |x|-x» A solution to this problem is given by (9) with and s(x) = A = L1-2 s + A sech 3 ( c - s J CO 2 126 (10) ( I D (In Chapter II, §3 we show these are the only solutions.) By comparing the solutions (10) and the observed pulses ZK concluded that they were the same. (see Zabusky [46])• As indicated, we will show that a similar procedure can be applied to (3). We observe solitary wave pulses in solutions to a partial differential equation corresponding to (8) and see how they explain the recurrence. We show they are special solutions to this PDE. We are 11. unable to give an exact expression for the solitons but we do give existence and uniqueness results (see Chapter II) and a method for computing approximations to these solitons (see Chapter III, §3). Before examining (3), we remark that considerable further work has been carried out on the KdV equation. These investigations are not directly concerned with the FPU problem but do involve the soliton or travelling wave solutions of (8) observed by ZK. We briefly summarize these results. In [9], Gardner, Greene, Kruskal and Miura (GGKM) present a method for solving the KdV equation (analytically). They assume solutions which approach a constant sufficiently rapidly as |x| -> 00 and can predict, from arbitrary i n i t i a l conditions, the solitons which will emerge as t ->- °°. The case of two solitons has also been treated by Lax [32]. In addition, GGKM and Su have published a series of papers [34,35,41,10,28] investigating-pro-perties of the KdV equation. A summary of their results appears on the first paper of the series [34], One of their principal results is that the eigen-values of the Schrodinger operator L, defined by Lv = v + au(x,0)v xx where a is a constant and u(x,0) is the solution to (8) at time 0, are constants of motion for (8). In fact Lax [32] gives a general method of associating a nonlinear evolution equation u = K(u) (K depending only on u and its x-derivatives) with an operator L so that the eigenvalues of L are constants of motion of the evolution equation. Using the particular relation between the KdV equation and the Schrodinger operator, Lax [32] 12. also demonstrates (see also [33]) that the solitons appearing in solutions to (8) have speeds c. related to the X. by c. = 4A.. Related results 3 3 3 3 for the KdV equation appear in Karpman [20 and 21], Shih [38], Benjamin [1] and Hirota [12]. Finally, Zabusky has compared solitons appearing in numerical solutions of (8), f i r s t , with the analytic results of GGKM [9] and, second, with observations made from actual shallow-water wave experiments. The first [see 48] establishes an adequate means of computing numerical solutions to (8) (with boundary conditions at infinity), and the second [see 50] shows that the solitons observed numerically can in fact be physically realized. 13. II. Equations of Evolution §1 Derivation of the Second Equation of Evolution. We have indicated in Chapter I that ZK's investigations of equation (1.1) led to the KdV equation. Here we will see how equation (1.3) leads to a similar evolution equation. In the remainder of this chapter we will give some results concerning existence and uniqueness of solutions to these, and more general, equations of evolution. Consider equation (1.3). Johnson [14] has shown that i t does not always have a global solution so, as discussed in Chapter I, we will add a fourth derivative term. As before, this does not guarantee a global smooth solution. (For a global existence theorem in a non-conservative environment, see Greenberg, MacCamy and Mizel [IT].) From (1.1) and (1.3) we get the general form yt t = *2 ( yx) yxx + g W ' t > 0 • X E [ 0 ' 1 ] ( 1 ) 2 where $ (w) H 1 + eg (w) , (2) with initial/boundary conditions given earlier, (1.2). Following ZK we now define w = y , v = y . Assuming the required smoothness for X t solutions to (1), we get the coupled equations w = v (3) t x 2 v„ = $ (w)w + gw • (4) t X X X X (1.1): Equation 1 of Chapter I. 14. We now introduce the Riemann Invariants of (1) with 3 = 0 . These are the functions constant along the characteristic directions 31 9 x . of that equation and are given by 1 W r ± = | ( ± v + J 0 $a)d?) Multiplying (3) by ±$ and adding to (4) we. obtain ±v + $w = ± $(±v + $w ) ± Bw , t t X X X X X ( r ± ) t = ± $ ( w ) ( r ± ) x ±| B w x x x . (5) We would like to obtain an equation involving only r+ and r_ . Now w an equation from which we cannot expect to get w explicitly as a function of r and r . In particular, consider the case $2(0 = 1 + e( l - ) (l+O* (that i s , $ as i t appears In Johnson's bounded sound speed wave equation (1.3)). By Liouvilles First Theorem (see Ritt [37], P- 20), i f $ has an 15, Integral- in terms of elementary functions, i t can be written as f<3> = v (w) + Y c. log v. (w) ; o ^ 1 i where the c^ are constants, are algebraic and r is a nonnegative integer. It can be shown that such an integral for $ must have r > 0 ; however we have been unable to find such an integral for $ . Nor could we find an integral for $ which was not in terms of elementary functions. In general, then, an explicit form for w is unlikely and we must settle for something less. We seek an approximation to w : r+ + r_ = io (1 + eg(O)* dE = / o (1+ 0(e)) % d£ = w + 0(e) . 2 We will now assume, as did ZK, that r+ = 0(e ) for a l l t . Then we have w = r_ + 0(e) Now consider the second of equations (5). Let g be continuous and satisfy g(r + 0(e)) = g(r) +.0(e).. Then, writing r for r_, and re-2 calling r = 0(e ), r = -$(r + 0(e))r - |(r + 0(e)) t x L xxx = -[1 + eg(r + 0(e))] %rx - | r ^ + 0(6e) - - [ 1 + f g(r + 0(e))]rx - | rx x x + 0(e2) + 0(Be) = - I I + f g(r)]rx - | rx x x + 0(e2) + 0(ge) . 16. 2 If we now ignore the small effects e and Be we get the equation r + (1 + | g(r))r + f r = 0 . (6). t 2 x 2 xxx Finally we wish to make a change of the independent variables, principally in order to include the nonlinearity e in the time scale. Computationally this will mean that we do not have to compute for long times in order for the effect of the nonlinearity to become evident. Let £ = x - t , n E — t and u(£,n) = r(x,t) . Then for notational convenience, replace E, ,n in the resulting equation by x,t and set e 2 2 As discussed earlier, 8 = 0(e ) so 6 is a small parameter, 0(e) The result is the generalized Korteweg de Vries, or GKdV equation, u„ + g(u)u + <S2u = 0 . (7) t X xxx In particular, the equations of interest are the KdV equation (1.8) (g(.u) = u) , and the bounded sound speed variation, or BKdV equation u + (1 1—)u + 62u = 0 . (8) t I — X xxx /l+u This derivation has followed ZK closely. However the following point, which was not considered by them, should be mentioned. We have made certain approximations in the above, and then disregarded terms of 17. 2 2 order e " or higher. However notice that g = 0(e ) and yet we did not-discard the term -jr gr . I n fact, since we certainly wish to retai-n-2 xxx this term, perhaps we would be more correct to also retain other terms 2 of order e ; in particular, those in the approximation of <J>(w)rx We consider the specific case of the FPU equation (1) where g('u) = u (which we will refer to as (II. 1 FPU)). It follows that w = r - | r2 + 0(e2) , (9) 2 so, retaining 0(e ) terms give a perturbed form of the KdV equation u + (u - f u2) u + 62 u =0 . (10) t 2 ' x xxx A natural question, then, is whether solutions to the KdV equation and (10) differ significantly. We will discuss this in Chapter III. We will also have cause in Chapter III to consider the approxi-mation to the g-term more closely. Computationally we will see that 2 the 6 term in the GKdV equation (7) is the most important regarding the soliton behavior, so i t is reasonable to ask i f inclusion of the neglected term 0(ge) would affect our results. Using (9) we see that inclusion of this term leads to the equation u + uu + 62(u - f u2) =0 . (11) t x 4 xxx We could also obtain modified versions of the BKdV equation (8). Because of their complexity, and because we do not gain anything from these 18, ihodifications., even for the KdV equation (as we will see in Chapter III,. §!•)>, we omit those derivations. §2 Existence Results for the Solitons We have seen how the wave equation (1) leads to the GKdV evolu-tion equation (7). Our motivation for the study of these equations is the recurrence phenomenon discussed in Chapter I. In this respect, §3 of that chapter indicates that our interest lies in solutions of (7) which are travelling waves. That i s , solutions u(x,t) = s(x-ct) where s is defined and has a continuous third derivative on (-00 ,00) and where s satisfies the boundary conditions (1.9) which, for convenience we restate: lim s(x) = s^ , a constant, (12a) x->-±o° lim s'(x) = 0 , (12b) X->±OD lim s"(x) = 0 . (12c) x->±<=° In the following, we will use the term "soliton" to denote such a solution, to the problem ( 7). (We will use the term "solitary wave pulse" to refer to the wave pulses which we will observe in approximate solutions to (7). See Chapter III.) In this section we show that for certain values of c , 19. no solitons exist; in the following section we show that when solitons exist, 2 they are uniquely determined by c (and 6 ) and the boundary conditions (12a) to (12c). (The exact meaning of the term "uniquely" will be defined in §3.) We distinguish two classes of solutions. A soliton will be called positive (or said to have positive amplitude) i f i t attains the value A = sup Is(x) - s I at some point x with s(x ) > s . A negative soliton is defined r o O oo ° similarly. For the following theorem we will assume in fact that s(x) - sm is of constant sign and nonzero for a l l x e (-00 ,00) . It w i l l be established in §3 that this assumption is justified. Other equations of evolution have been studied and results given concerning travelling wave solutions. In particular, Cole [ 5] has found explicit travelling wave solutions to Burger's equation u + uu - vu =0 , . (13) t x xx v ' V a constant, and Canosa [ 4] has shown how solutions of Fischer's 20, population propagation equation 2 u^ + k(u - u) - nu = 0 , t xx X) and k constants, degenerate into solitons after long times. With such observations for similar evolution equations and knowing the results of ZK, one might expect to find solitons for the BKdV equation (8). We will look at numerical results regarding this in Chapter III. However we now wish to show that for certain values of c , the GKdV equation (7) has no solitons moving with speed c , For comparison, recall that in the evolution equation u = cu , c is the speed of signal propagation. A similar behavior is found in (7). 3 Theorem 2.1: Let s be a (real valued) function in C (-00 , 00) which satisfies the boundary conditions (12). Suppose that s(x) - s^ is of constant sign and nonzero for x e (-«,,°°). Let g be continuous. (i) If m •< g(G) < M for a l l 6 e(-°°,«) , then only for m < c < M can the GKdV equation (7) have solutions u(x,t) = s(x-ct) ; (ii) If mQ « g(8) for 8 > s^ , then only for c > can the GKdV equation (7) have solutions u(x,t) = s(x-ct) which are positive solitons. Broof: Suppose u(x,t) = s(x-ct) is a solution. It is immediately clear from the boundary conditions (12) that there must be some finite point x for which s"(x ) = 0 . For, i f not, then by continuity s" o o is of constant sign and s' must be monotone. This is not consistent with, s' .being continuous and having limit 0 as |x| ->- °° 21. On the other hand, we now show that for c i (m,M) , s" can never be zero. From (7) s must satisfy the ordinary differential equation (ODE) fiV' - cs' - g(s)s' . (14) Integrating (14) we get, for some constant K , <52s" = cs - /g(s)ds + K = f(s) . (This has the intuitive interpretation that the convexity of a soliton s is simply a function of its height.) Now |f(s) = c-g(s) , so i f c 4 m or c M , is of constant sign (and possibly zero) . In any case, f is monotone. Note that (by continuity) f(s0 0) = 0 Now suppose that for some (finite) x^ , s"(x^) = 0 ; with s^ = s (x^) , f(s^) = 0 . Recall that by hypothesis s^ f s^ , so by the monotonicity of f , for all x with s(x) between s, and s , s"(x) = 0 . But i 00 s -»- s as x oo ; by' continuity then there is a finite point x„ such 00 Z that x > x ' implies s(x) is between and s . Thus, for Z I 00 x > x^ , s"(x) = 0 . However this is consistent with the boundary con-ditions only i f s(x) = s for x > x^ . Thus by the hypothesis that s(x) - s ^ 0 , (i) follows. 22. -To show (ii) we repeat the above argument, in the case c < m , to show that there cannot be any finite point x., with o ± s"(x^) = 0 and s(x^) > . q.e.d. The proof of the theorem can be simplified i f we consider c i [m,M] (rather than c i (m,M)) when we show there can be no point x with s"(x) =0 . For then f cannot be zero and f must be strictly monotone. Then (since s(x) - s ^ 0 for any x ) , i t is immediate oo that for a l l x , s"(x) 4- 0 . It is of interest to note that this theorem was motivated by numerical work. Attempts to compute solutions to the BKdV equation (8) with c > 1 became rapidly unbounded. As we have noted, the theorem applies to solitons. Notice that the conclusions of the theorem agree with ZK's results. In particular, motivated by numerical results, they sought positive solitons. Theorem 2.1 says that their speed must be more than s^ , which agrees with ZK's result, (I.11), that A= 3(c-S r a) > 0 . The theorem makes no statement that solitons exist i f c lies within the bounds on g . In fact, in the following section we will see that solitons with speed c < .s do not exist for the KdV equation, CO even though g has no lower bound. At this point one consequence of the theorem should be pointed out. In the case of the BKdV equation soliton speeds are bounded by 1. This is not so for the KdV equation (since g(u) = u is unbounded), and this difference will be reflected in the length of time i t takes the solitary wave pulses (which are in fact solitons) to regroup to cause the recurrence (as described in Chapter I.) 23. To conclude this section we note that Theorem 2.1 may.be genera-lized to other evolution equations, including a generalization of Burger's equation (13). ! Theorem 2.2: Consider the evolution equation u + g(u)u + =0 , n > 2 , (15) Z X 3xn where v is a non-zero constant. Suppose s is a (real valued) function in Cn(-o°,K>) which satisfies the boundary conditions s -> s^ , s', s", s^n_1^ •> 0 as |x| -*• 00 . Suppose s(x) - s^ is of constant sign and nonzero. Then the conclusions (i) and (ii) of Theorem 2.1 apply to the solutions of equation (15). Proof: As for Theorem 2.1. §3 Uniqueness of the Solitons We now consider the GKdV equation (7) and show that travelling wave solutions u(x,t) = s(x-ct) which satisfy the boundary conditions (12) (that i s , a l l solitons) are uniquely determined by c and sM We first clarify what we mean by uniqueness. Recall that s = s(x) must satisfy the ODE (14): 62s'" = cs' - g(s)s' . (14) Clearly, since x does not explicitly appear in (14), solutions to (14) are translation invariant. That i s , for any solution s(x) and for any constant p , s(x+p) is also a solution. Furthermore, (14) has the 24. trivial solution s(x) = s^ . We will also see later that there may be both positive and negative solutions.with .the same speed c and background s^ Thus we say that the solution to (14) is unique fro mean that a l l nontrivial positive (negative) solutions to (14) which satisfy the boundary conditions (12) are translates of some particular positive (negative) solution. We recall at this point that a l l solitons are required to possess three continuous derivatives at al l points x e (-°°,°°) value problem (12) and (14) to a suitable i n i t i a l value problem, and then to apply well known uniqueness results from the theory of i n i t i a l value problems for ODE's. For this reason we require, the following preliminaries. Consider a general second order ODE We let y = y,y~ = y' and we write (26) as a system of two first order The fundamental idea of our proof is to relate the boundary y" = f(x,y,y ') (16) ODE's: 07) y^ = f(x,y1,y2) Define the vector y = and the vector valued function f = 1 f \ J Then (17) is a special case of the general where f ^ x . y ^ ) = y2 vector ODE y' = f(x,y) (18) 25. We let | • | denote the Euclidean vector norm. Definition: We say the function f satisfies a Lipschitz condition (with respect to y ) in a region R of (x,y)-space i f there exists a finite Lipschitz constant L such that |f(x,y) - f(x,z)| <L|(y-5)| (19) for a l l (x,y) , (x,z) e R . (We also say that the ODE (18) satisfies the Lipschitz condition.) We now state the basic uniqueness result for i n i t i a l value problems (see for example Birkhoff and Rota [ 3 ] , p. 103). Theorem 2.3: Consider the system of ordinary differential equations (T8) . Suppose i t satisfies a Lipschitz condition in a region R . Then there is at most one solution to (18) which satisfies a given i n i t i a l condition y(x ) = y in R o o Now suppose the scalar valued function f in (16) does not depend on y1 ; i.e. f(x,y,y') = f(x,y) . Assume f satisfies a Lipschitz condition with respect to the second variable y . Then with f^(x,y^,y2) y^ , and f^ = f it is clear that f = must also satisfy a Lipschitz condition. Thus we have Corollary 2.4: Consider the second order ODE (16) and suppose f(x>y»y') = f(y) satisfies a Lipschitz condition in a region R = {(x,y) ; |y[ < b} . Then (16) has at most one solution satisfying the i n i t i a l conditions y(x ) = y and y'(x ) = y. , where |y I < b o o o 1 '•'o1 We also require the following two lemmas. 26. Lemma 2.5: Suppose y is a function defined for a l l x e (-03,00) and which is symmetric about two distinct lines x = x^ and x = x^ Then y must be periodic and have period- dividing Sjx^-x^l Definition: We say that the ODE (16) is reflection invariant i f for every function y(x) satisfying (16), z(x) = y(-x) also satisfies (1 ) Lemma 2.6: Consider the second order ODE (16) and suppose that y is a solution satisfying the i n i t i a l conditions (20, y'(xo) =0 for some finite point Xq . Suppose that the ODE (16) is reflection invariant and that f = f(y) satisfies a Lipschitz condition in a region R containing (x0>yo) • Then the function y is the only solution of (16) satisfying (20), and in fact y must be symmetric about the line x = x o Proof: By Corollary 2.4, y is uniquely determined by the i n i t i a l condi-tions. Let yg(x) = y(x+xQ) . Then, since f = f(y) , solutions to (16) are translation invariant and y is the unique solution to (16) satisfying s the i n i t i a l conditions y (0) = y , y'(0) = 0 . But the ODE (16) S O S is reflection invariant by hypothesis. Therefore z(x) = y (-x) also satisfies (16). But z satisfies the same i n i t i a l conditions as yg , so z = yg . Thus yg is symmetric about 0 , and therefore y is symmetric about x . q.e.d. 27. We now show the uniqueness of solutions to (14) which satisfy the boundary conditions (12). We will assume that g is a given continuous function defined on (-00 , 00) Let G(s) = /g(s)ds . Then integration of (14) gives 62s" = cs - G(s) + K (21) where, by continuity and the boundary conditions (12) , K = G(s ) - cs CO oo Equation (21) is now of the form (16) with f(x,s,s') = f(s) = cs - G(s) + K We now multiply (21) by s' and integrate again. Then %62(s')2 = j s2 - H(s) + Ks + L (22) where H(s) = /G(s)ds c 2 and L = H(s ) - TT S - K S CO / CO 00 Using (22) we define the function <}.(t) = | t2 - H(t) + Kt + L . (23) We now give the uniqueness result. 28. Theorem 2.7: Assume g is continuous on (-co,co) and suppose G satisfies a Lipschitz condition (19) on any region R = {(x,s): |s| < b} . Let the roots of <j> satisfy the condition that there is a unique root of tj) which is the least root of <f> greater than s . Then there is at most OO one nontrivial positive solution to the ODE (14) which satisfies the boundary conditions (12). Furthermore this solution has the following properties: (i) there is a point x at which s'(x ) = 0 ; o o (ii) s(x ) = v, ; o 1 (i i i ) s is symmetric about XQ ; (iv) s decreases monotonically from Xq to °° A similar statement may be made for the negative solitons. • Proof: Let s be such a solution. Since s is continuous, at the point x at which s attains its maximum we must have s'(x ) = 0. r o o By (22), SQ = S(Xq) must be a root of <f> which is bigger than (Note that s is always a root of <j> ; this corresponds to the trivial CO solution). It may be easily shown that the ODE (21) is reflection invariant. Therefore, by Lemma 2.6, s is uniquely determined by the value s and s is symmetric about x . Thus the number of solitons o J o is limited by the number of roots of <j> Now suppose that Sq = v for some zero v of <j> with v > Since s -+ s as Ixl -> »>' there must be some point x. > x for which co r 1 O Sl ~ s(xi ) = vi • By C 2 2 ) » s ' C x ^ = 0 • B u t then, by Lemma 2.6 the solution s is uniquely determined by s^ and s' (x^) and is symmetric about x^ . Thus s is symmetric about Xq and x^ ; by Lemma 2.5 s is periodic. Clearly then s cannot satisfy the asymptotic boundary con-ditions. We conclude therefore that s = v, and the (positive) solution to o 1 (21) is unique (to within'translation) and symmetric about x 29. Finally consider (iv). By an argument similar to the above, there can be no points x > x for which s(x) g s . thus for a l l o o x > x , s(x) < s(x ) . But d> has no. roots, between v-., and s , O O 1 oo » so there can be no points with s1(x) = 0 except x^ . Therefore s decreases monotonically for x > Xq , unless perhaps there is some point x„ > x with s(x~) = s . In this case (as before), the solution must 2 o 2 oo be symmetric about Xq and x^ and hence periodic. Thus this case is not possible and (iv) holds, q.e.d. We now discuss this result as i t pertains to the KdV and BKdV equations. As well, we use the function. <j> to extend the (non-) existence results of §2. First we state Corollary 2.8: Suppose there exists a positive (negative) solution to (14) which satisfies the boundary conditions (12). Let (J> be given by (23). Then for a l l s^ < t < (u^ < t < s^ , where denotes the biggest zero of rj> less than s^ ) <J>(t) must be positive. In particular,if any such solutions to (14) exist, we must have <j>"(s ) 5- 0 Proof: From (22) we have %62(s'(x))2 = <Ks(x)) (24) for a l l x . Since the left hand side of (24) is positive, the rip;ht hand side must be also.But by (iv) of Theorem-2.7, the values taken on by s(x) lie in (s^.v^) (in the positive case). Thus the first claim follows. The second is clear, for i f ^"(s^) < 0 then, since <f>(sco) = 0 , there is a neighborhood about s for which <j>(s)<0 . q.e.d. 30. Now consider the KdV equation. In this case g(s) = s in (14), and <j> is a cubic with roots s , s and 3c - 2s' . Thus, since there CO CO co is only one nontrivial root of <j> , for given c and. s^ , there cannot be both positive and negative solitons.. Now recall ZK's result. In the case that c > s^ , they gave an explicit form for solutions of (14), namely s(x) = s + A sech2 T" , A ( 2 5 ) 2 where A = 3(c-s ) , A = — r — . (see 1.10). Notice s'(0) = 0 and s(0) = 3c-2s , in accordance with CO Theorem 2.7 and the roots of A . In the case c < s , (25) is not T CO 2 applicable since A < 0 . However we know from Theorem 2.7 that in this case a l l solitons must have negative amplitude; a reasonable sugges-tion for such a solution is a function s of the form (25) with 2 2 126 A = r~ • However such s do not satisfy the ODE (14). In fact, -A solitons of negative amplitude do not exist. This follows from Corollary 2.8 since c < s implies <t(t) < 0 for 3c-2s= u < t < s . CO CO ^1 CO Now consider the BKdV equation. The situation here is somewhat different and less simple. In this case we get 4>(t) = "M" t2 + (t sinh"1 t - /i+t2) + Kt + L (26) where K = (l-c)s - sinh ^ s 31, Now, in addition to s^ , <j> may have 2 or more roots, and may have roots on both sides of s^ . This suggests that there m"ay now be solitons of eq.ual speed, with the same background,, andi with, both-positive, and negative amplitude. This can in fact occur (supposing of course that any solitons exist at a l l ) . For example, suppose = 0 and s satisfies the ODE (14) and the boundary conditions (12). Let s^(x) = -s(x) . Then s^ , satisfies (14) and the same boundary conditions and has amplitude of opposite sign to s On the other hand, for certain values of c and s , even ' CO ' with 0 < c < 1 , there may be no solitons. From (26) we get D>"(S ) = — (1-c) . (27) CO Thus by Corollary 2.8, i f ^ < 1-c , there are no solitons with speed 00 c satisfying the boundary conditions (12) specified by s^ . As well, (27) gives a simple way to test i f both positive and negative solitons can exist (with the same c and s^ ) .. For i f both exist, <j>(t) must be positive for t on both sides of s^ (and near), so ^"(s^) > 0 We summarize with Corollary 2.9: Let g(s) = 1 ~ — and let <J> be given by (26). A+s2 Then (14) has no solution satisfying (12). when A ? < l-c , 32. and (14) can have both positive and negative solutions of the same speed c and background s which satisfy (12) only i f Finally we note that in Chapter III, §2, negative solitons will arise naturally in solutions to the BKdV equation (8). §4 A Uniqueness Result for the Evolution Equation We now wish to consider the uniqueness of solutions to the GKdV equation (7). First we will change the problem slightly and consider solu-tions of (7) which are periodic on some interval [0,L] and which satisfy the periodic i n i t i a l data u(x,0) = f(x) , where . f(x) = f(x+L) for a l l x . For computational purposes this problem is actually more appropriate than considering solutions defined on (-03,00) , as in §§2,3, for we can compute solutions only on a finite interval. We will discuss this point further in Chapter III. For functions f with period L we use the norm given by ||f||2 = J f2(x)dx . o Now consider uniqueness. Sjoberg [40] has shown that smooth periodic i n i t i a l conditions uniquely determine periodic solutions to the KdV equation which exist for a l l t >0 . (See Benjamin, Bona and Mahoney [2] for a criticism of the existence part of Sjoberg's proof. However . existence of smooth solutions has also been shown by Temam [43] and in a 33. series of papers by Kametake [15,16,17,18].) The uniqueness part of the proof has also been given by Lax [32]; however Lax fails to show the exi-stence of a time-independent a priori bound on max |u (x,t)|. This X E [ 0 , L ] X bound is required in the proof. For the GKdV equation we were unable to obtain such a bound, but, under the assumption that i t exists, we will extend the uniqueness proof of Lax and Sjoberg to equation (7). Before we do, we look briefly at the bound 8. Sjoberg obtained 8 from the first 4 of an infinite series of conservation laws for the KdV equation given by GGKM [35]. It is because we could not generalize more than three of these laws that we are unable to obtain an a priori bound 8 for max | | for the GKdV equation. In fact, for g(u) = , p a positive integer, i t has been shown in [28] that only three "polynomial conservation laws" exist, except in the case p = 1 - the KdV equation - or p = 2. Thus it is quite possible that no more than three conservation laws can be found for (7). For completeness we state the three conservation laws we could generalize, and we use one of these to show that a bound 8 can be obtained, given a bound on }|u (x,t)|J. The conservation laws are: u,. = -(e(u) + s2u ) t xx x (%u2) = -(<Ku) + <52(uu - %u2)) t X X * x"x ' 2 2 o (tf(u) -6 u ) = -U(u)- + 6^'(u)u - 2g(u)u2) X t T X X B X 4 2 -26 (u u - kxi ) 1 x xxx d xx7 Jx 9'(u) = g(u) , <j,'(u) = ug(u) , f (u) = 2g(u) and £« (u) = i|/*(u)g(u) 34. Using the periodicity these give the conserved quantities ( u(x,t)dx = a , ' o o 2. L 2 |u|f = Jo u (x,t)dx = , (28) / (Y(u) - 62u 2) = a X for a l l t. We require the following Sobolev lemma. Lemma 2.10: Let 0 ^ m < n be integers and let E > 0 be a given constant. Then there is a constant c such that for a l l functions y, sufficiently differentiable on [0,L], ,m 2 max I i l l ; „ m' XE[0,L] 9X' max I 1 $ E xe[0,L] 3xm 3x n + c 3 y 3x A bound on max |u | now follows. For assume that x u (x,t) < B 1 XX o for a l l t. Then, using Lemma 2.10 and (28), 35, max |u | ^ ellu || + c||u| xe[0,L] X o 1 = B. We now give the uniqueness result. Theorem 2.11: Let u and v be smooth periodic solutions to where g has a continuous derivative. Suppose for a l l t and some B max |u | < B . If u(x,0) =v(x,Q) then for a l l t , u(x,t) = v(x,t) X E [ 0 , L ] X Proof: Let w = u - v . Then w satisfies 2 w = u - v = -g(u)u + g(v)v - 6 w t t t X X . xxx = -g(u)u + g(v) uv - g(v)u + g(v)v - «S2w 2 = "(g(u) - g(v))u - g(v)w - 6 w X x xxx Applying the Mean Value Theorem to g , we have for some £ = £(x,t) % J w dx = / ww dx a t Q J t o L = -/ (wg'(C)wu + g(v)ww + 6 ww )dx o Consider the terms separately and integrate by parts. Remember w has period L L 2 L L f 6 ww dx = 6 {ww I - J ww dx} J xxx xx1 ; x xx o o o 9 9 L = 6Z {0 - %w I } = 0 . o L 2 L L 2 / g(v)ww dx = %g(v)w | - %/ g'(v)v w dx x x o o L X 2 = -¥/ g' (v)v w dx o Thus we have 9 I I , J12 • ^ 2 % ^Ilwll =~/ w (g'(5)ux - %g'(v)v )dx o x L 9 4 / v A i g ' C O u J + %|g'<v)v |)dx o L < B./ w/(|g'a)[ + %|g'(v)|)dx . o But, applying Lemma 2.10 again we have max" - | u| < e ||u || + c][.u[.J xe[0,L] X < cL2 .max ' |u | +--ca* X E [ 0 , L ] X " L 4 eSL% + ca* so u and v must be bounded. Clearly then £ is also. Since continuous, we have |g'(£)| , |g'(v)| bounded, say by • 37, Then % |[w||2 -4 y BB /L w2dx at " " " 2 ""2 o | 8 B 2 |MI 2 Since |[w(x,0)j| =0 we must therefore have |[w(x,t)|| = 0 for a l l t and hence w = 0 . q.e.d. The proof of the above theorem applies with virtually no change to the solutions of (7) defined on the real line (-00 , 0 0) . In this case it can be shown that smooth solutions which approach a constant sufficiently rapidly as |x| -»- 00 are uniquely determined by their i n i t i a l data (again, assuming the existence of the bound B ) . 38. I l l Numerical Results §1 Numerical Simulation of the KdV Equation In this chapter we present numerical results (in picture form) which demonstrate that the solitons observed by Zabusky and Kruskal also occur in approximate solutions to the BKdV equation and that they provide an explanation for the recurrence observed in approximate solutions to Johnson's bounded sound speed wave equation (II.1) (which we refer to as equation II.IB). In this section we firs t discuss some necessary prelimi-naries and then, for comparison, we look at the solitons observed by ZK in solutions to the KdV equation u + uu + 62u =0 . (1) t X xxx To begin we replace the fixed boundary conditions in (1.2) by periodic boundary conditions. In particular, we require u(x,t) = u(x+2,t) for a l l t , with solutions of (1) to be computed on [0,2] . In this way the soliton behavior will become more evident. We must ask, however, i f periodic boundary conditions are suitable for studying soliton behavior, for recall that a soliton is defined on ( - c o } c o ) and is asymptotically constant as |x| -*• 00 , and hence certainly not periodic (unless i t is identically equal to s^ . We are not interested in this trivial case.) In fact, because of the rather special behavior of the solitons, they may be studied in this way. To explain, suppose there is a solution to (1) in which more than one soliton 39. is present. Then the peak of each soliton s = s(x) will move with constant velocity c and constant amplitude (except where i t is nearly coincident with another peak), just as i f the solution to (1) were u(x,t) = s(x-ct). Thus solitons act independently. Now the effect of periodic boundary conditions is to create, from a single wave pulse on [0,2] , an infinite family of well separated wave pulses on (-co,co) . Since the pulses are well separated, the original peak on [0,2] behaves as i f i t alone were present. We illustrate with an example. Using (I.10),let u(x,0) = s(x) = A sech2 ~ - for x E [0,2] , (2) 2 where A = 1, 62 =.000484 , A2 = The approximate solution to (1) computed with i n i t i a l conditions (2) and using the periodic boundary conditions is exhibited in Figure 6. (In the figure, the computed solution to (1) is displayed at the sequence of times t = 0,.2,.4,... . Each horizontal line represents the interval [0,2] of the x-axis. This is the manner in which a l l computations involving the KdV or BKdV equations will be displayed.) It is easily seen that the peak lying on [0,2] moves with constant velocity -^ ; this agrees with its predicted velocity, as a soliton, of -j (see I . l l ) . Next we must choose a difference scheme to compute approximate solutions to (1). (That i s , we now specify the scheme used to obtain Figure 6.) Define xj = Jh j = 0,1,...,J t =nk n = 0 ,1,2 ,. .. n \ (3) 4.0. with Jh = 2 . First consider the scheme used by Zabusky [46] (see his Bibliography, #6): n+1 n k.^ , n,2 , n^ n. .2, ^ n ... v. = v. - — (D (v.) + v.D v.) - 6 kD.D D v. . (4) 3 j 3 o y 2 o 3 + o - j Here, as throughout this chapter, v denotes a solution to a difference scheme. (The solution to the differential equation is denoted by u .) v is defined on the lattice (3) and v? = v(x.,t ) . This scheme is 3 3 n explicit - that i s , vj " * ^ c a n ^e c o mPu ted directly from the values v), v] i and v^ „• - and hence is simple to apply. Unfortunately, to j j ± l j+2 ensure that (4) is stable (we will define exactly what we mean by stability 3 in Chapter IV), k must be 0(h ) (See Kreiss [23]). Since values of h like .01 and smaller will be used frequently (in order to resolve individual wave pulses), this imposes much too severe a restriction on the size of k (since we must compute for long times to observe the recurrence). Consequently we do not use Zabusky's scheme. (We remark that Zabusky does not discuss stability of his scheme; nor does he indicate the size of k required for his computations.) It is of interest to note Zabusky's choice of difference quo-tient to approximate the term uu . In particular, by writing 1 2 uu = — ((u ) + uu ) x 3 x x it is clear that the natural second order centered difference approxima-tion 1/^ r n. 2 , n„ ns tr.\ 3(Do(v.) +v.Dov.) (5) 41. gives the second term on the right in Zabusky's scheme (4). It is not clear why Zabusky chose (5), rather than the obvious approximation vUD v11 , but i t is interesting because we will see below that exactly. J o 3 the same choice was made by Sjoberg. We now examine three similar difference schemes. (The motivation for these schemes appears in Chapter IV, §3). Using the splitting (5), we obtain Sjoberg's scheme [39]: On the other hand, one frequently writes equations such as (1) in their conservation form: Using the "natural" form (1) of the KdV equation, we get n-1 (7) 1 2 2 u + — (u ) + 6 u = t 2 x xxx 0 This leads naturally to the scheme n-1 (8) Each of these is implicit; that i s , a system of linear equations must be solved in order to obtain the values v\+^~ . Now, for his computations, 42. Sjoberg used the rather special scheme (6). This suggests that, in some way, the scheme (6) is to be preferred to (7) or (8). This does not seenv to be the case. We illustrate with an example. An approximate solution to (1) was computed using each of the schemes (6), (7) and (8) (the compu-tations were otherwise identical); these are displayed in Figures (7), (8) and (9), resp. It is evident that the results are indistinguishable. 2 Other comparisons (not displayed), using several values of h,k and 6 also indicate no reason for preferring (6). We therefore compute a l l further approximate solutions to (1) using (7), since (it is easily shown that) i t requires the least arithmetic. We remark that (7) is also more readily adapted to equations (11.10) and (11.11). to solve (1). Our object is to compare (approximate) solutions to (1) to computed solutions to (II. 1FPU). For a meaningful comparison, the i n i t i a l data must correspond. Let the i n i t i a l data for (II.1FPU) be specified by As a final preliminary we determine the i n i t i a l data required y(x,0) = a sin T T X x e [0,1] (9) From the derivation of u (Chapter II, §1), y (x,o) u(x,0) = |(-yt(x,0) + J 0 Using the definition of $ (II.2), and (9) > 43. u(x,0) = — cos T T X + 0(e) = C(x) + 0(e) . We ignore the 0(e) term and use as starting data for (7) vj H c (V 5 ?j j=0,l,...,J . (10) Since (7) is a 3-level scheme, we require starting data at t = k as well. Using Taylor's Theorem and truncating, u(x,k) = u(x,0) + ku (x,0) + 0(k2) = c ~ k(cc' + « V ) , from which we get vj " cj " k ( cj?j + &2ty j = 0,1,...,J . (11) We now discuss the relationship between ZK's solitons and the energy recurrence observed in solutions to (II.1FPU). We demonstrate the correspondence using three examples. Each example consists of (i) an approximate solution to (II.1FPU) obtained using particular values - of e, h (the spatial mesh width) and a (see (9)); and (ii) an approximate solution to (1) obtained using 7 , h2 6. (= i s e e- Chanter II, §1) and in i t i a l data ( see (10) and (11) ) corresponding to ( i ) . -These examples are exhibited as figures 10 and 8 , ; Figures 11 and 12 and Figures 5 and 13. The exact values of e, h, a and 2 6 used appear with the figures. (Three examples are included to show the variation in the number of solitons formed, and to support observation 44. 3 ° , ( i ) , below.) ZK claim that soliton interaction explains the energy recurrence. The following observations support this claim. 1° The energy recurrence: Each of Figures 10, 11 and 5 (abbreviated as F 10, 11, 5) exhibits a recurrence of the energy into the first mode at some time t ; that i s , there is a recurrence of the i n i t i a l , R state of the system. In each of F 8, 12, 13, the i n i t i a l state is also reproduced at some time t . (Except for a phase shift). See for example t = 6.0 in Figure 8 where- the original cosine wave is clearly reproduced. 2° The recurrence mechanism: In F 8, 12, 13 we can see how the recurrence takes place. Take Figure 8 for example. At (a) there is 2 the i n i t i a l data cos TTX . From (a) to (b) , the 6 term of (1) has no effect and the classical steepening takes place where the slope is 2 negative. At (c), 6 u becomes significant and prevents formation of a discontinuity. From there to (d) several wave pulses are created. Then, after they are completely formed at (d), each pulse moves with constant velocity and amplitude , losing its identity only when i t is nearly coincident with another pulse. At (e) the phase relationship of the pulses is reversed from (d). The pulses then interact and reproduce the original waveform. ZK confirm (see [46]) that each of these wave pulses is a soliton (as defined in Chapter II, §2. See as well 1.9, I.10). 3° The detailed modal histories: The interaction of the solitons not only causes a recurrence of the i n i t i a l state, but also explains details of higher mode energy histories. We note the following points. 45. (i) In F 10, 11, 5 the energy remains concentrated in a finite number of modes. In F 8, 12, 13 a finite number N of solitons are generated. In each case = . That i s , the number of solitons appearing in an approximate solution to (1) coincides with the number of energy modes which are excited in the corresponding approximate solutions to (II.1FPU). (ii) In F 10, 11, 5 the recurrence interval [0,t ] is divided R into i equal subintervals by the maxima of the energy in the i * *1 mode (i = 2,...,N ) . For example, at t = ~-t , the energy is predominantly. X Z R in the 2nc* mode. Correspondingly, in F 8, 12, 13, at t = y t the solitons have formed into two distinct groups and the approximate solution to (1) appears as the 2nc* harmonic of the i n i t i a l data; that i s , for example, as cos 2TTX in Figure 8. A similar correspondence holds for the higher modes. On the basis of these observations we conclude, as did ZK, that there can be l i t t l e doubt that the soliton phenomenon does provide the mechanism by which the recurrence takes place. We finish this section with a brief discussion of a discrepancy not mentioned by ZK: the times t and t are not in exact argument. In the three examples at hand the magnitude of the error is about 19%, 22% and 15% (Recall that the time scale in F 8, 12, 13 is E 2 reduced by a factor of -r- . By t* we mean — t ) . We cannot explain 2 ^ r e r this discrepancy, but we indicate at least that two of the most immediate possible explanations are not valid. 46. . Recall that in Chapter II we motivated the derivation of a modified form of the KdV equation, namely (I.10): e 2 2 u + (u - ^ u )u + 6 u =0 . t 2 x xxx e 2 One might suggest that inclusion of the previously neglected term — u u X could lead to solutions which s t i l l exhibit the soliton behavior, but in which the time t * agrees with t„ . An example shows this is not r • R the case. Consider Figure 14, which displays an approximate solution to 2 (11.10) obtained using the same 6 and i n i t i a l data as in Figure 8. (The approximate solution was obtained using a difference scheme formed from (7) by using the difference quotient (y\ - 2"(Vj)2)D0Vj t C > s^m u-'-a t e the nonlinear term.) Clearly Figure 14 displays the same soliton behavior as Figure 8; the recurrence time t is however unchanged. To see why this term has no effect, and to motivate our second suggestion, we look statistically at the relative size of the terms I u2u I and 16 2u I : 1 2 x1 1 xxx1 2 (i) From Table I we see that <5 u attains much larger xxx e 2 values than — u ux and, as well, i s , on the average, an order of magnitude larger; E 2 (ii) From Table II we see that |— u u^| is seldom larger 2 than 16 u | , ind that this occurs only when both terms are quite small. X X X 2 We see then that 6 u becomes quite large and significantly X X X e 2 affects the solution; — u u remains comparatively small and does not. Z. X Furthermore, since i t is dominant, one might suggest that an adjustment 2 to the 6 term might correct the recurrence time discrepancy. Such 47. a modification appears in equation (11.11): 2 e 2 u + u u + 6 ( u - - r u ) = 0 . t x 4 xxx 2 However in Figure 15 5 and the i n i t i a l data again correspond to Figure 8, but the same soliton behavior is observed with no change in t r 48. TABLE I FIGURE 8 | f u 2 u I I 6 2 u [ '2 x11 xxx' FIGURE 12 |fn2u||62u | 12 x1' xxx1 FIGURE 13 [fu2u||52u | 1 2 x11 xxx1 MINIMUM MAXIMUM MEAN STANDARD DEVIATION 0. 5.21 69.6 .16 4.3 .52 8.5 SKEWNESS 370.4 254.5 0. 0. 0. 0. 2.16 88.1 11.5 127.3 .09 6.8 .28 7.0 .27 12.3 .89 12.9 132.7 98.0 558.3 367.4' TABLE I: SIZE COMPARISON OF THE TWO TERMS ^u2u AND <S2u 2 x xxx The statistics were obtained from values for each of the two terms at each grid point (x.,t ) (computed using difference quotient approximations - J n from the solutions to (7)) where 0 £ i < J and t is one of the times J n appearing in the appropriate figure (8, 12 or 13). 49. TABLE II FIGURE 8 | - |u 2 u I I 6 2u | 12 x11 xxx1 FIGURE 12 | - |u 2 u I I 6 2u I '2 x'1 xxx1 FIGURE 13 | § u 2 u I I 6 2u I '2 x11 xxx1 PERCENTAGE MAXIMUM MEAN STANDARD DEVIATION 3.7% .7% .37 .03 .03 ,11 .02 .02 .063 .015 .016 .038 .008 .009 4.9% .81 .47 .03 .01 .04 .03 TABLE II: SIZE COMPARISON OF THE TWO TERMS fu u AND 6 u IN THE 2 x xxx CASE THAT THE FIRST TERM IS LARGER (IN ABSOLUTE VALUE). The statistics were computed as in Table I, but values were only used i f 2 2 |yu u I > 16 u | . The percentage figure indicates how often the first £. X X X X term was larger. 50. §2 Solitons and the BKdV Evolution Equation .The BKdV equation is given by u + g(u)u + 62u = 0 , (12) t x xxx where g(u) = 1 1 9 , (1+u f We now demonstrate that the solitons observed in §1 also appear in (approximate) solutions to (12) and we relate them to energy recurrences observed in approximate solutions to (II.IB). The preliminaries follow §1. We consider periodic boundary conditions on [0,2] and we compute approximate solutions on the grid (3) We use a difference scheme motivated in the same way as those in §1: (I+k<$2D 1) D )v?+ 1 = 2kg(v?)D vn + (I-k62D D D )v" 1 . (13) + o- j j o j + o - j The difference quotient discretization of the nonlinear term is the natural second order one. The i n i t i a l data for the scheme (13) is given at time t = 0 by (10) and at time t = k from Taylor's Theorem (obtained in the same way as (11)): We now exhibit four examples which demonstrate the energy recurrence - soliton relationship. (The cost of each computation limits the number of examples given.) Each example consists of a pair of 51, corresponding figures, as in §1. These are Figures 16 and 18, 17 and 19, 20 and 22, and 21 and 23. The first two examples illustrate the soliton behavior when i t is particularly simple; only 2 (positive; see 2° below) solitons are present. In fact Figures 18 and 19 are nearly identical but both are included to show that the energy recurrence - soliton relationship holds with different values of e (see Figures 16 and 17). The last two pairs show the interaction with more solitons present. The following observations support the claim that solitons occur and explain the energy recurrence. 1° The energy recurrence: F 16, 17, 20, 21 clearly show a recurrence (at time t ) of the energy into the first mode; that i s , a return to the i n i t i a l state of the system. The computed solutions to the BKdV equation (12) (F 18, 19, 22, 23) also exhibit a return to the in i t i a l state. In F 18, 19 this is clear (at t = 28, 27.5 resp.). In F 22, 23 the recurrence is less obvious, although at t = 26, 15.5 resp., the solution appears to be primarily composed of the first mode. To make this more precise we introduce Figures 24 and 25. These were obtained by decomposing the approximate solution appearing at each time step in Figures 22, 23 resp., into a (finite) Fourier series and then plotting each Fourier coefficient as a function of time. It is clear from these figures that the first mode is again predominant at some time t^ > 0 (a 26, 15.5 resp.) 2° The wave pulses: These are evident in F 18, 19, 22, 23. As in §1, the i n i t i a l data appears at t = 0 . Then the classical 2 steepening takes place, the <5 term being negligible. In constrast to §1, 52. however, this takes place in regions of both positive and negative slope. 2 Then, before discontinuities can form, 6 u becomes important and xxx wave pulses are created. Notice that both positive and negative pulses, appear. They are equal in number and have corresponding amplitudes. Careful observation of these pulses verifies that each pulse, once formed, moves with constant velocity (except of course during interaction with other pulses). The recurrence then occurs when the phase relationship of the pulses is such that they interact to reform the i n i t i a l waveform. We also note: (i) Consider Figure 18. We claim that two positive solitons are present, as indicated by the arrows at t = 6 . For, between t = 20 and t = 28 i t is the mutual interaction of the two positive (and the two negative) pulses which causes the recurrence. By inference then we conclude the smaller pulse is a travelling wave pulse, even though i t is so obscurred by the other pulses that its velocity cannot be determined. Similarly in Figure 22 there are three wave pulses. (ii) The wave pulses move slowly, in comparison to those observed in §1 (as can be determined by direct measurement). This is reflected in the magnitude of the time t . For example, in Figure 8, with e = .068 and 8 solitons present, t £ 9.7. But in Figure 18, with e -• .1 , there are only 2 positive solitons. Yet t = 28. This agrees with our prediction following Theorem 2.1 about the bounded soliton speeds. 3° The solitons: Recall that a soliton is defined as a travelling wave solution u(x,t) = s(x-ct) to (12). We have noted that the observed wave pulses have the characteristics of travelling wave 53. solutions: they move with constant velocity and amplitude. In the following section we establish (numerically) that the observed wave pulses are solitons, just as ZK did for the KdV equation. 4° The higher mode energy: Certain features of the higher mode energy histories can be determined from observation of the soliton interactions: (i) F 16, 17, 20, 21 show that no energy is fed into the even modes. Examining F 18, 19, 22, 23 we see that solitons are generated in pairs, one positive and one negative, each with the same velocity. (Recall from Chapter II, §3 that u(x,t) = s(x-ct) and u(x,t) = -s(x-ct) are solitons with the same speed and whose amplitudes are equal in magnitude but of opposite sign.) This results in no even modes being present; see for example Figure 24 and 25. (ii) Only a finite number of the odd modes ever attain significant energy (see F 16, 17, 20, 21). As in §1, there is a correspondence with the finite number of solitons generated in the (approximate) solutions to (12). In this case, however, the number of odd modes excited agrees with the number of pairs (i.e. a positive soliton and the equal negative soliton) of solitons appearing. On the basis of these observations we conclude that ZK's explanation can be extended to Johnson's wave equation (II.IB). We complete this section by remarking that the explanation, although basically correct, i s , as in §1, not completely accurate. We do not present any explanation for these inaccuracies, first because we have already seen in §1 that the obvious explanations are not correct, and secondly because the modifications 54. to the BKdV equation (12) corresponding to (11.10) and (11.11) lead to difference equations whose solutions are prohibitively expensive to compute. The failings of our procedure are: (i) As in §1, the times t * and t do not agree. However ' r R ° the agreement is as good as in §1, being about 11%, 5%, 7% and 25% in the four examples at hand. (ii) The solitons do not smoothly reform the i n i t i a l waveform, particularly in Figures 22 and 23. This is caused by the phase relation-ship of the solitons not being correct at that time at which they interacted to move the energy back into the first mode. To a lesser degree this occured in §1; see Figure 13 where at the recurrence time t =9.7 there is a small amount of the higher modes superimposed on the solution which consists primarily of mode 1. §3 Computation of Individual Solitons We have two goals in this section. First, we seek a satisfactory method of numerically solving the boundary value problem defining the solitons for (12); for convenience we repeat (11.14) and (11.12): i-2 i n f , \ . . 1 • 6 s = (c-l)s + •- s with boundary conditions given by lim s(x) = s^ , a constant |x(-*=° (14a) (14b) lim s'(x) = | x|-*>° lim s"(x) = 0 |x|-x» (14c,d) 55. This enables us to obtain an approximation to a particular soliton 2 (described by {' , its speed c and its background s^ ) . Secondly we wish to show that the wave pulses observed in §2 are in fact solitons; we do this by comparing them to particular computed soliton approximations. Our basic approach is to treat problem (14) as an i n i t i a l value problem, or IVP, (in much the same way as in Chapter II, §3); however we first indicate why we feel standard numerical techniques for treating boundary value problems are unsuitable. The principal difficulty with (14) is that the boundary conditions are specified at infinity. We examine three common procedures. The first is the method of finite differences, which involves discretizing (-00 , 00) using a finite number N of points, and then replacing derivatives by difference quotients at the grid points (just as we have done for the partial differential equations). Unfortunately this leads to a nonlinear system of equations;since the system will have N equations and N must be large (to achieve the resolution necessary), we do not choose this method. Secondly, consider the so called "shooting" methods (see Keller [22]). Here a boundary value problem is handled by solving a sequence of IVP's whose solutions (hopefully) converge to the solution of the boundary value problem. Since we will shortly obtain an approximate solution to (14) by solving a single IVP, there is l i t t l e reason to pursue shooting methods. Finally there is Galerkin's method (see Kantorovich and Krylov [19], p.262). It attempts to approximate the solution s by a linear combination of preselected functions. (Since an approximation to s is thereby obtained which is defined at every point on (-00 , 00) , and not just at some finite collection of 56. points, i t is often called a global method.) There are two principal reasons why we discard this method. First, i t is difficult to choose suitable approximating functions which match the boundary conditions. Secondly, recall that (14) always has the trivial solution s(x) = s^ There is a tendency of the method (in this case) to find this solution. We now turn to the problem of solving (14) as an IVP. For numerical purposes we assume existence of the desired solution (specified 2 by 6 , c and s ). We restrict our attention to positive solitons. CO Since solutions to (14) are translation invariant, we can assume the solution s has a peak at 0 . Then s'(0) = 0 . Let S Q 5 s(0) Using the symmetry of the solution (Theorem 2.7), we only compute the solution for x > 0 Recall that (14) may be integrated to second order. From equation (11.21), 62s" = (c-l)s + sinh"1s + K (15) where K = (l-c)s - sinh ^s 2 Thus, given 6 , c, s^ and Sq , we have a completely specified IVP with s(0) = s , s'(0) = 0 . Now, we have seen in §3 of Chapter II how o 2 2 Sq may be determined from 6 , c and . However, given <5 and c , equation (11.22) equally well determines s^ as a function of SQ In particular, s is a root of 57, i(t) = ~- t 2 + A+t 2 + (l-c)s t - s sinh- 1t + M where M = s2 + s sinh - Al+s2 2 o o o o (In most cases s^ is uniquely determined as the only root of G less than Sq . However i t may happen that 8 has two such roots, in which case the one nearer s should be chosen, or none, in which case there is no o 2 soliton for the given c and s^ .) So, given only 6 , c and S Q , our IVP is completely specified. This is important, since i t is precisely 2 6 , c and s^ (but not s^ ) which may be determined for the travelling waves observed in §2. We remark that by solving (14) as an IVP we do somewhat restrict cur computed solutions. For (once is known) we have seen (Theorem 2.7) that only when Sq is the particular root v-^ of <J> (<J> depends on s^ ; see 11.23) does the solution of (15) asymptotically approach s^ . Now, in every algorithm for computing solutions to an IVP, errors are incurred (by truncation and roundoff). The effect of these errors may be considered as a perturbation of the i n i t i a l data. Thus (to an IVP algorithm) i t appears that SQ ^ v^ > therefore that the solution sought is not asymptotically equal to s^ . We conclude then that we should not expect any IVP algorithm to give approximate solutions to (15) which become arbitrarily close to s^ for values of x becoming arbi-trarily large. Fortunately, however, as will be seen below, this limi-tation does not prevent us from obtaining satisfactory approximate solutions to (15). 58. Recall our first objective. We seek a computed approximation to a travelling wave solution to (12). Now suppose we were to compute (using the difference scheme (13) of §2) an approximate solution to the BKdV equation (12) with i n i t i a l data taken as the true solution s to (15). The result would be a travelling wave of course (as was found -see Figure 6 - for the KdV equation (1)). It would therefore be reasonable to consider our approximate solution to (15) as adequate i f the same behavior was observed when the i n i t i a l data was obtained from the approxi-mation to s . An example shows that, in this respect, our IVP approach provides satisfactory solutions (15). See Figure 26, where a travelling wave is clearly obtained by using as i n i t i a l data for the scheme (13) 2 the solution to (15) computed with 6 = .0003, c = .22 and Sq = 3.03. The wave travels with speed .21, in close agreement with its theoretical value. We have thus achieved our first objective. We conclude by looking at the second. We illustrate with three examples. Consider: (i) the largest (positive) pulse in Figure 21. (ii) the second largest (positive) pulse in Figure 21. (i i i ) the large pulse in Figure 16. (We are unable to look at any of the small pulses since, as we mentioned in point 2° of §2, they are so obscured by the other pulses that their velocity may not be determined.) 1° use the displayed solution of (12) (Figure 21 or 16) to 2 measure the velocity c and amplitude s of each pulse (6 is given). 59, 2 2° Use 5 , c and Sq (from 1°) to compute an approximate solution to (15) (on some finite interval). 3° Pick a particular time at which the pulse in 1° is least affected by its neighbours. (This seems to be just after a l l the pulses are ini t i a l l y formed.) 4° Compare the pulse from 3° with the soliton from 2° by simply plotting them on a single set of axes. See Figures 27, 28 and 29 for the comparisons for the pulses ( i ) , (ii) and ( i i i ) , resp. The comparison is only of interest near the peaks, for two reasons. First, the asymptotic part of the observed pulse (in 3°) is completely obscured by its neighbours. Secondly, the computed solutions to (15) (see 2°) cannot be determined accurately in the asympto-tic region, as we have indicated earlier. Clearly, however, each pulse, near its peak, bears a close relationship to the soliton of the same velocity and amplitude. We therefore conclude, with l i t t l e doubt, that the observed pulses are solitons; that i s , special travelling wave solutions to (12). 60. IV. Stability §1 Frozen Coefficient Stability for the Wave Equations We now consider the stability of the difference scheme used to compute approximate solutions to the general wave equation (II.1). We let v denote the approximate solution defined on the grid given in §2 of Chapter I and we let v1} = v(x.,t ) . The difference scheme is given 3 3 n by vn + 1 = 2vn - v ^1 + k2{(l + cg(D vn) )D D vn + eD2D2vn} , (1) which explicitly defines v at time from its values at times t and t , , given the starting values n n-1 v5" f<V 1 v2 v. = f(xj) + ^ [ ( i + £ g(f,(xj)))f"(xj) + ef-'Vj)] for j = 0,1,...,J and the boundary conditions v^ = v^ = 0 for a l l n (See Cliapter I, §2). We assume k is given as some function of h such that k ->• 0 as h -»- 0 We now define what we mean by stability of (1). In order to do this we require the following. As written, (1) defines v only at the grid points (x^ . ,t^) . However i f we write (1) as v(x,t+k) = 2v(x,t) - v(x,t~k) +k2{[!+ e g f D (v(x,t))] o x D D v(x,t) + BD2D2v(x,t)} , then i t is clear that,at each time level t , we may consider v to n be a function defined on the entire interval [0,1] (provided of course we define the i n i t i a l data on the whole interval).. We let vn(x) = v(x,tn) for x e [0,1] . We assume v" is square integrable on [0,1] and use |vnlj to denote the 1, = I^fOjl] norm given by 2 1 2 vn|| = / (vn(x))Zdx . Then we have the Definition: (1) is I^-stable i f there exist constants K and y , independent of h and k , such that ||vn|I < Ke^11 | v ° l f o r a l l n . (2) By stable we will mean I^-stable. We say (1) is strongly stable i f the inequality (2) holds for some K with y = 0 , and conditionally stable i f i t is stable for some restricted class of functions defining k We would like to obtain a practical criterion for the conditional stability of (1); that i s , an easily verified condition on h and k which will ensure stability of (1). Since (1) is nonlinear we are unable to present such a condition which is globally valid; rather, we will follow the commonly accepted practice of examining the corresponding constant 62, (or "frozen") coefficient problem to obtain a condition for stability which will be close, locally, to the true condition for stability of the nonlinear scheme. We use the technique introduced by von Neumann to study the stability of constant coefficient difference schemes. Let a be a constant and consider n+1 „ n n-1 . 2, 2^ n ,,^ 2^ 2 n. v. = 2v . - v. + k {a D.D v. + 3D.D v.} , 3 3 3 + ~ 3 + ~ 3 (3) which is the constant coefficient scheme corresponding to (1). As.for scheme (1), we say (3) is stable i f the inequality (2) can be satisfied. To investigate a 3 level (or 2-step) scheme like (3), we first write i t as a 1-step scheme by letting wj 2 vj ^ ' f ) v w v. J n+1 3 2 2 2 2 2 2 + k a D+D_ + k 3D+D_ -1 0 In v w I )3 (4) In the same way as we defined vn(x) on the entire interval [0,1] above, we can now define a vector valued function v_ (x) = v (x) w (x) on [0,1] . We now use some basic concepts from the theory of Banach spaces to study (4). Let B be the Banach space of a l l vector valued functions f v(x) = v(x) w(x) with V,WE L2[0,l] and with I|V«2H /0 |v|2dx • 2 2 ,(v (x) + w (x))dx v||2+||w||2 (5) €3. We can then consider that (4) defines an operator P on B and write (4) as 1 1 + 1 T) 1 1 fC\ y_ = P v (6) We define the norm of P as the induced operator norm: Pvf P|| = sup {; v v i 0, v e B} To study (6) we introduce the Fourier series of a function v e B : ijx v(x) = Y C. e 3 3 1 where C. = / v(x)e 1;)Xdx J o is a 2-vector. The Fourier coefficients C. of v define a sequence 3 - r C = {C^ .} . For any sequence of 2-vectors, C = {C^ } , let C_. = a. 3 b. and let |[C.|L be the Euclidean vector norm. Then the collection of 3 a l l sequences C with llcll2 = Illc ||2 3 = Y(a2 + b2) < oo • 3 3 3 forms a Banach space B with the norm [|*[] . Notice that by the Parseval equality ||v|| = ||c|l i f C = {C^ .} is the sequence of Fourier coefficients of v 64. Since P is an operator defined on B , i t induces an operator P on B which can be defined by its action on each Fourier coefficient: {PC}j = G (h)C , where G^(h) , called the amplification matrix, is a 2 x 2 matrix given by G = G^(h) r2-2a± + 6B1 + (a1-4g1)-2 cos £ + g *2 cos 2 £. -1 1 (7) with an a2k2 gk2 2 > 0 » 3-L - — A 5 = jh . (G is also known as a h h "symbol" of P ; see Lax and .Nirenberg [31].) Thus i f yj*" has Fourier coefficients C. , vn~*"^will have Fourier coefficients 3 -(G (h))" C For stability of (3) we require that none of the Fourier coef-ficients of be amplified. Thus for a l l j and a l l h , we insist that G..(h)||2 < 1 , (8) where || G|2 denotes the Euclidean or spectral norm of the matrix G To see that this implies stability of (3) according to our definition, choose v E B with Fourier coefficients C = {C.} . Then 65. I I P C I I ^ I I I G . C . H 2 ^ ! ^ ! ! 2 H e . | f < I l l c J 2 = H e | | 2 . Thus ||p|i < 1 (where again we mean the induced operator norm on P ) , and i t follows from the Parseval equality that ||P|| < 1 . Then, since the induced operator norm on P is compatible with the norm on B , from (6), l l y l 2 -flpv^1!2 < ( I I P M I Z " 1 ! ! ) 2 n-li|2 Using (5), llvnP + ||wn||2 < llv11"1!2 + |w n"Y = Hw n|| 2 + | | v n - 2 | | 2 or, I v l 2 < | | v n - 2 H 2 « . . . v"*"|2 , i f n odd , Oil 2 ., v II , i f n even. Letting K = max {1, V~ } , we see that the uniform boundedness condition llv°|| (8) implies that (3) will be strongly L^-stable. 66. Now let us examine the amplification matrices G (h) . A necessary condition for (8) to hold is that for a l l £ , both eigenvalues of G l i e on or within the unit circle. This follows easily, since for any matrix norm 1| 'H , the spectral radius p (G) of G , defined as the maximum (in modulus) eigenvalue of G , satisfies (p(G))n< «Gn|! < || G« n for any integer n > 0 , From (7) the eigenvalues of G satisfy A2 - 2{1 - a ; L( l - cos O + g (3 - 4 cos g + cos 2 £)}A + 1 = 0 2 It is easily shown that the quadratic \ - 2bA + 1 has both roots bounded by 1 in magnitude i f and only i f b is real and fb| < 1 Thus (8) can hold only i f 0 < a^l - cos O - g (3 - 4 cos £ + cos 2 O < 2 2 or, 0 < a^l - cos g) - 2g1(l - cos £) < 2 2 for a l l £ . Set x = 1 - cos £ and let S(T) = a^ x - 2g^ x Then, for stability, we must have 0 ^ S(T) < 2 for 0 4- x < 2 Notice s(0) =0 . B y differentiating i t follows that s has a single al c ritical point at c = . (We assume g , and hence g^ , is nonzero since the case g = 0 is well known. See, for example, Isaacson and Keller [13],, Chapter 9, §3.) Because s is a quadratic, we can distinguish the following cases: 67. (i) B < 0 : here c { [0,2] , so (8) holds provided s(2) < 2 . This will be true i f 2 46 \ < 1 . (9) h (ii) B > 0 and a2 < ^ f" : h e r e c e [°>2) 3X1(1 (8) h o l d s i f h < a2k < 4/B • (10) h 2 86 ( i i i ) B > 0 and a £ 2 : again c £ [0,2) and (8) holds i f h (9) is satisfied. 2 We now have the desired stability criterion. That i s , given a and B w e can, knowing h , find k to ensure stability (of the constant coefficient scheme). However we would like to discuss these conditions in more detail. First consider very small B (small that is with respect to h and k ). Then we have case (i) or case ( i i i ) and the stability condition is given by ( 9 ) . For these small g's , (9) is just a slight alteration of the well known wave equation (i.e. g = 0) stability condi-ctk tion — < 1 , as we would hope. For g > 0 we have a slightly relaxed condition (that i s , for a given h we may choose a somewhat larger k ) , but for g < 0 we must in fact have < p < 1 for some p 2 Secondly, suppose a and g are given and h has been chosen 2 86 so that a < 2 > s° we have case ( i i ) . Then the left hand inequality h 83 ' 2 of (10) gives the stability condition -y < 2a . This has the effect h of placing a lower bound on the value of h we may use: in this case h must satisfy 68, i h 2 < 2M 2 a 46 for stability. In particular, i f is not small, no value of h a which is small enough for practical purposes may be permissible. Finally, when 3 5. 1 , (10) places a strong restriction on k . For although the right hand inequality in (10) is independent of h , we have, assuming h satisfies the left hand inequality, . 4/3~ 43 ,2 k < — N < " 2 s < h a a 2 Thus a large a (relative to 3 ) demands a small k (relative to h) for stability. Finally we would like to discuss the situation we considered in Chapter III, where 3 = — > 0 . Then the stability conditions become: (ii) i f a2 < -j , (8) holds i f J 7 < T/3 and a2 » \ \ (11) n J 2 2 ( i i i ) i f a > j , (8) holds i f k2 2 1 \ (a - f> < 1 . (12) h 2 To determine the values of a of interest, we must examine the nonlinear scheme (1). We should require the constant coefficient scheme (3) to be 2 stable for those values of a attained by the nonlinear term in (1). 69. (That i s , we require (1) to be locally stable at a l l points (x,t) .) When g(u) = 1 - — ^ •-. , the nonlinearity is bounded by 1 and 1 + e , (l+O* 2 so for the scheme (3) we require stability for • 1 < a < 1 + e . When g(u) = u , however, the nonlinear term 1 + ey^ is unbounded and the situation is less clear. In particular, i f for some point (x,t) , 1 + eyx(x,t) < 0 , (1) becomes an illposed elliptic problem and a stable scheme is not possible. Moreover, we can see from (11) that a 2 1 stable constant coefficient scheme is not possible for a < 3" • Thus, to obtain a stability criterion, we assume that 1 + ey (x,t) > — X j for a l l points (x,t) . (This is not an unreasonable assumption, since for small e , a physical string should maintain 1 + ey^ near 1 .) Thus requiring (local) stability we get the conditions k rr , 2 2 1 — <J /3 for -r- > a 5- — h 3 3 k2 , 2 1. . , 2 2 (a - —) <i 1 for a > j h The second condition is more restrictive than the f i r s t , so i t is sufficient to choose k and h satisfying the second condition. For the nonlinear scheme (1), and for either g , we can summarize this condition by requiring ^2 ( m a x l-l + egj - y) < 1 h 70. For the string equation (II.IB) this can be written more simply as ^ (e + h 4 1 • h J This is the condition used to choose appropriate values of h and k for the computations performed in Chapter III. Recall that i t does not guarantee stability, first because i t was derived as an approximate condition by considering the constant coefficient problem, and secondly because i t follows from only a necessary condition for stability. 71. §2 Numerical Stability of Gaussian Elimination for the Evolution Equations In the remainder of this chapter we consider the difference schemes used to compute approximate solutions to the GKdV equation u + g(u)u + 62u = 0 . (13)-t x xxx As explained in Chapter III, §1, we assume periodic i n i t i a l data on [0,2] given by u(x,0) = <; (x) , and we require u(x,t) = u(x+2,t) for al l x and t . The approximate solutions are computed on the grid (III.3) We restrict our attention to implicit schemes of the form vn+l = vn - l _ n _ y2 n ( vn + 1 + v""1) . (14) J J 3 + o - j j The motivation for such schemes will be discussed in §3. As in §1, v = v(x,t) denotes the solution to the difference scheme (14), and by we mean v(x.,t ) . For the present we are not concerned with the 3 3 ri-te rm (j)1^ . It is assumed to be some approximation to the term g(u)u 3 x at (x^.t^) involving v only at times t < t^ , so that the system of equations obtained below will be linear. We look at this term more 2 closely in the following section. Now let Q = -6 D + D Q D and write (14) as (I-kQ)vJ+1 = -2k(j>*J + (I+kQ)Vj~ 1 Defining n , n n nN T V = ( VQ, V ]_ , . . . . V j ) 72. the scheme .can be written in matrix form as Qvn+1 = r (15) for some matrix Q and vector r . It is easily seen that Q has the form (here illustrated for J = 6, with dots representing nonzero elements) /•• • ' A .... •» * Q = V " * * * / The nonzero entries in the corners result from the periodic boundary conditions. In order to minimize the number of arithmetic operations required to solve this system numerically, i t would be advantageous to perform Gaussian elimination without pivotivig (GEWP) , for this preserves the form of Q . However i t is important to ensure that this process is well defined (that i s , a l l the natural pivot elements are nonzero), and that the lack of pivoting does not result in numerical instability of the algorithm. To make this more precise, recall (see for example Isaacson and Keller [13], Chapter II, §1) that solving the system (15) by GEWP is equivalent to decomposing Q into a product Q = LU , where L is lower triangular with diagonal elements 1 and U is upper triangular, and then solving the systems Lr = r and Uvn+"'' = r .We then say o o GEWP is (numerically) stable i f the elements of L and U are bounded 73. and the diagonal elements u_„ of U (which are the natural pivot elements) are a l l nonzero. We now give a restricted result. We show that the pivot elements u.. are a l l nonzero and bounded. Notice this proves that the vectors v +^ are a l l well defined, for the result implies that Q ^ exists. (A different proof of this fact has been given by Kreiss [25]; see Theorem 4.4 in the following section.) To study the system (15) we introduce the more general differential equation ufc + G(u) + P(x,t,3)-u = 0 , (16) where G(u) is some (possibly) nonlinear function of u and its x-ueriva-tives, and P is a linear operator of the form P(x,t,9)u = I ' {3j+1<eyu) + 8 j(B.9j+1u)} , j=0 J 3 3 with 3 E — and each 3. = 3.(x,t) real valued and 2-periodic in x 3X 3 3 We consider periodic solutions of (16), as before. Notice that the GKdV 62 equation occurs as the special case m = 1, 3Q = 0 > = 2 G(u) = g(u)u^ . Define the difference operator S = S(x,t) by m S = T {D-?D g .DJ + D?3.DJD } >0 + o j - + J - o S is clearly a discretization of P . Choose 0 < c x < l (a = y in (14)) and consider the difference scheme 74. n+1 n-1 „n . • ,' n+1 . . n-1, v. = v. - 2kG. - 2kS{av + (l-a)v. } , 3 3 3 3 3 where Gn. is a discretization of G at (x.,t ) which involves no J 3 n values of v at times t > t . A s before this ensures a linear system n of equations. The scheme can be written (I+2kaS)v^+1 = -2kG*J + (l-2k(l^)S)v"_ 1 , or in matrix form Qvn+1 = r (17) for some matrix Q and vector r . Since (15) is a special case of (17) we use the same notation. Clearly Q is "almost" a banded matrix, having a band of 2m + 3 nonzero diagonals along with elements in each corner corresponding to the periodic boundary conditions. Exactly the same comments regarding the solution of (15) by GEWP apply to the system (17).' We begin our stability proof with Lemma 4.1: Let the real matrix H of order q be skew symmetric. Then the eigenvalues of H are given by ± iX . j = 1,2,..., y i f q is even, or 0,± iX. j = 1,2 Q-^- i f q is odd, 3 *• 75. where each A is a nonnegative eigenvalue of the hermetian matrix iH Proof: Since H* = -H, (iH)* = (-i)(-H) = iH , or, iH is hermetian. Hence iH has real eigenvalues; let A be any such.. Then -iA is an eigenvalue of H . Since H is real its eigenvalues must occur in complex conjugate pairs. q.e.d. We now give the stability theorem. Stability for the system (17) will follow as a corollary. Theorem 4,2: Let the q x q matrix A = I + H where H is real and skew symmetric. Then A is non singular and the decomposition A = LU by GEWP is stable in the sense that the pivot elements u.. of U are bounded and nonzero: 2 1 4 u.. < 1 + A for a l l i , i i where A is the largest nonnegative eigenvalue of the hermetian matrix iH Proof: Let A^ denote the k x k matrix obtained from A as the upper left hand corner of k rows and columns. Do similarly for H, L and U . Recalling that L is lower triangular with diagonal elements 1; and that U is upper triangular with diagonal elements u.. , det A^ = det • det U, = 1 • II u. i=l i i Hence d £ t \ \k det Ak_] 76. Consider the case where k is even and write p = — . (The case where k is odd is similar). The determinant of A^ is the product of its eigenvalues. Let 0 < X, < . . . < X denote the nonnegative eigenvalues of 1 p 111^ and let 0 = U q < < . . . < y^ ^ denote those of IH^.i • Then by Lemma 4.1, P 2 det A, = n (1+X .) , * j=l J P"1 2 det A, = n (1+y.) . By a well known theorem (see e.g. Noble [36], p. 415) concerning eigen-values of hermetian matrices, i t follows that ° <. * i < ^ < \ < < V i v< y P - i * X P H e n c e n (l+x2) . 1 + 12 i=l 3 2 P _ 1 \k = ^ - ^ * ~ i p-l " 3 = 1 1+y n (i+y:) 3 j=i J 2 < 1 + xp , 2 P-l 1+AL and (1+X.) n J — kk 1 . .. ., , 2 3=1 l+y. J > 1 + X 2 > 1 The assertions of the theorem now follow. Corollary 4.3: The conclusions of Theorem 4.2 apply to the system (17) ~ n+1 Qv = r Proof: It need only be shown that the difference operator S defined on v n+1 defines a skew symmetric matrix operator S acting on n+1 First note that D,, D and D define matrix operators + o -M + , and M_ , respectively. Illustrating for J = 4 and , M 1_ 2h ( 0 1 -1 0 -1 -1 1 0 1 -1 0 -1 The following facts are easily verified. (i) M, , M and M commute (since they are circulant; see + o -Lancaster [29], p. 267); (ii) M£ = -M_ ; T ( i i i ) M = -M (M is skew symmetric). o o o Now define the diagonal matrices = ^..(t) with (Bj).^ = 3j(x^,t) . Note a l l B are of course real and symmetric. Then the operator S is j given by m . . . S = J (M^M B.M J + M^B.M^M ) . >0 + o j - + j - o To see that S is skew symmetric, consider i t termwise. 78. (M?M B.M^ + MIB.M^M )T ^ + 0 3 - + J - o = (MT)JB.MT(Mf)d + MT(M T )V ( M V - 3 0 + o - 3 + = (-l)2j+1{M?B .M + M M?B.M^ } + 3 0 - o + 3 -= -(M^ B.M^ M + M?M B.MJ) + J - o + 0 3 - ' = -(M?M B.M^ + M?B.MJM ) + 0 3 - + 3 - 0 Skew symmetry now follows and the proof is complete. Sjoberg [39] has considered a general implicit difference equation similar to the system (17). For such difference equations he showed that Gaussian elimination without pivoting is well defined in that no zeroes appear as pivot elements, and that the matrix operator (corresponding to n+1 Q ) acting on v was nonsingular (so that v was well defined for a l l times t ). In the case in which his a.'s are a l l zero, Corollary 4.3 n 3 gives an alternate proof of his results, and shows as well that the pivot elements are bounded. Care should be taken when reading his paper, for he does not clarify the role of the nonlinearity in ensuring that a linear system of equations is obtained. 79. §3 Kreiss' Stability Theorem: The Energy Method We now examine the stability of the difference schemes used to compute solutions to the GKdV equation (13). We begin with some preli-minaries . We assume approximate solutions to (13) are computed on the grid (III.3). As in §1, the approximate solution v to (13) obtained at times t (by solving some appropriate difference scheme) will be taken to be defined as a function vn(x) on the whole interval [0,2] . Then we use the = L^tO^] inner product 2 (v,w) = J vwdx 0 2 andnorm |[v|j = (v,v) . We say a difference scheme which gives an approximate solution to the GKdV.equation (13) is stable (or strongly stable in the sense of §1) i f there exists a constant K (independent of h and k) such that | | v nI < K « v°| | for a l l n . We will insist on this strong form of stability since we have shown (see equation 11.28) that, for a solution u of (13), [|u(x,tn)| = ||u(x,0)|| for a l l n , and so we would not expect an approximate solution which was allowed to grow (as the stability definition Jv^'fl < Ke"^nk||v* I would allow) to simulate the true solution adequately. We also say that 80, n+1 n-1 01 , / nN , .2. n . , n+1 n-1. v. = v. - 2kA(v.) - kS D,D D (v. +v. ) J J J + o - j j is a consistent approximation to the GKdV equation provided <}> (u(x^ . ,t^)) -g(u(x.,t))u(x.,t)-»-0 as h 0 for a l l points (x.,t ) . We use &j n x j n 1 j n the following theorem (see Kreiss and Widlund[25], p. 42) (which we state for the case of a scalar real-valued function v of a scalar variable x ) Theorem 4.4: Suppose n is a constant with 0 < n < 1 Consider a difference scheme of the form (I-kQ^) v(x,t+k) = (I+kQ^ v(x,t-k) + 2kQQv(x,t) , (18) where QQ and are operators acting on L such that Q^ v e for a l l v e and i = 0,1 . Suppose that for a l l v e (i) (v,Qov) = 0 , (ii) ( v . Q j V ) < 0 , ( i i i ) k|Q v|| « (l-n)|v|| . Then, given any v^ (x) , v"*"(x) e there is a unique solution v of the difference scheme (18) with n{|[vn|2 + | v n + 1 S 2 } < (2 -n ) { | |v° l 2 + fv 1 ! 2 } for a l l n . Notice that the conclusion of the theorem ensures stability according to our definition. For 81. Bvn||2<llv"|2+I|v"tl|2 We remark that Theorem 4.4 may only be applied to linear operators Q Q , for, although no assumption of linearity is made in the theorem (as stated by Kreiss),Kreiss' proof requires a lemma which is only valid (in general) for linear operators Qq . (See Kreiss and Widlund [25], Lemma 9.1, p. 43). We now examine the stability of the difference schemes. In Chapter III, §1, we indicated that explicit schemes are not suitable, and so we- restrict our attention to implicit schemes, and in particular, to schemes of the form (18). Notice that the schemes (III.6) - (III.8) and (III.13) are a l l of this form. It is of course precisely this fact which motivated the choice of these schemes. We first show stability of the frozen coefficient (or linearized) version of these schemes. It is given by the difference equation (I+y 2D,D D )vn + 1= ( I - M V D D ) vn _ 1 - 2kaD vn , (19) + o j + o - j o j ' where a is a constant. (Note that (19) may also be taken as a difference scheme giving approximate solutions to a particular example of the GKdV equation, namely, 2 u + au + 6 u = 0 , t x xxx 82, i.e., g(u) = a .) We require Lemma 4 . 5 : Let v,w be any 2-periodic (real valued) functions in L^. Then (i) (v,Dow) = -(DQV,W) (in particular (y,.D v) =0 ) (ii) (v,D+w) = -(D_v,w) Proof: The results follow easily from the definition of the inner product and the periodicity condition, q.e.d. We now prove stability of (19). Theorem 4 . 6 : The difference scheme (19) is stable provided kL£L< ! h < 2 Proof: (19) is of the form (18) with Q =r-<5 D,D D and 1 + o -QQ = ~3l>Q • By Lemma 4 . 5 , for any v e L^ (v,Q v) = (v,-aD v) = -a(v,D v) = 0 o o o and (v.Q-jV) =-62(V,D+DQD_V) =62(D_V,DO(D v)) = 0 . Finally, using the periodicity and Schwarz' inequality, 2 2 2 li r, t|2 r / ™ N2 , I t rv(x+h)-v(x-h), , |QQvfl = / (aDQv) dx = a / H> 3 d X o o . 2 2+h 2-h 2 = ^-jij v (x)dx + / v (x)dx - 2 / v(x+h)v(x-h)dx} 4 h h -h o 2 2 2 « ^2 + l v l 2 + 2t / v2(x+h)dx]%[J v2(x-h)dx]%} 4 h o o =4M2 > h 83. and the result follows by Theorem 4.4. q.e.d. Theorem 4.6 may also be proven using the von Neumann technique introduced in §1, although the proof is considerably longer. Now consider the nonlinear schemes themselves. For each of these schemes (i.e. III.6, 7,8 and 13) Qq and may be determined appropriately. In each case is given as in Theorem 4.6 (and so condi-tion (ii) of Theorem 4.4 is satisfied). But in each case Qq is non-linear. Therefore Theorem 4.4 may not be applied to these schemes and so we cannot present a proof of stability. It is possible however to state Theorem 4.4 in a form which may be applied to certain nonlinear Qq . We now present this form, even though we were unable to apply i t to prove stability of any of the schemes at hand. Theorem 4.7: Theorem 4.4 remains true i f condition (i) is replaced by the requirement that ( i ) ' (V,QQW) + (QQV,W) < 0 for a l l v,w e ; Qq may be a nonlinear operator. Proof: As given by Kreiss and Widlund [25]. Note however that Lemma 9.1 ([25], p. 43) is not required. (It should be noted that, in the original statement of Kreiss' theorem which allows v to be a complex vector valued function of a vector variable x , the proper form of ( i ) ' i s Re(v,Q w) + Re(Q v,w) ^ 0 '^ o o 84. for a l l v,w e .) Theorem 4.7 is stronger than Theorem 4.4 since for linear Q o (v,Q v) = 0 for a l l v e L„ implies o Z (v,Q w) + (Q v,w) = 0 for a l l v,w e L0 o o / To conclude we wish to pursue a suggestion made by Kreiss regarding Sjoberg's scheme (III.6) for solving the KdV equation (III.l). For this 1 2 scheme, Q v = - —(D (v ) + vD v) . Notice that o 3 o o (v,Q v) = --|{(v,Dov2) + (v,vD v)} = - | { . ( v , D o V 2 ) + ( v2, D Q v ) } = - "|{(v,Dov2) - (v,D o V 2 ) } 0 . That i s , Qq satisfies the condition (i) of Theorem 4.4. For this reason, Kreiss[24] stated that Sjoberg's scheme (III.6) should be used to compute approximate solutions to the KdV equation, in preference, say, to (III.7) or (III.8). However the reason for this suggestion is not clear since, as we have already mentioned, Theorem 4.4 does not apply to nonlinear QQ On the other hand, let us assume that Kreiss had another reason (which we have been unable to ascertain) for choosing Qq to satisfy ( i ) . In the event that such a reason exists, we present a scheme for the GKdV 85. equation which satisfies ( i ) . In fact we consider the more general equa-tion u + g(u)u + ( £ 6 3 2 + 1) u = 0 , (20) j=0 J 8x J where each 6^ is a constant (not a l l of which are zero), and we compute approximate solutions to (20) using a difference scheme (I-kQ ) v(x,t+k) = 2kQQv(x,t) + (I+kQ ) v(x,t-k) ,., (21) where Qq is some operator (to be specified below) and Q^ is given by m Q = - y 6 .D-JD Dj . (22) 1 , L n i + o -j=0 J (It follows from §2 that this scheme is well defined; that i s , i f (21) is written as a system of linear equations, the matrix corresponding to I-kQ^ is non-singular. That the scheme is well defined also follows from Theorem 4.4, since the existence part of the proof does not involve the nonlinear operator QQ .) Theorem '4.8: Consider the equation (20). Let g be continuous, Define d(u) = j ug(u)du . u Then with Q1 defined by (22) and Qq by -QQV = DQ(d(v)v) + dDQV , (23) 86, the difference scheme (21) is a consistent approximation to (20) and the conditions (i) - ( i i i ) of Theorem 4.4 are satisfied provided 2 - max | d| < 1 . (24) h x,t Proof: We consider the derivation of d(u) ; i t will then be obvious that (21) is a consistent approximation. We wish to split g(u)u and write g(u)uv = (d.(u)u) + d„(u)u X X X £~ X and then define -Q v H D (dn(v)v) + d_(v)D v . o o 1 Z o To satisfy condition (i) of Theorem 4.4, we must have 0 = (v,QQv) = - {(v,Do(d1(v)v)) + (v,d2(v)DQv)} = - {(v,DQ(d1(v)v)) + (vd2(v),DQv)} = - {(v,Do(d1(v)v)) - (v,Do(vd2(v)))} . This will be true provided d^ = d2 = d . Then g(u)u = (ud(u)) + d(u)u X X X = (ud'(u) + 2d(u))ux and d must satisfy the ordinary differential equation 87. ud'(u) + 2d(u) = g(u) A solution of this equation is given by d(u) = -y / ug(u)du . It is u clear then that Qq defined in (23) gives a consistent approximation to (20) satisfying condition ( i ) . Consider ( i i ) . - ( v . Q ^ v ) = (v,(£ 6 D J D D ^ ) V ) = Y 6 . ( V , D ! D Djv) u j + o -= £ ( - I ) V ( A , D O ( A ) ) = 0 by Lemma 4.5, and (ii) holds. Finally assume max |d(v)| exists and xeI0,2] t>0 denote i t by M . Write w(x,t) = d(v(x,t)) . Then, suppressing t in our notation, and using the periodicity and Schwarz' inequality }ov|2 -||Do(wv) + wDov|2 1 2 = — j I (w(x+h)v(x+h) - w(x-h)v(x-h) 4h o + w(x)[v(x+h) - v(x-h)]}2dx « J' { [| w(x+h^ | + | w(x) | ] [ v(x+h) 1 4h o + [|w(x)[ + [w(x-h)| ][ v(x-h)[ }2dx iM2 2 ^2 / {|v(x+h)[ + |v(x-h)| }Zdx 4h o 88. Thus ( i i i ) holds i f < 1 . q.e.d. For the linear GKdV equation (i.e. g(u) E a ) , Theorem 4.8 gives the scheme (19), and hence Theorem 4.6 follows as a Corollary of Theorem 4.8. For the KdV equation, Theorem 4.8 gives Sjoberg's scheme (III.6), as we would expect. - But consider the BKdV equation (II.8). Then d(u) = Now d is unbounded at zero, the maximum M does not exist and the stability condition (24) cannot be satisfied for any choice of k except 0 89. V Conclusions We wish to briefly summarize our results. These can be grouped into those concerned with the solitons, and those concerned with stability of the difference equations. Consider the former. We have demonstrated that, from smooth in i t i a l data, solitons, that is solitary wave pulses, do appear in solutions to the BKdV equation. We have also shown that the interaction of these wave pulses is responsible for the energy recurrence observed in solutions to the nonlinear wave equation (1.3). We have established (numerically) that these solitary wave pulses are in fact particular solutions to the BKdV equation. These results are similar to those obtained by ZK for the KdV equation. We have determined certain properties of the solitons. We showed that they are uniquely determined (when they exist) by their speed and background and the sign of their amplitude, and that they are symmetric unimodal travelling waves. We note that we did not analytically establish the role of solitons in general solutions to the BKdV equation. Such results have been given for the KdV equation by Lax [32] and GGKM [9]j see the discussion at the end of §3 in Chapter I. These results could not be extended since they rely on the association of the KdV equation with the Schrddinger operator (see [32]). In the general case, no tech-nique is known for associating a suitable operator with the GKdV equation, and so we are unable to predict the size, speed or number of solitons which 90. will arise from arbitrary i n i t i a l conditions. For the difference schemes our results concern two different equations. First we have considered the generalization 2 = a y + By , (1) tt xx xxxx a and B constants, of the well known wave equation yt t= a yx x ' ( 2 ) The difference scheme y(x,t+k) - 2y(x,t) + y(x,t-k) = n + D y ( X j t ) . ( 3 ) is well known to give stable approximations to (2) provided Ictlk « 1 . (4) h We have taken the natural generalization of (3) to the equation (1) obtained 2 2 by adding the term BD+ D y(x,t) to (3) and have given a necessary condition for stability for this scheme. We note that for small B it is just a per-turbation of the well known condition (4). Secondly we examined the BKdV equation (IV.13). Following the work of Kreiss we computed solutions to (IV.13) using schemes of the form (I + k62D D D )v.n + 1 = 2kg(v.n)D v.n + (I - k62n D D )v.n _ 1 (5) + o - J J ° J + o - j We showed that the constant coefficient problem corresponding to (5) is < 1. stable provided 1. We also noted that, although there seem 91. to be some theoretical reasons for preferring simulations for the nonlinear term of a more complicated nature than the natural g (v_.n)Do (v^.n) , there is no practical evidence to support the theory. Finally we note that schemes of the form (5) lead to a system of linear equations (I + H)x = b where H is skew-symmetric and I is the identity matrix. We demonstrated that such systems are nonsingular and that Gaussian elimination without pivoting is well defined and leads to an LU decomposition in which the pivot elements are bounded both above and away from zero. 92. • TIME Figure 1 MODAL ENERGY DECOMPOSITION OF AN APPROXIMATE SOLUTION TO WAVE EQUATION 1.1. THE ORIGINAL FPU RECURRENCE. Parameters: e = 1/64, h = .03125(J =32) , k = .01104854 (k/h = 1//8) , Initial data sin irx TIME Figure 2 MODAL ENERGY DECOMPOSITION OF AN APPROXIMATE SOLUTION TO JOHNSON'S WAVE EQUATION 1.3. RECURRENCE IS OBSERVED FOR JOHNSON'S EQUATION. Parameters: e = 1/8, h = .03125(J = 32), k = .01104854, Initial data sin T T X 93. Figure 3 MODAL ENERGY DECOMPOSITION OF AN APPROXIMATE SOLUTION TO WAVE EQUATION 1.1 COMPUTED WITH. NON-UNIFORM SPATIAL DISCRETIZATION. Parameters: e = 1/64, J = 33, k = .0013681, Initial data sin TTX . Figure 4a MODAL ENERGY DECOMPOSITION OF AN APPROXIMATE SOLUTION TO WAVE EQUATION 1.1 COMPUTED WITH J NOT A POWER#OF 2. Parameters: e = 1/64, h = 0.3(J = 30), k = .01104854, Initial data sin T T X . Note: J even. cf. Figure 4b. o r~. ~ Figure 4b MODAL ENERGY DECOMPOSITION OF AN APPROXIMATE SOLUTION TO WAVE EQUATION 1.1 COMPUTED WITH J NOT A POWER OF 2. Parameters: e = 1/64, h = .0285714(J = 35), k = .01104854, Initial data sin T T X . Note: J odd. cf. Figure 4a. 9 5 , Figure 5 MODAL ENERGY DECOMPOSITION OF AN APPROXIMATE SOLUTION TO WAVE EQUATION 1.7 WITH g-TERM, B = h2/12. Parameters: e = 1/64, h = .03125 (J = 32), k = .025, Initial data sin TTX 0.40 J\ .1.00 A .1.60 J-L2.20 L12.80 L J.3.40 J_ 4.00 0.20 h 0.80 1.40 .2.00 L2.60 _ 3.20 JL 3.80 11 0.00 0.60 1.20 1.80 2.40 3.00 ± 3.60 Figure 6 APPROXIMATE SOLUTION OF KDV EQUATION III.l SHOWING SINGLE TRAVELLING WAVE PULSE. Initial data sech2((x-1)/A) , A 2 = 1262 . Parameters <52 = .000484, h = .01 (J = 200), k .005 Figure 7 APPROXIMATE SOLUTION OF KDV EQUATION III.l COMPUTED USING SPLIT SCHEME III.6. •Initial data cos T T X .. Parameters: 62 .0.01302083, (= 1/768), h = .02 (J = 100), k = .01 . Figure 8 APPROXIMATE SOLUTION OF KDV EQUATION III.l COMPUTED USING NATURAL SCHEME III.7. Initial data cos T T X . Parameters: 62 = .001302083 (= 1/768), h = .02 (J = 100), k = .01 . 98. Figure 9 APPROXIMATE SOLUTION OF KDV EQUATION III.l COMPUTED USING CONSERVATION SCHEME i n . 8. Initial data cos T T X . Parameters: <52 = .001302083 (= 1/768), h = .02 ( J = 100), k = .01 . 99. Figure 10a MODAL ENERGY DECOMPOSITION OF AN APPROXIMATE SOLUTION TO WAVE EQUATION II.1 FPU, MODES 1-5. (See Figure 10b for Modes 6-8) Parameters: e = 1/16, h = .03125 (J = 32) , k = .025, B = h2/12 . Initial data L sin T T X . TIME Figure 10b MODAL ENERGY DECOMPOSITION OF AN APPROXIMATE SOLUTION TO WAVE EQUATION II.1 FPU, MODES 6-8. (See Figure 10a for Modes 1-5). Parameters: e = 1/16, h = .03125 (J = 32), k = .025, B = h2/12 . Initial data £ sin T T X . TT 100. Figure 11a MODAL ENERGY DECOMPOSITION OF AN APPROXIMATE SOLUTION TO WAVE EQUATION II.1 FPU, MODES 1-5. (See Figure lib for Modes 6-8). Parameters: e = .0688705 (L/(3 x 4.84)) , h = .02 (J = 50), k = .016, 3 = h2/12, Initial data — sin T T X IT Figure lib MODAL ENERGY DECOMPOSITION OF AN APPROXIMATE SOLUTION TO WAVE EQUATION II.1 FPU. MODES 6-8. (See Figure 11a for Modes 1-5). Parameters: e = .0688705, h = .03125 (J = 32), k = .025, 3 = h2/12 . Initial data 1 s i n ^ TT Figure 12 APPROXIMATE SOLUTION OF KDV EQUATION II I . l . Initial data j cos Parameters: S2 = .0052085 (= 1/192), h = .02 (J = 100), k 1.30 M L , » l 2 . 7 0 102. 4.10 1101*15.50 H IbdJL6 .90 IUU/U,8.30 lUMlU 9• 70 1.20 MJLjdl2.60 4.00 5.40 6.80 8.20 U I A J U U 9 . 6 0 1 0 U 2 - 5 0 /U43-90 \kM5-3° AA6'10 8.10 U U l i i / 9 . 5 0 A ; •°° E42-40 ft/iA/t3-80 \/IJW5-20 Jut6-60 ww8-00 9.40 °-9° I)U42:30 IvM^ 3-70 & M 5 A 0 AA6-™ WAA90 XW9-30 ^ .00 / ^ y ^ 6 . 4 0 y ^ | ^ 7 . 8 0 0.80 0.70 '•6° fkh lL^ 2-10 i lk 3 - 5 0 A l i i 4 - 9 0 AA 6- 3 0 7.70 9.20 9.10 0.60 1)^2.00 ^ i^3-40 ^l44-80 T^A 6.20 0.50 A , X 1 . 9 0 0.40 21__/1.80 1JUUJL3.20 A A fl 4.60 3-3° AAA*-10 h h 6 A 0 MA 7.50 JLUflWB.90 6.00 U.vyJ, 7.40 8.80 i i i k . 0.30 \ X\ .70 HLMi3 .10 A M K 4.50 7.30 AJULH 8.70 0.20 IX X- 1 -60 5.90 3.00 JJU/t4.40 JUU_iliW_ 5.60 JUliLi7 . 2 0 OJLUJU8.60 o. io _ \ - ^ i . 5 o f ^ ^ 2 . m j y ^ 4 . 3 o 5 . 7 0 y y ^ 7 . i o l y y y ^ a . s o ^rA 2 - 6 5 i|/44-20 lU 5^-60 Afl/A7-00 ^ iv6-40 0.00 A /_ 1 .40 Figure 13 APPROXIMATE SOLUTION OF KDV EQUATION II I . l . Initial data cos T T X Parameters: Sz = .000484, h = .008 (J = 250), k = .0025 . Figure 14 APPROXIMATE SOLUTION OF THE MODIFIED KDV EQUATION 11.10. Initial data cos T T X . Parameters: 62 = .001302083 (= 1/768), h = .02 CJ = 100), k = .01 . cf. Figure 8. Figure 15 APPROXIMATE SOLUTION OF THE SECOND MODIFICATION OF THE KDV EQUATION, 11.11. Initial data cos T T X . Parameters: <S2 = .001302083 (= 1/768), h = .02 ' (J = 100), k = .01. cf. Figure 8. 105. Figure 16 MODAL.ENERGY DECOMPOSITION OF AN APPROXIMATE SOLUTION TO WAVE EQUATION II.1 B. Parameters: re = .1, h = .02 (J = 50), k = .016, 6 = h2/12, Initial data — sin T T X Figure 17 MODAL ENERGY DECOMPOSITION OF AN APPROXIMATE SOLUTION TO WAVE EQUATION II 1 B Parameters: e = 1/4, h = .03125 (J = 32), k = .025, 3 = h2/12, Initial data s in TT x TT 106. 5.00 n nru. / \ . _/-s n nn / , ya.nn 29.00 A x r 7 V A 4.00 A. 10.00. i R nn / y ^ _ 22.00 28.00 3 ^ 3.00 \ q nn ~ lA i5.oo v / \ ^ r 2 1 . 0 0 27.00 2.00 A^^X^ R nn ~/ y /-^ 1 4 nn. / l ^ ^ ?n nn / ^ -"\ 26.00 T l .00 \ S 1 nn ^ \ /•>• i3.00^LX-^w= 19.00 ,/ 25.00 3.00 A ' 6.00 J V 1 ? nnv_/ y ~ r I R nn / V^ -N ? 4 nn \ Figure 18 1 APPROXIMATE SOLUTION TO THE BKDV EQUATION III. 12. Initial data j cos TTX Parameters: S2 = .0003, h = .04 (J = 50), k = .01 . 107, Figure 19 -APPROXIMATE SOLUTION TO THE BKDV EQUATION III. 12. Initial data y cos TTX . Parameters: S2 = .00032552083, h = .04 (J = 50), k = .01 . 108. TIME Figure 20 MODAL ENERGY DECOMPOSITION OF AN APPROXIMATE SOLUTION TO WAVE EQUATION II. I B . Parameters: e = .1, h = .015625 (J = 64), k = .014, B = h2/12, Initial data J _ sin T T X 2TT o TIME Figure 21 MODAL ENERGY DECOMPOSITION OF AN APPROXIMATE SOLUTION TO WAVE EQUATION II.IB. Parameters: e = .1, h = .02 (J = 50), k = .016, B = h2/12, Initial data sin TTX 109. . 4.00 ^^J^ 8"5° -AA^y 13.00,1^ ^.SO.ApA 22.00A^ __7 ZB.SO^ .^. 31. 00 3.50 M ^8.00 1.00 4 ^ 7 . 5 0 1 2 . 0 0 , ^ 1 6 . 5 0 / ^ 2 1 . 0 0 1 ^ 2.50 AJl^ll.OO^U^lS.SOT^^ 2.00 I I ^6.50 1.50 TV ^6.00 / V L , , A 10-50 A . / \ P . 15.00A^K 19.50\Arv J '24.00,/^^ 28.50 T 25.00,/RV\ 29.50 24.50, A \ 29.00 1.00 r\_^5.50 10.0014.50/1,.A/U 19.0Q.AXA /I 23.50 0.00 N^/* 4-50 0.50 "\ X 5 . 0 0 A / l . ,,^9.50 J1A-A^ 1 4 . 0 0 ^ 1 U 18.50.AA 7123.00 1 3 . 5 0 J i A A J L . 18 .00 Figure 22 3 APPROXIMATE SOLUTION TO THE BKDV EQUATION III. 12. Initial data . -j cos TTX . Parameters: 62 = .00020345052083, h = .02 (J = 100), k = .005 . 110, Figure 23 APPROXIMATE SOLUTION TO THE BKDV EQUATION III. 12. Initial data j cos TTX . Parameters: <52 = .0003, h = .013 (J = 150), k = .004 . 111. Figure 25 FOURIER SERIES DECOMPOSITION OF THE APPROXIMATE SOLUTION TO THE BKDV EQUATION APPEARING IN FIGURE 23 . 0.60 0.40 0.20 0.00 Figure 26 APPROXIMATE SOLUTION OF BKDV EQUATION III..12. A TRAVELLING WAVE. Initial data: an approximate soliton. Parameters: 62 = .0003, h=.02 (J=100), k=.001. Figure 27 COMPARISON OF SOLITARY WAVE PULSE (i) AND APPROXIMATE SOLITON. Parameters: 62 = .0003, c=.22, s = 3.03. Dotted line is soliton. o to CT ^ i O — i 1 Figure 28 COMPARISON OF SOLITARY WAVE PULSE (ii) AND APPROXIMATE SOLITON. Parameters: 62 = .0003, c = .14, sQ = 1.92. Dotted line is soliton. Figure 29 COMPARISON OF SOLITARY WAVE PULSE ( i i i ) AND APPROXIMATE SOLITON. Parameters: 62 = .0003, c = .032, sn = .78. Dotted line is soliton. 114. BIBLIOGRAPHY 1. BENJAMIN, T.B., 1972, The Stability of Solitary Waves, Proc. Roy. Soc, V.A328, pp.153-183. 2. BENJAMIN, R.B., BONA, J.L., and MAHONEY, J.J., 1972, Model Equations for long waves in nonlinear dispersive systems, Phil. Trans. Royal Society of London, v. 272, #1220, pp.47-78 3. BIRKHOFF, G., and ROTA, G.C., 1962, Ordinary Differential Equations, 1st ed. (Ginn and Company, Waltham, Mass., 1962). 4. CANOSA, J., 1971, On a Nonlinear Equation of Evolution, IBM Scientific Center Report, Palo Alto, Calif. 5. COLE, J.D., 1951, On a Quasi-linear Parabolic Equation Occurring •in Aerodynamics, Q. Appl. Math., v.9, p.225. 6. COOLEY, J.W., and TUKEY, J.W., 1965, An Algorithm for the Machine Calculation of Complex Fourier Series, Math, of Computation, v.19, #90, pp. 297-301. 7. FERMI, E., PASTA, J.R., and ULAM, S.M., 1955, Studies of Nonlinear Problems, I, Los Alamos Report #La-1940 (May, 1955) (unpublished). 8. FORD, J., 1961, Equipartition of energy for Nonlinear Systems, J. Math. Phys., v.2, p.787. 9. GARDNER, C.S., GREENE, J.M., KRUSKAL, M.D. and MIURA, R.M., 1967, Method for solving the Korteweg-de Vries Equation, Phys. Rev. Letters, v. 19, #19, p.1095. 10. GARDNER, C.S., 1971, Korteweg-de Vries Equation and Generalizations,IV. The Korteweg-de Vries Equation as a Hamiltonian System, J. Math. Phys., v.12, p.1548. 11. GREENBERG, J.M., MACCAMY, R.C., and MIZEL, V.J., 1968, On the Existence, Uniqueness and Stability of Solutions of the Equation o'(u )u + Au _ = p u , J. Math and Mech., v.17, #7, 1968. X X X xtx o tt 12. HIROTA, R., 1971, Exact Solution of the Korteweg-de Vries Equation for Multiple Collisions of Solitons, Phys, Rev. Letters, v.27, pp.1192-1194. 13. ISAACSON, E. and KELLER, H.B., 1966, Analysis of Numerical Methods, (John Wiley Sons, New York, 1966). 115. 14. JOHNSON, G.D., 1967, On a Nonlinear Vibrating String, Ph.D. Dissertation, U.C.L.A., 1967. 15. KAMETAKE, Y., 1969, Korteweg-de Vries Equation I. Global Existence of Smooth Solutions, P. Jap. Acad., v.45, p.552. 16. KAMETAKE, Y., 1969, Korteweg-de Vries Equation II. Finite Difference Approximation, P. Jap. Acad., v.45, p.556. 17. KAMETAKE, Y., 1969, Korteweg-de Vries Equation III. Global Existence of Asymptotically Periodic Solutions, P. Jap. Acad., v.45, p.656. 18. KAMETAKE, Y., 1969, Korteweg-de Vries Equation IV. Simplest Generalization, P. Jap. Acad., v.45, p.661. 19. KANTOROVICH, L.V. and KRYLOV, V.I., 1964, Approximate Methods of Higher Analysis (trans, by Curtis D. Benster, Interscience Publishers, New York, cl964). 20. KARPMAN, V.I., 1968, Some Asymptotic Relations for Solutions of the Korteweg-de Vries Equation, Phys. Letters, V.A26, p.619. 21. KARPMAN, V.I. and SOKOLEV, V.P., 1968, On Solitons and the Eigenvalues of the Schrodinger Equation, Sov. Ph. JER, v.27, p.839. 22. KELLER, H.B., 1968, Numerical Methods for Two-Point Boundary Value Problems, (Blaisdell, Waltham, Mass.', 1968). 23. KREISS, H.O., 1960, Uber die Losung von Angangsrandwertaufgaben fur Partielle Differentialgleichungen mit Hilfe von Differenzengleichungen, Trans. Royal Inst, of Technology, Stockholm, #166. 24. KREISS, H.O., 1969, Course Notes, Spring 1970, Cal. Inst, of Tech. 25. KREISS, H.O. and WIDLUND, 0., 1967, Difference Approximations for Initial Value Problems for Partial Differential Equations, Report NR 7, (Sept. 1967) Department of Computer Science, Uppsala University, Sweden. 26. KRUSKAL, M.D., 1965, Asymptotology in Numerical Computation: Progress and Plans on the Fermi-Pasta-Ulam Problem, in Proceedings of the IBM Scientific Computing Symposium - Large Scale Physical Problems in Physics, Dec. 9-11, 1963 at IBM Thomas J. Watson Research Lab. 27. KRUSKAL, M.D. and ZABUSKY, N.J., 1964, Stroboscopic Perturbation Procedure for Treating a Class of Nonlinear Wave Equations, J. Math. Phys. v.5, pp. 231-244. 116. 28. KRUSKAL, M.D., MIURA, R.M. , GARDNER, CS. and ZABUSKY, N.J. , 1970, Korteweg-de Vries Equation and Generalizations, V. Uniqueness and Nonexistence of Polynomial Conservation Laws, J. Math. Phys., v.11, #3, p.952. 29. LANCASTER, P., 1969, Theory of Matrices (Academic Press, New York). 30. LAX, P.D., 1964, Development of Singularities of Solutions of Nonlinear Hyperbolic Partial Differential Equations, J. Math. Phys., v.5, p.611. 31. LAX, P.D. and NIRENBERG, L., 1966, On Stability for Difference Schemes; a Sharp Form of Garding's Inequality, CP.A.M., v. 19, #4, pp.473-492. 32. LAX, P.D., Integrals of Nonlinear Equations of Evolution and Solitary Waves, CP.A.M., v.21, pp.467-490. 33. LAX, P.D., 1969, Nonlinear Partial Differential Equations and Computing, SIAM Review, v. 11, #1, p.7. 34. MIURA, R.M., 1968, Korteweg-de Vries Equation and Generalizations I. A Remarkable Nonlinear Transformation, J. Math. Phys., v.9, #8, p.1202. 35. MIURA, R.M., GARDNER, CS. and KRUSKAL, M.D., 1968, Korteweg-de Vries Equation and Generalizations II. Existence of Conservation Laws and Constants of Motion, J. Math. Phys., v.9, #8, p.1204. 36. NOBLE, B., 1969, Applied Linear Algebra, (Prentice-Hall, Inc., Englewood Cl i f f s , New Jersey, 1969). 37. RITT, J., 1948, Integration in Finite Terms; Liouville's Theory of Elementary Methods (Columbia University Press, New York, 1948). 38. SHIH, L.Y., 1971, The Korteweg-de Vries Equation: Transition Period, J. Math. Phys., v.12, #10, pp.2124-2126. 39. SJOBERG, A., 1969, Numerical Solution of the Korteweg-de Vries Equation, Report NR 25, (Oct. 1969), Department of Computing Science, Uppsala University, Sweden. 40. SJOBERG, A., 1969, On the Korteweg-de Vries Equation; Existence and Uniqueness, Department of Computing Science Report, Uppsala University, Sweden. 41. SU, CH. and GARDNER, C.S., 1969, Korteweg-de Vries Equation and Generalizatijns, III. Derivation of the Korteweg-de Vries and Burger's Equation, J. Math. Phys., v.10, #3, p.536. 117. 42. TUKEY, J.W., 1958, The Estimation of (Power) Spectra and Related Quantities, in On Numerical Approximation. (Proceedings of a Symposium Conducted by the Mathematics Research Center, United States Army, at the University of Wisconsin, Madison, April 21-23, 1958), ed. R. Langer. 43. TEMAM, R., 1969, Sur un Probleme Non Lineaire, Jour, de Math. P.A., v.48, p.159. 44. ZABUSKY, N.J., 1963, Phenomena Associated with the Oscillations of a Nonlinear Model String (The Problem of Fermi, Pasta and Ulam), in Proceedings of the Conference on Mathematical Models in the Physical Sciences, Stefan Drobot ed., pp. 99-133 (Prentice-Hall, Englewood Cl i f f s , New Jersey, 1963). 45. ZABUSKY , N.J., 1962, Exact Solution for the Vibrations of a Nonlinear Continuous Model String, J. Math. Phys., v.3, pp. 1028-1039. 46. ZABUSKY, N.J. and KRUSKAL, M.D., 1965, Interaction of Solitons in a Collisionless Plasma and the Recurrence of Initial States, Phys. Rev. Letters, v.15, pp. 240-243. 47. ZABUSKY, N.J., 1967, A Synergetic Approach to Problems of Nonlinear Dispersive Wave Propogation and Interaction, in Nonlinear Partial Differential Equations: A Symposium on Methods of Solution. (Academic Press, 1967, ed. W. Ames). 48. ZABUSKY, N.J., 1968, Solitons and Bound States of the Time-Independent Schrodinger Equation, Phys. Rev., v.168, #1, p.124. 49. ZABUSKY, N.J., 1969, Nonlinear Lattice Dynamics and Energy Sharing, J. Phys. Jap., v.26 (5), p.196. 50. ZABUSKY, N.J. and GALVIN, C.J., 1971, Shallow-water Waves, the Korteweg-de Vries Equation and Solitons, J. F l . Mech., v.47, pp.811-824.
- Library Home /
- Search Collections /
- Open Collections /
- Browse Collections /
- UBC Theses and Dissertations /
- Numerical simulation of a nonlinear wave equation and...
Open Collections
UBC Theses and Dissertations
Featured Collection
UBC Theses and Dissertations
Numerical simulation of a nonlinear wave equation and recurrence of initial states Buckley, Albert Grant 1972
pdf
Notice for Google Chrome users:
If you are having trouble viewing or searching the PDF with Google Chrome, please download it here instead.
If you are having trouble viewing or searching the PDF with Google Chrome, please download it here instead.
Page Metadata
Item Metadata
Title | Numerical simulation of a nonlinear wave equation and recurrence of initial states |
Creator |
Buckley, Albert Grant |
Publisher | University of British Columbia |
Date Issued | 1972 |
Description | In 1955 Fermi, Pasta and Ulam (FPU) [7] observed an unusual recurrence to initial state in numerical solutions of a nonlinear wave equation. Zabusky and Kruskal (ZK) [47] have subsequently found an explanation for this phenomenon based on special travelling wave solutions ("solitons") of the (nonlinear) Korteweg de Vries (KdV) equation. In this thesis we extend ZK's explanation to a similar nonlinear wave equation given by Johnson[14]. We investigate existence and uniqueness of solitons for a (nonlinear) generalization of the KdV equation. (Chapter II) and present computational results to illustrate ZK's soliton explanation of the recurrence, both for FPU's equation and Johnson's equation (Chapter III). In Chapter IV we give some results concerning the stability of the difference schemes used to obtain solutions to the nonlinear partial differential equations. |
Genre |
Thesis/Dissertation |
Type |
Text |
Language | eng |
Date Available | 2011-02-25 |
Provider | Vancouver : University of British Columbia Library |
Rights | For non-commercial purposes only, such as research, private study and education. Additional conditions apply, see Terms of Use https://open.library.ubc.ca/terms_of_use. |
DOI | 10.14288/1.0080359 |
URI | http://hdl.handle.net/2429/31823 |
Degree |
Doctor of Philosophy - PhD |
Program |
Mathematics |
Affiliation |
Science, Faculty of Mathematics, Department of |
Degree Grantor | University of British Columbia |
Campus |
UBCV |
Scholarly Level | Graduate |
AggregatedSourceRepository | DSpace |
Download
- Media
- 831-UBC_1973_A1 B82.pdf [ 4.89MB ]
- Metadata
- JSON: 831-1.0080359.json
- JSON-LD: 831-1.0080359-ld.json
- RDF/XML (Pretty): 831-1.0080359-rdf.xml
- RDF/JSON: 831-1.0080359-rdf.json
- Turtle: 831-1.0080359-turtle.txt
- N-Triples: 831-1.0080359-rdf-ntriples.txt
- Original Record: 831-1.0080359-source.json
- Full Text
- 831-1.0080359-fulltext.txt
- Citation
- 831-1.0080359.ris
Full Text
Cite
Citation Scheme:
Usage Statistics
Share
Embed
Customize your widget with the following options, then copy and paste the code below into the HTML
of your page to embed this item in your website.
<div id="ubcOpenCollectionsWidgetDisplay">
<script id="ubcOpenCollectionsWidget"
src="{[{embed.src}]}"
data-item="{[{embed.item}]}"
data-collection="{[{embed.collection}]}"
data-metadata="{[{embed.showMetadata}]}"
data-width="{[{embed.width}]}"
data-media="{[{embed.selectedMedia}]}"
async >
</script>
</div>
Our image viewer uses the IIIF 2.0 standard.
To load this item in other compatible viewers, use this url:
https://iiif.library.ubc.ca/presentation/dsp.831.1-0080359/manifest