This de¬nition is an algorithm for computing the ¬rst order Taylor term. We merely

multiply (x) by a dummy variable and then take the partial derivative of N (u, x; )

with respect to in the usual way, that is, with u and x kept ¬xed. The dummy variable

takes us away from the complicated world of functional analysis into the far simpler realm

of ordinary calculus.

The need for the “direction (x)” in the de¬nition of the Frechet differential can be

seen by expanding u(x) in a spectral series. The operator equation N (u) = 0 becomes a

system of nonlinear algebraic equations. When a function depends upon more than one

unknown, we cannot specify “the” derivative: we have to specify the direction. The Frechet

differential is the “directional derivative” of ordinary multi-variable calculus.

What is important about the Frechet derivative is that it gives the linear term in our

generalized Taylor expansion: for any operator N (u),

N (u + ) = N (u) + Nu (u) 2

) [Generalized Taylor Series] (C.13)

+ O(

For the special case in which N (u) is a system of algebraic equations,

2

(C.14)

F a+ = F (a) + J +O

the Frechet derivative is the Jacobian matrix J of the system of functions F (a).

With this de¬nition in hand, we can easily create a Newton-Kantorovich iteration for

almost any equation, and we shall give many examples later. An obvious question is: Why

C.2. EXAMPLES 529

bother? After all, when we apply the pseudospectral method (or ¬nite difference or ¬nite

element algorithms) to the linear differential equation that we must solve at each step of

the Newton-Kantorovich iteration, we are implicitly computing the Jacobian of the system

of equations that we would have obtained by applying the same numerical method ¬rst,

and then Newton™s method afterwards.

There are two answers. First, it is often easier to write down the matrix elements in

the Newton-Kantorovich formalism. If we already have a subroutine for solving linear

differential equations, for example, we merely embed this in an iteration loop and write

subroutines de¬ning F (x, u) and Fu (x, u). It is easy to change parameters or adapt the

code to quite different equations merely by altering the subroutines for F (x, u), etc.

The second answer is that posing the problem as a sequence of linear differential equa-

tions implies that we also have the theoretical machinery for such equations at our disposal.

For example, we can apply the iterative methods of Chapter 15, which solve the matrix

problems generated by differential equations much more cheaply than Gaussian elimina-

tion.

In principle, we could use all these same tricks even if we applied Newton™s method

second instead of ¬rst ” after all, we are always in the business of solving differential

equations, and reversing the order cannot change that. Conceptually, however, it is much

simpler to have the differential equation clearly in view.

C.2 Examples

The second order ordinary differential equation

(C.15)

uxx = F (x, u, ux )

generates the iteration

’ Fu x ’ Fu = F ’ uxx

(i)

(C.16)

xx x

(i)

where F and all its derivatives are evaluated with u = u(i) and ux = ux , and where (x),

as in the previous section, is the correction to the i-th iterate:

u(i+1) (x) ≡ u(i) (x) + (C.17)

(x)

The Frechet derivative is more complicated here than for a ¬rst order equation because u(x)

appears as both the second and third arguments of the three-variable function F (x, u, ux ),

but the principle is unchanged.

For instance, let

F (x, u, ux ) ≡ ±(ux )2 + exp(x u[x]) (C.18)

The needed derivatives are

(C.19)

Fu (x, u, ux ) = x exp(x u)

(C.20)

Fux (x, u, ux ) = 2 ± ux

We can also extend the idea to equations in which the highest derivative appears non-

linearly. Norton (1964) gives the example

1

u uxx + A(ux )2 + B u = B 20 + (C.21)

sin(πx)

12

APPENDIX C. NEWTON ITERATION

530

which, despite its rather bizarre appearance, arose in analyzing ocean waves. We ¬rst

de¬ne N (u) to equal all the terms in (C.21). Next, make the substitution u ’ u + and

collect powers of . We ¬nd

N (u + ) = u uxx + A(ux )2 + B u ’ B 20 + 1

(C.22)

sin(πx)

12

2

+ {u }+ 2

+ uxx + 2 A ux +B + A( x)

xx x xx

To evaluate the Frechet derivative, we subtract the ¬rst line of the R. H. S. of (C.22) ” the

term independent of ” and then divide by . When we take the limit ’ 0, the term in

2

disappears as we could have anticipated directly from its form: it is nonlinear in , and

the whole point of the iteration is to create a linear differential equation for . Thus, the

Frechet derivative is equal to the linear term in the expansion of N (u + ) in powers of ,

and this is a general theorem. Our generalized Taylor series is then

N (u + ) ≈ N (u) + Nu (u) (C.23)

and equating this to zero gives

Nu (u) = ’N (u) (C.24)

which is

u(i) (i)

+ (u(i) + B) (C.25)

+ 2 A ux

xx x xx

1

= ’ u(i) u(i) + A(u(i) )2 + B u(i) ’ B 20 + sin(πx)

xx x

12

Inspecting (C.25), we see that the R. H. S. is simply the -independent term in the power

series expansion of N (u + ). The L. H. S. is simply the term in the same series which

is linear in . Once we have recognized this, we can bypass some of the intermediate steps

and proceed immediately from (C.22) to (C.25). The intermediate equations, however, are

still conceptually important because they remind us that the iteration is simply a generalized

form of Newton™s method.

We have not discussed boundary conditions, but these are usually obvious. If the ¬rst

guess satis¬es the boundary condition, then we can safely impose homogeneous boundary

conditions on (x) at all iterations. However, we must be careful to impose conditions that

are consistent with those on u(x). If we have Dirichlet conditions, u(’1) = ±, u(1) = β,

then we would set (±1) = 0. With Neuman conditions on u(x), du/dx(’1) = ± and

du/dx(1) = β, then we would impose d /dx(±1) = 0. Finally, periodic boundary condi-

tions, i. e. using a Fourier series, should be imposed on (x) if the problem speci¬es such

conditions for u(x).

The Newton-Kantorovich method applies equally well to partial differential equations

and to integral equations. For example,

u2 + u2 = 0 (C.26)

uxx + uyy + x y

generates the iteration

ux + uy

x y

(C.27)

+ +

xx yy

2 2

(i) (i)

ux + uy

2 2

(i) (i)

’ uxx + uyy +

(i) (i)

ux + uy

C.3. EIGENVALUE PROBLEMS 531

We can create an iteration even for the most bizarre equations, but the issues of whether

solutions exist for nonlinear equations, or whether they are smooth even if they do exist,

are very dif¬cult ” as dif¬cult as the Newton-Kantorovich method is elementary.

C.3 Eigenvalue Problems

When the equation is nonlinear, the distinction between eigenvalue and boundary value

problems may be negligible. A linear boundary value problem has a solution for arbitrary

forcing, but a linear eigenvalue problem has a solution only for particular, discrete values

of the eigenparameter. It is rather different when the eigenvalue problem is nonlinear.

For example, the periodic solutions of the Korteweg-deVries equation (“cnoidal waves”)

and the solitary waves (also called “solitons”) of the same equation satisfy the differential

equation

uxxx + (u ’ c) ux = 0 [Korteweg-deVries Eq.] (C.28)

where c is the eigenvalue. The KDV equation can be explicitly solved in terms of the so-

called elliptic cosine function, and this makes it possible to prove that solutions exist for a

continuous range of c:

Solitary waves exist for any c > 0 (C.29)

u(’∞) = u(∞) = 0 :

Cnoidal waves exist for any c > ’1. (C.30)

u(x) = u(x + 2π) :

What makes it possible for the eigenparameter to have a continuous range is that for a

nonlinear differential equation, the relative magnitude of different terms changes with the

amplitude of the solution. For the Korteweg-deVries problem, the wave amplitude for either

of the two cases is a unique, monotonically increasing function of c.

For other problems in this class, the wave amplitude need not be a monotonic function

of a nor is the solution necessarily unique. For example, if we replace the third degree term

in (C.28) by a ¬fth derivative, one can show (Boyd, 1986b) that in addition to the usual

solitons with a single peak, there are also solitons with two peaks, three peaks, etc. (The

multi-peaked solitons are actually bound states of the single-peaked solitons.)

The existence of solutions for a continuous range of the eigenparameter, instead of merely

a set of discrete values of c, is a common property of nonlinear equations. It follows that

to compute a solution to (C.28), for example, one can simply choose an almost-arbitrary

value of c ” “almost-arbitrary” means it must be larger than the cutoffs shown in (C.29)

and (C.30) ” and then solve it by the Newton-Kantorovich method as one would any other

boundary value problem. The equation for the correction (x) to u(i) (x) is

+ (u(i) ’ c) = ’u(i) + (u(i) ’ c) u(i)

(i)

(C.31)

+ ux

xxx x xxx x

The only possible complication in solving (C.31) is that the operator on the L. H. S.

has two eigensolutions with zero eigenvalue. In a word, (C.31) does not have a unique

solution. The physical reason for this is that the solitons of the Korteweg-deVries equation

form a two-parameter family.

The ¬rst degree of freedom is translational invariance: if u(x) is a solution, then so is

u(x + φ) for arbitrary constant φ. If we Taylor expand u(x + φ) and keep only the ¬rst term

so as to be consistent with the linearization inherent in (C.31), then

du

+ O(φ2 ) (C.32)

u(x + φ) = u(x) + φ

dx

APPENDIX C. NEWTON ITERATION

532

must be a solution of the linearized Korteweg-deVries equation for arbitrary φ. It follows

that the solution of (C.31) is not unique, but is determined only to within addition of an

arbitrary amount of

du

≡ (C.33)

1 (x)

dx

which therefore is an eigenfunction (with zero eigenvalue) of the linear operator of (C.31).

Fortunately, the solitons are symmetric about their peak. The derivative of a symmetric

function is antisymmetric. It follows that the column vector containing the spectral coef-

¬cients of du/dx is not a linear eigenfunction of the discretized form of (C.31) when the

basis set is restricted to only symmetric basis functions (either cosines, in the periodic case,

or even degree rational Chebyshev functions on the in¬nite interval).

The second degree of freedom is that if u(x) is a solution of (C.28), then so is

√

v(x) = ’ + (1 + ) u 1 + x (C.34)

for arbitrary . (In technical terms, the KdV equation has a continuous Lie group symmetry,

in this case, a “dilational” symmetry.) Again, we must Taylor expand v(x) so as to remain

consistent with the linearization inherent in (C.31). We ¬nd that the other eigenfunction

[the O( ) term in (C.34)] is

x du

≡ ’1 + u(x) + (C.35)

2 (x)

2 dx

We can exclude 2 (x) on the in¬nite interval by imposing the boundary condition

u(±∞) = 0. Adding a term proportional to 2 (x), although consistent with the differ-

ential equation (C.31), would violate the boundary condition because 2 (x) asymptotes to

a non-zero constant at in¬nity. In a similar way, we can remove the non-uniqueness for the

periodic case by imposing the condition that the average of u(x) over the interval is 0. The

simplest way to implement this condition is to omit the constant from the basis set, using

only the set {cos(nx)} with n > 0 to approximate the solution.

Once we have computed the soliton or cnoidal wave for a given c, we can then take that

single solution and generate a two-parameter family of solutions by exploiting these two

Lie group symmetries. Thus, the translational and dilational symmetries are both curse

and a blessing. The symmetries complicate the calculations by generating two eigenmodes

with zero eigenvalues. However, the symmetries are also a blessing because they reduce