Eight definitions of the slow manifold: seiches, pseudoseiches and exponential smallness

, .dynamics

or at ,mospneres arl£1 o c e a n s

ELSEVIER

Dynamics of Atmospheres and Oceans 22 (1995) 49-75

Eight definitions of the slow manifold: seiches, pseudoseiches and exponential smallness J o h n P. B o y d Department of Atmospheric, Oceanic and Space Science, University of Michigan, 2455 Hayward Avenue, Ann Arbor, MI 48109, USA Received 20 August 1993; revised 6 June 1994; accepted 21 July 1994

Abstract

The slow manifold is the Holy Grail of initialization schemes for weather forecasting, a hypothetical subspace of the model's phase space which is free of high-frequency (and forecast-wrecking) gravity waves. Using new machinery, we analyze seven previous definitions of the slow manifold and one definition original to this work. One new tool is the mathematics of nonlocal solitary wave theory. The traditional definition of a slow manifold is too restrictive, but a generalization in the spirit of nonlocal soliton theory is computable. Another novelty is the conceptual distinction between 'seiche' and 'pseudoseiche'. The former are free, unforced gravity waves, whereas pseudoseiches are forced gravity mode oscillations which mimic seiches over a finite time interval. We show through a linear model that it is possible to remove seiches by initialization, even though the usual algorithms are divergent. However, intermittent gravity waves, i.e. pseudoseiches, appears as inexorable as death and the tides.

1. Introduction

R i c h a r d s o n (1922) discovered the hard way that insertion of unaltered observations into a forecast model is disastrous. H e predicted a c h a n g e of 145 m b a r for a 6 h forecast for a day in 1910 when the actual pressure change was less than 1 m b a r (Platzman, 1969; Lynch, 1992). The p r o b l e m is that m e a s u r e m e n t errors, being random, are m u c h m o r e ageostrophic than the atmosphere. This generates large-amplitude gravity waves which wreck the forecast. T h e goal of initialization schemes is to adjust the observations to a nearly geostrophic state like the real atmosphere. This goal was formalized w h e n Leith (1980) introduced the c o n c e p t o f the 'slow 0377-0265/95/$09.50 © 1995 Elsevier Science B.V. All rights reserved

SSDI 0 3 7 7 - 0 2 6 5 ( 9 4 ) 0 0 3 9 6 - 3

50

J.P. Boyd / Dynamics of Atmospheres and Oceans 22 (1995) 49-75

Table 1 Eight flavors of the slow manifold Name

Definition

References

Bandlimited

Fourier transform F(w) of f(t) vanishes for I w I > W No seiches or pseudoseiches

Warn (unpublished data, ca. 1982), Lynch and Huang (1992) Lorenz and Krishnamurthy (1987), Lorenz (1992)

Fast modes are unique functions of Rossby amplitudes dku/dt k .~ O(1), O(E k)

Vautard and Legras (1986), Jacobs (1991), Lorenz (1994) Kreiss (1980), Browning et al. (1980)

No seiches

Baer (1977), Baer and Tribbia (1977), Tribbia (1984)

Gravity mode ordinary differential equations replaced by ruthorder Baer-Tribbia Minimizes norm of time derivatives of the gravitational modes over a finite interval Seiches and pseudoseiches are exponentially small

Daley (1980), Debussche and Temam (1994)

Pseudoseichefree (hysteresis free) Gravity mode servitude Bounded derivative Seiche-free multiple scales (Baer-Tribbia) Differentialalgebraic (DA(m)) Optimal control

Fuzzy

Courtier (1987), Lorenc (1988), Courtier and Talagrand (1990) T. Warn (unpublished data, ca. 1982)

manifold': a hypothetical N-dimensional manifold embedded in the 3N-dimensional phase space of a primitive equations model which is devoid of gravity wave activity. Long before his paper, however, the purpose of initialization schemes was to adjust the raw data onto a state of minimal gravity waves. Unfortunately, there is no straightforward, direct way to minimize the highfrequency ripples. Practical initialization schemes have therefore attacked a variety of mathematical proxies. Most of the literature refers to the slow manifold, as if the state of minimal gravity waves was uniquely defined. In reality, different algorithms for minimizing different proxies for the gravity waves generally compute distinct N-dimensional manifolds. Table 1 is a precis of eight definitions of the slow manifold. The 'pseudoseichefree' and 'fuzzy' manifolds are not computational procedures but rather attitudes towards how much 'slowness' can be achieved. The other six manifolds are defined by specific algorithms. The next section is a brief introduction to the normal mode (Hough basis) formalism, which is the most natural framework for discussing the slow manifold. Sections 3 - 7 analyze linear slow manifolds, that is, the application of standard initialization methods to a forced linear oscillator. As the hydrodynamic equations are nonlinear, this linear analysis might seem an irrelevant diversion. However, to a first approximation, the gravitational modes are passively forced by the interactions of the slowly evolving Rossby modes. The reward for neglecting the weak modification of the Rossby forcing by the gravity

J.P. Boyd / Dynamics of Atmospheres and Oceans 22 (1995) 49-75

51

waves and the equally weak self-interaction of the gravitational modes is that the linear analysis is straightforward. In particular, it is possible to make a clear distinction between free gravity waves or seiches, which can be suppressed by the proper choice of initial condition, and pseudoseiches, which are forced solutions with small, high-frequency oscillations superimposed on slower variations. For the forced, linear equation, we can show that (1) the oscillations in the pseudoseiche cannot be eliminated but only minimized and (2) the minimum amplitude of the high-frequency oscillations is exponentially small in l / e , where E is the Rossby number. In Section 8, we discuss how the physics of the slow manifold is modified from that of the forced, linear manifold by the addition of full gravity mode with gravity mode nonlinearity. Although rigorous proof is not possible, we can offer some good arguments that these modifications are small and quantitative rather than qualitative. Section 9 is a compendium of brief miscellaneous comments on the various definitions of the slow manifold. Section 10 is a catalogue of unresolved issues.

2. The Hough function representation The natural framework for slow manifold theory is to represent the spatial dependence of the flow as a series of Hough functions, that is, as a generalized Fourier series in which the basis functions are the normal modes of a resting atmosphere. The forecast model becomes a set of coupled ordinary differential equations for the time dependence of the modal amplitudes. Each of the modes is the product of a height-dependent factor, which satisfies an ordinary differential equation in height ('vertical structure equation') with a latitude-and-longitude factor which is a normal mode of the shallow water wave equations (Kasahara, 1977, 1978; Daley, 1991). We omit the gory details because the numerical coupling coefficients and the spatial structure factors are irrelevant to slow manifold theory. What matters is that the dynamics has been transformed into a system of ordinary differential equations in time, and that the unknowns may be partitioned into two sets: geostrophic mode amplitudes and gravity mode amplitudes. The simplest model which captures the mathematical essence of the slow manifold and also shows most clearly its connection with nonlocal solitons is that introduced by Lorenz (1986) and Lorenz and Krishnamurthy (1987). The model is obtained by truncating the Hough series to just five modes: three Rossby modes and two gravity waves. A similar philosophy of maximum truncation of a spectral series for fluid convection gave the three equations now known as the Lorenz system, which are a star attraction of any text on chaos and dynamical systems theory. The five-mode set, henceforth labeled the 'LK Quintet', plays an analogous role in discussions of the slow manifold. After nondimensionalization and rescaling of coordinates, which eliminate all

52

J.P. Boyd ~Dynamics o f Atmospheres and Oceans 22 (1995) 49-75

explicit parameters except the nondimensional parameters shown below, the model is the system of five ordinary differential equations Ut = - V W ; b ' V z - a U F } V t = U W - b Uz - a V + Rossby triad Wt = - UV-

aW

,LK Quintet

(1)

X t = - - Z -- a x

z t = b U V + x - az I gravity diad

where a is the damping coefficient, b and b' are nonlinear coupling coefficients, and F is the forcing, a constant independent of time. We adopt the convention of denoting the amplitudes of the three Rossby waves by upper-case letters (U, V, W ) and the gravity waves, which are much smaller in amplitude, by lower-case ( x , z). Although (1) is derived from fluid mechanics, the unknowns are the amplitudes of terms in a generalized Fourier series, not coordinates or velocities. The two real-valued first-order equations which form the gravity dyad may be written in two alternative forms. First, setting a = 0 for simplicity, the dyad is equivalent to the single, real-valued second-order equation Ztt -~ Z ~

-bl'Vtt

(2)

The other option is the single, first-order complex-valued equation u, + io~u = f ( t )

(3)

where u = z +/x, oJ = 1 + ia, and f ( t ) is the nonlinear forcing. The reason for introducing the generic symbols oJ and f into (3) is that these three options--a coupled pair of first-order equations, a single second-order equation, and a single first-order but complex-valued equation--are available for more realistic models, too. There is one equation of the form of (3) for each complex-valued gravity mode which is kept in the truncation.

3. Linear slow manifolds

One important characteristic of initialization algorithms is that they can be applied to linear systems, too. This observation would be insignificant were it not for the weakness of the nonlinear interactions of the gravity waves with each other. The forcing of the gravity waves, f ( t ) in (3), is dominated by the Rossby modes. Because of this, 'It appears to be valid for a wide variety of models where the fast variables may be regarded as constituting a forced damped linear oscillator, with the forcing...[f(t) in (3)] supplied entirely or mainly by the slow variables' (Lorenz and Krishnamurthy, 1987, p. 2949). It follows that it is a good first approximation to pretend that (3) is a forced, linear equation and ask: what is the slow manifold, as computed by various algorithms, for this linear model? Below, we shall address the complications when the model is fully nonlinear.

J.P. Boyd/ Dynamicsof Atmospheresand Oceans22 (1995)49-75

53

It should be noted that (3) actually retains the dominant nonlinear terms in the sense that Rossby-Rossby interactions are not directly approximated in any way. Only the gravity mode equations are linear, and linear only in the sense of being linear in the gravity mode amplitudes; these equations contain quadratic products of Rossby amplitudes in f(t).

4. Seiches and pseudoseiches The general solution to (3) is the sum of a homogeneous solution plus a particular solution p(t): u =Aexp(-ioot)

+p(t)

(4)

where A is an arbitrary constant determined by the initial condition, i.e. A = u(0) - p ( 0 ) . The homogeneous solution is a free, unforced gravity wave. For brevity, we shall borrow the oceanographic term for free oscillations in a lake and refer to the homogeneous solution as a 'seiche'. Constructing the linear slow manifold then seems easy: we choose a particular integral which is free of wiggles and then set u(0) = p(0) so that A = 0. The difficulty is that it is impossible, except in special cases, for the particular integral to be wiggle-free. In other words, a linear slow manifold does not exist in the strict, classical sense. However, a generalized slow manifold does exist in the sense that if we cannot eliminate the gravity waves completely, we can at least minimize them. For example, let us define

r( t ) = e x p ( - R o t )

exp(

-st 2)

(5)

Then the solution to

Pt + itop = r( t )

(6)

1 ~1/2 . p(t) = -~( ~ ) erfl(s)l/2t] exp(-iwt)

(7)

is

as shown in Fig. 1 for a typical choice of the width parameter s. Asymptotically,

l(rrll/Zlexp(-itot), p(t)~-~ s] ~-exp(-iwt),

t~ t~-~

} (8)

Although the forcing is strictly localized in the sense that the amplitude of r(t) decays exponentially with t, the solution is nonlocal: the flow oscillates forever at frequency to. If we add a seiche to (8), that is, a multiple of exp(-ioJt), we can reduce the amplitude of the oscillations for either large negative or large positive t, but not both. The particular integral always has seiche-like oscillations. This motivates the following.

54

J.P. Boyd / Dynamics of Atmospheres and Oceans 22 (1995) 49-75 Solution: rbD(t)=G.5 sqr(pi/2~) er[(t/5)

~(

i t)

r h c )

10

1,~ t ,,p:

F©rc~:

a

0

~

I-

1 ( r )

r(t)=oxD(t't/25)

oN)(-t

t]

VIV ]r]

/~k

1

-1 i

i

-I0

3

i0


Fig. 1. Upper panel: real part of the solution p(t) = (1/2X~-/s) 1/2 erf(sl/2t) e x p ( - i~ot) for s = 1/25 and o)= 1 to the equation Pt + i~op = r(t), (The graph was scaled by dividing by the constant (1/2Xrr/s)l/2.) Lower panel: the real part of the forcing r(t) -=--e x p ( - st 2) e x p ( - ioot).

Definition 1. Pseudoseiche A pseudoseiche is a solution which includes oscillations which intermittently mimic a sine wave, i.e. a free gravity wave or seiche. This definition is useful in making a distinction between what is and is not possible for a linear slow manifold. The exemplary solution p(t) shows that a pseudoseiche-free slow manifold, that is, a manifold which is completely free of oscillations with apparent frequency w, is not always possible--at least for the linear model (3). Fig. 2 illustrates a second solution family that is seiche-free but not pseudoseiche-free. The forcing m(t) and solution ~ ( t ) are

m(t;s,P) =s e x p ( - i w t )

•

(-1) nsech2[s(t-nP)]

"~= -~ i~(t;s,P) = e x p ( - i o J t ) ~ ( - 1 ) n t a n h [ s ( t - n P ) ] n=

(9)

--0o

This family of solutions is temporally periodic with period P and therefore can be expanded as a Fourier series. If P is not an integral multiple of 2~-/oa, none of the frequencies in the Fourier series is equal to ~o. The function ~(t,s P) is obviously the slow manifold in the sense that it lacks a component equal to w. Nevertheless, the response is oscillating with roughly unit frequency for all t including the gaps between the peaks of the forcing. Does full nonlinearity help? A third example suggests that it does not. Boyd (1994) explicitly computed temporally periodic nonlinear solutions to the LK Quintet (1) with a = F = 0, i.e. no damping or forcing. The periodic solutions for b' = 0, which are strictly linear in the gravity waves and constitute a linear slow

J.P. Boyd/ Dynamics of Atmospheres and Oceans 22 (1995) 49-75 !~zlutiTn: R ~ I

i~a:~ of :a;(t)

(

15

F~,rcL~{4:9~1

Z~st of m(t;s

-1 O

55

0

i'

[,~ { Z :

i'3

Fig. 2. Upper panel: real part of the periodic solution ~z(t, s = 1, P = 6~'). Lower panel: the forcing re(t; 1, 6~'). manifold, differ only in the fourth decimal place from the fully nonlinear manifold even for fairly strong nonlinearity. This third example suggests that allowing the Rossby-driven gravity waves to self-interact (and also in turn to alter the Rossby dynamics) does not eliminate the high-frequency oscillations of the slow manifold. In the next section, we shall describe a general method for solving (3) and demonstrate through a 'regularization' procedure that the generic linear slow manifold is a pseudoseiche, like these three examples.

5. Regularization and the genericity of pseudoseiches in the linear slow manifold To understand the regularization procedure about to be described, it is helpful to rewrite forced linear oscillator Eq. (3) as

u t + i~ou = E f ( e t )

(10)

so as to display explicitly its dependence on the Rossby number e, which both specifies the amplitude of the forcing of the gravity waves by the Rossby waves and also measures the slow time scale on which the Rossby modes evolve. Then, whenever f ( t ) has a well-behaved Fourier transform where the transform is defined by 1 f ( t ) = (277.)1/2 f_ e - i w t f ( w ) d w (11) (10) has the formal solution u(t)=

1

f

~

(2.n-) '/2 - =

F(w/e) e -iwt

i(w - w )

dw

(12)

J.P. Boyd/Dynamics of Atmospheres and Oceans22 (1995) 49-75

56

Solid: abs(U(w)) Dashed:abs(F(w))

10~ 100 10.~ 10-5 10-3 10-, 10-~ 10-¢ 10-~ 10.,

012

014

016

018

112

114

116

1'.8

w

Fig. 3. Schematic graph of the absolute value of the Fourier transform of the solution to (10), IU(w)L (continuous curve) and of the corresponding forcing function, IF(w) l (dashed curve), to = 1. Different choices of integration path around the pole at w = to are equivalent to adding multiples of e x p ( - i t o t ) to the solution and thus implicitly represent the usual constant of integration multiplying the homogeneous solution. The pole in the integrand also amplifies components of the forcing for which w = to as illustrated in Fig. 3. The spike in the Fourier transform of u(t) is not the spectrum of a true seiche, e x p ( - i t o t ) , which would be a Dirac t%function centered on w = to. However, the pole-induced p e a k is sufficient to produce oscillations with apparent frequency to. These oscillations make (12) a pseudoseiche. To demonstrate this assertion, we can eliminate the pole in the integrand through a procedure we shall call 'regularization' (of the Fourier integral). The key idea is to add and subtract a multiple of the explicit, analytic pseudoseiche p(et). We choose the multiplier q so that the subtraction cancels F ( w / e ) at w = to. As the Fourier transform of r(t) ( - e x p ( - s t 2) e x p ( - i t o t ) ) is

R(w)

(2s)1/2 exp

(13)

this multiplier is

q = (2s)l/2F(to//e)

(14)

We can write, without approximation,

u ( t ) = o'(t) + q p ( t )

(15)

J.P. Boyd / Dynamics of Atmospheres and Oceans 22 (1995) 49-75

57

.4

o

~,,

~,

A~,,~

~,

~,v~'~v,"~ ~ a

0

~3

0

~

: 0

4O

Z~

n

g~

i ,

i ,

40

.. i

6O

60

i L

i i

o ,

60

Fig. 4. The solution to u t + i u = - • sech(et), its 'nonlocal part', q p ( t ) , and their difference or(t)=- u The imaginary parts o f u, q p, and o- are shown as the top, middle, and bottom panels, respectively, e = 1 / 3 , s = 1/100. q p ( t ) is also shown as the dotted curve in the top panel. For this forcing a n d with to = 1, qp(t) = ( ~ ' / 2 ) s e c h [ r r / ( 2 e ) ] s i n ( t ) erf(tsl/2).

qp(t).

where

2 )1%/o, e It s l' l and ~r(t) is the Fourier transform of a function which does not have a pole at W=to:

1 ~(t)

o~

(2rr) 1/2 f _ - e -iw'

F(w/e) -F(to/e) exp[-(w-to)2/(4s)] dw i(to-w) (17)

Standard Fourier transform theory (Boyd, 1991a) implies that tr(t) decays exponentially fast with It]. Thus, the gravity mode oscillations for large t are contained entirely within the p(t) term and the amplitude of these oscillations is proportional to F(to/e). tr(t) is the pseudoseiche-free function which is the classical ideal of what the slow manifold should be. The additive pseudoseiche, qp(t), is the messy reality: the smooth, slowly varying function has oscillations superimposed. Although these oscillations are the superposition of many frequencies via the Fourier transform, collectively, the frequency band centered on the pole at w = to mimics a sine wave of the single frequency to. Fig. 4 illustrates the decomposition of u(t) into its 'local' and 'nonlocal' parts, or(t) and qp(t). The shapes vary with the width parameter s, which determines how rapidly the error function in p(t) rises away from the origin, but the qualitative nature of the decomposition is independent of s.

J.P. Boyd / Dynamics of Atmospheres and Oceans 22 (1995) 49-75

58

6. Exponential smallness

Theorem 1 If F ( w ) is the Fourier transform of f ( t ) , then f(t)

analytic for all [ I m ( t ) l < A = l F ( w ) l < B e x p [ - A R e ( w ) ]

(18)

In other words, F(to/e), which is the amplitude of the high-frequency oscillations, is O ( e x p ( - A t o / e ) ) , where A is the width of the strip of analyticity for f ( t ) in the complex t-plane. Under the mild restrictions that (i) f ( t ) is analytic for real t and (ii) f ( t ) has a well-behaved Fourier transform, it is a generic property of the solutions to the forced, linear oscillator that the amplitude of the pseudoseiche is exponentially small in l / e , where e is the Rossby number. As an illustration, to = 1 for the LK Quintet (1). For the special case a = F = b' = 0, Lorenz and Krishnamurthy (1987) have shown that the Rossby modes have the family of solutions U= • sech(et),

V=• tanh(et),

W= -• sech(et)

(19)

These hyperbolic functions of time have simple poles at +_irr/2, so A = 7r/2 and the amplitude of the nonlocal, high-frequency oscillations should be proportional to e x p [ - T r / ( 2 •)]. Lorenz and Krishnamurthy (1987) have shown that it is:

z(t)=-b~rexp-~e

sin(Itl),

Itl--'~

(20)

just as predicted by Theorem 1. This exponential smallness is a great blessing. When • = O(1/10), say, the amplitude of the nonlocal oscillations is only O(10-6). It is for this reason that standard initialization schemes work as well as they do. This exponential smallness is also a great complication because e x p ( - A / • ) is not analytic as • ~ 0. This implies that most initialization schemes are asymptotic but divergent, as explained in the next section.

7. The canonical initialization and divergence For (3), that is, u t + itou = f ( t ) , at least four different algorithms, including multiple scales and the bounded derivative principle, all define the same manifold via the series

io~ n=0k to } which we shall refer to as the 'canonical' initialization. Table 2 lists four different derivations of (21) with references. The series (21) is really only the special case t = 0 of a multiple scales perturbation theory which approximates u(t) on the slow manifold for all t.

J.P. Boyd / Dynamics of Atmospheres and Oceans 22 (1995) 49-75

59

Table 2 Four derivations of the canonical initialization Name

Method or condition

References

Multiple scales

u(t) is a function of the 'slow' time

Baer (1977), Baer and Tribbia (1977), Tribbia (1984) Machenhauer (1977), Kreiss (1980), Lorenz (1980), Browning et al. (1980) Tribbia (1984), Daley (1991, p. 287) This work

T = ~t Bounded derivative

ui(O)= 0 where u i = jth time derivative

Power series

Expand f(t), u(t) in power series and re-sum the homogeneous solution Expand 1/(~o - w) in integrand of Fourier transform as geometric series in w

Fourier transform

The Fourier transform derivation, which is given here for the first time, substitutes the geometric series 1 ¢O--W

O) n = 0 \ O )

into the Fourier solution

u(t)

1 (2rr)

1/z f

~ -~

e-iwt

r(w)

i(w - w )

dw

(23)

The infinite series of Fourier integrals is transformed by the identity dnf _

dt"

1

(2~r)T]2£ J -iwl"e-iwtF(w)dw

(24)

into a series of derivatives which, at t = 0, is simply (21). This derivation also suggests that the initialization series is asymptotic rather than convergent. The range of the integration variable in the Fourier transform is infinite, but the geometric series converges only for Iwl < ~o. The divergence of the geometric series on most of the integration interval forces the initialization series (21) to diverge, too. Three other arguments for the divergence of (21) are summarized in Table 3. The empirical argument is illustrated by Fig. 5, which shows the terms of even degree in the initialization series (21) for the unforced, inviscid LK Quintet when the Rossby modes are given by the hyperbolic solutions (19). (As noted above, the complex first-order Eq. (3) can be written alternatively as a second-order equation for z = Re(u); the even n terms in (21) are equivalent to perturbatively solving zt, + z = - b e 3 [ 2 sech3(et) - sech(et)l.) The terms show the typical behavior of asymptotic but divergent series: for fixed e, the terms decrease, reach a minimum, and then grow without bound with increasing degree. The minimum defines the optimal truncation order nopt for a given Rossby number e.

60

J.P. Boyd / Dynamics of Atmospheres and Oceans 22 (1995) 49-75

Table 3 Four arguments for divergence Name

Argument

References

Empirical

Numerical observation

Exponential smallness

e x p ( - ,~/~) is not analytic at e = 0, so the radius of an e-power series for this function is zero If the power series of f ( t ) has a finite radius of convergence d, then terms of the initialization series must be O(n! d - " ) Geometric series for 1~(to - w) is applied beyond its radius of convergence

T. Warn (unpublished data, 1983), Lorenz (1980), Lorenz (1986), Warn and Menard (1986) Warn and Menard (1986), Lorenz and Krishnamurthy (1987), Boyd (1991b) Vautard and Legras (1986)

Forcing power series

Fourier transform

Another exp(-A/e) implies that scales series

This work

a r g u m e n t f o r d i v e r g e n c e is e x p o n e n t i a l s m a l l n e s s . T h e f u n c t i o n is n o t a n a l y t i c a t • = 0, b u t r a t h e r h a s a n e s s e n t i a l s i n g u l a r i t y . T h i s t h e r a d i u s o f c o n v e r g e n c e o f a p o w e r s e r i e s i n e, s u c h as t h e m u l t i p l e w h i c h g i v e s (21), m u s t h a v e a z e r o r a d i u s o f c o n v e r g e n c e . L o r e n z a n d

le_3

+

e~s=i/3

E

e~s:I/6

o

6~s:i/9 , x

t:~:~:i/12

I00

i0

1

1

•O l

•001

le-4 "~(. le-5 le-6

2¢" "'w --*('"

i

i i0

i

i 20

order n

Fig. 5. The absolute error, scaled by dividing by e 3, in the multiple scales approximation to z(0) for the simplified inviscid LK Quintet where ztt + z = - - [ - - 6 sech(et)]tt, for four different values of e. The optimal truncation order nop t is that perturbation order which gives the minimum error for given e, i.e. the trough in each of the curves shown• The first ten terms for u(0) (and by conjecture, all) are the absolute values of the (2n + 2 ) - d Euler numbers, which asymptote to (2/~') 2n +2(2n + 2)!

J.P. Boyd / Dynamics of Atmospheres and Oceans 22 (1995) 49-75

61

z(0)

11

+~/n error

.i j

Omin term

J

.01 u n

I,~litude c~_non]ocal

%

.001 le-4

~

"/'t.

s

c le-5 a le-61 e1 le-7

~ •

de le-91e-8

.... '~ '-~

r le-10

O~'

~.'.?,,

o le-ll r le-12

•

le-13,

~

%

le-14 le-!5 I

__

i

,

0

i0 i/e~ilon

i

L

20

LnfiF~ty hbrm I F .i ! .01

u

[

+min error

-

.001 ]

-

e

-

'a

_~

d le-7 |

Omin term

-

n s le-4 ! c le-5 ~ a . 1 le-6

~amplitude of ncxdocal

Q

-

le-8 r le-9 : r =. o le-10 I

-

r le-ll !

-

le-12

!

le~13 ~

-

-

le-14 l_

~

, 0

--I~ i/epsilon

20

Fig. 6. Dashed curve with asterisks: a = ~ - e x p [ - ~ ' / ( 2 e ) ] , which is the amplitude of the nonlocal oscillations in time for the hyperbolic solution to the LK Quintet. Continuous line with plus signs: the error in that term in the perturbation theory which has the smallest error for the given e,n = hopt, Dotted line with circles: that term in the perturbation series, of whatever order n, which has the smallest magnitude. Panels differ only in the definition of 'smallest error'. Top panel: 'Smallest' means values at t = 0 only. Bottom panel: 'Smallest' is the L~ norm, that is, the error and the magnitude of the perturbation term are evaluated for t e [ - ~ , ~ ] and the largest value is taken before comparing different orders n. The error norm of the optimally truncated series is equal to a, the amplitude of the nonlocal oscillations, so the dashed and continuous curves are indistinguishable.

62

J.P. Boyd / Dynamics of Atmospheres and Oceans 22 (1995) 49-75

Krishnamurthy (1987), the numerical observations of Warn and Menard (1986), and the more general Fourier transform argument (Sections 5 and 6) all suggest that an exponential dependence on 1 / c is a generic feature of the slow manifold. The E-power series and other equivalent initialization schemes must then diverge. The exponentially small part of the slow manifold not only creates the divergence of the multiple scales series, but also controls the error in the optimally truncated series. Let a denote the amplitude of the gravitational oscillations of the pseudoseiche. For the special solution of Lorenz and Krishnamurthy (1987) to (1), a is the amplitude of the oscillations of unit frequency which persist long after the Rossby modes have returned to a non-interacting steady state. Eq. (20) shows = 7rb e x p [ - T r / ( 2 E ) ]

(25)

Fig. 6(a) compares a(E) with (i) the error in the optimally truncated initialization series (21) and (ii) the magnitude of the smallest term in that series. For all E, both the minimum error and the minimum term are O(a). Indeed, if we repeat the analysis for the multiple scales series for general t, not just t = 0, by defining 'minimum error' and 'smallest term' in terms of the L~ norm, the minimum error simply is a as shown in Fig. 6(b). The optimally truncated series misses the high-frequency oscillations completely because ~(E) is a transcendental in I / e , which an e-power series cannot represent. However, for n = n o p t , the series mimics the true oscillatory manifold as closely as possible except for omitting the oscillations. When the initialization series has converged to within O(a), adding more terms only increases the error. Thus, (21) diverges as n ~ ~ for fixed E. Perturbation series, exponential smallness, and forced linear equations have been discussed in depth by Boyd (1991a). Exponential smallness is not a weird and special property of the slow manifold, but rather is a generic property of nonlocal solitary waves and a wide range of other important physical phenomena including metastable quantum molecules and ions, dendrite formation at solid-liquid interfaces, propagation in optical fibers and many others. Boyd (1987, 1989a, 1991a, 1995a) gave references to the large and rapidly increasing literature on exponential smallness in the reciprocal of parameter. Vautard and Legras (1986) gave a simple argument which specifies the rate of divergence: the terms in (21) must ultimately grow as n! if the forcing f(t) is singular anywhere in the complex t-plane except at infinity. This factorial rate of divergence is typical of asymptotic series.

8. Nonlinear slow manifolds

Up to now, the form of the forcing f(t) has been unspecified. In particular, there was nothing to restrict f(t) to be the product of the Rossby amplitudes only. Because f(t) was purely a symbol for the function on the right-hand side of the gravity mode equation, the analysis still applies when f(t) is adjusted to include the self-interaction of the gravitational modes.

J.P. Boyd/ Dynamicsof Atmospheres and Oceans22 (1995) 49-75

63

For example, in the LK Quintet, the Rossby amplitude W(t) depends weakly on x(0) and z(0), the initial amplitudes of the gravity waves. Thus, an initialization of a given order, such as u(0) = ioJ n=0k ~o ] dr" (0)

(26)

where u = z(t) + ix(t), requires an iteration (Tribbia, 1984) because

f(t)-bW,[t;

,, U(O), V(O), W(O), u(0)]

(27)

In other words, u(0) appears on both sides of (26) in a fully nonlinear treatment. However, the numerical evidence strongly suggests that the dependence of the Rossby mode W(t) on u(0) is very weak (Tribbia, 1984; Boyd, 1994). Similarly, the rapid convergence of all reasonable iteration schemes for operational forecasting models implies that the gravity waves have little effect on f(t) even for models with millions of degrees of freedom. We must confess, however, that we cannot yet prove that the weak feedback of the gravitational modes on the Rossby modes and on the gravitational modes themselves does not somehow change f(t) in a significant way. After all, Boyd (1991a) has convincingly demonstrated that the exponentially small high-frequency oscillations, precisely because they are so small, are very sensitive to changes in f(t). Do the linear slow manifolds really teach us anything at all about the nonlinear atmosphere? Fortunately, our analysis of linear manifolds depended only on the analytical properties of f(t), and not upon its precise numerical values. For example, the Vautard-Legras (1986) proof of the divergence of the initialization series (21) requires only that f(t) have singularities somewhere in the finite complex t-plane. That is, the series diverges unless f(t) is an entire function and perhaps even then, too. This argument falls short of a proof of divergence for a fully nonlinear model because we cannot prove that f(t) is not an entire function. However, it is believed that the atmosphere is chaotic. Chang et al. (1983) and Takaoka (1989) have shown that the solutions to certain simple chaotic flows have an infinite number of singularities distributed on fractal curves in the complex t-plane. In other words, these functions are the exact opposite of entire with lots of poles instead of none. The second important analytical property is that the Fourier transform F(w) must vanish at the seiche frequency oJ if u(t) is to lack high-frequency oscillations. Otherwise, we can use the regularization procedure of Section 5 to add and subtract the product of F(~o) and the regularizing function p(t) to demonstrate the existence of oscillations of apparent frequency ~o in the forced solution. In principle, there could be a feedback between the gravity and Rossby waves such that in a fully nonlinear model, the feedback deforms f(t) so that F(og) is compelled to be zero. However, no rationale for this zero-making feedback has ever been offered. Even when F ( w ) = 0, however, the initialization series diverges! As explained by Boyd (1991a), the analyticity of f(t) controls the convergence or divergence of

64

J.P. Boyd / Dynamics of Atmospheres and Oceans 22 (1995) 49-75

the series. Multiple scales series are almost always divergent. It is difficult with hindsight to understand why meteorologists ever thought it likely that (21) converged. Nonlinearity also muddles the distinction between seiche and pseudoseiche. There are exceptions; in Boyd (1994), the slow manifold is periodic in time whereas the general solution is doubly periodic with an O(e) frequency and an O(~o) frequency. The seiche is the sum of the components oscillating at the second, O(w) frequency. When the flow is chaotic, distinguishing seiche and pseudoseiche is much harder. The frequency power spectra will be broadband, and it is not possible to look at a single such spectrum and unambigously assign part to a seiche. However, the arguments offered above suggest that with or without an initialization, the power spectra will resemble that of Fig. 3 with a large spike at w---w. For any given solution of a numerical model, it is possible in principle to identify the seiche-like portion of the spectra, on a mode-by-mode basis, by subtracting that part of u(t) which is directly excited by F(w), the spectra of the sum of the nonlinear terms, f ( t ) ; the seiche is the remainder after the subtraction. Such a decomposition, as for the linear problem, is useful only in a theoretical sense. The practical problem is simply to minimize the high-frequency wiggles. We contend that even for fully nonlinear models the distinction between seiche and pseudoseiche is still a useful theoretical concept to allow us to appreciate what an initialization scheme can and cannot do.

9. S o m e c o m m e n t s on various slow m a n i f o l d s

9.1. The bandlimited slow manifold A function u(t) is 'bandlimited' with bandwidth l-I if its Fourier transform U(w) has the property

U(w) = 0 for all I wl > 12

(28)

Warn observed a decade ago, in an unpublished note, that the bandlimited manifold was the only definition which could be unambiguously labeled slow. In the absence of condition (28), u(t) contains components of arbitrarily high frequency. The power spectra of observations and numerical models and analytical solutions like those of Lorenz and Krishnamurthy (1987) and Boyd (1994) all suggest that it is impossible to adjust the initial values of gravity modes so as to produce a bandlimited solution. Another independent body of evidence which argues for the same conclusion is based on the Fourier integral derivation of (21). When 12 < w , the geometric series (22) is used only within its radius of convergence, as the Fourier integral is now evaluated only on the finite interval w ~ [ - 1 2 , 12]. The integrals (24) are bounded by 212 n max(I F(w) l), implying that the initialization series (21) must converge if f(t) is bandlimited. The observed

J.P. Boyd / Dynamics of Atmospheres and Oceans 22 (1995) 49-75

65

divergence of the initialization scheme as observed by Lorenz (1980) and others is then evidence that the solution is not bandlimited. Nevertheless, the goal of 'dynamic initialization' (Miyakoda and Moyer, 1968; Nitta and Hovermale, 1969; Daley, 1991) and also 'time filtering' (Lynch, 1985; Lynch and Huang, 1992) is to adjust the flow onto bandlimited slow manifold. Both families of algorithms strive to strongly damp or reduce high frequencies without making any distinction between free and forced gravity waves or between seiche and pseudoseiche. Both are as effective as any other strategy, in spite of the fact that the bandlimited manifold technically does not exist, because the irreducible minimum of high-frequency motion is exponentially small in l / e . To put it another way, it seems likely that an optimally initialized flow fails to satisfy the bandlimiting condition (28) because

U(w) ~ e x p ( - c o n s t a n t / E ) 4:0

I w [ >t 1~

(29)

that is, the slow manifold is bandlimited to within exponentially small components of the flow. To the extent that these small high-frequency components cannot be neglected, we are forced away from the simplicity of a manifold that is truly slow in the sense of the complete absence of high-frequency components to a mostly-slow-but-a-littie-fast manifold.

9.2. The pseudoseiche-free manifold This manifold is not an algorithm but an attitude. For example, Lorenz and Krishnamurthy (1987, p. 2943) asserted: 'for t > 10, there are unmistakable oscillations of apparent period 27r, which must be gravity waves and which imply that there is no universal slow manifold'. This assertion is true only if the slow manifold is defined in the restrictive sense of a manifold free of pseudoseiches as well as free gravity waves. This definition is defendable, as the primary aim of initialization schemes has always been to suppress high-frequency oscillations, and pseudoseiches have them. However, the nonexistence of the pseudoseiche-free slow manifold, though true, is only half the truth. The other half is there is a well-defined manifold which minimizes high-frequency gravity mode activity.

9.3. The gravity mode servitude manifold This definition is based on the notion that if the seiches are suppressed so that all gravity wave activity is forced directly by the slow modes, then the fast modes will be analytic functions of the slow modes. That is, if x(t) and z(t) are the fast modes and U, V, and W are the slow modes, then 'gravity mode servitude' requires

x=f(U, V, W),

z=g(U, V, W)

(30)

The key restriction is that f and g are independent of time so that the fast modes vary with t only through the time variations in U, V, and W.

66

J.P. Boyd / Dynamics of Atmospheres and Oceans 22 (1995)49-75

Kopell (1985), Vautard and Legras (1986) and Jacobs (1991) have all proved the existence of a manifold satisfying the conditions (30) for various simplified models in narrow parameter ranges. However, the analyticity condition is only a proxy. No one has yet offered a rigorous proof that (30) unfailingly suppresses high-frequency oscillations. Instead, Lorenz (1992) has shown that trajectories on Jacobs' manifold eventually escape from the convergence domain of Jacobs' power series, fold back upon themselves and develop high-frequency corrugations. Our comment is a second counterexample to the proposition that enslaving the fast modes to the slow as in (28) is a sure-fire wiggle-suppressor. For the inviscid LK Quintet (1), Lorenz and Krishnamurthy (1987) have given the analytic solution U = e sech(et) ~ t = ( l / E ) sech-l(U/E)

(31)

As there is a unique, one-to-one relationship between t and the slow mode U(t), it follows that, whatever the dependence of the fast modes on time, they will necessarily be analytic functions of U. As shown by Lorenz and Krishnamurthy (1987), however, the gravitational modes of the LK Quintet are pseudoseiches. Analyticity in the slow variables does not guarantee freedom from high-frequency oscillations. 9.4. The bounded derivative and seiche-free multiple scales manifolds

The definitions of these two manifolds are seemingly very different. The jth-order bounded derivative manifold is defined by the condition that for each gravity mode , =0 diu dt j t=0

(32)

(Machenhauer, 1977; Kreiss, 1980; Browning et al., 1980). This was dubbed the 'superbalance' condition by Lorenz (1980). This condition can be relaxed by replacing the right-hand side of (32) by O(e0, as is sometimes convenient. Kopell (1985) showed the equivalence of these two versions of the bounded derivative condition. The multiple scales manifold is derived by assuming the forcing f ( t ) in (4) is a function of a slow time variable T = et, rewriting (4) in terms of the slow variable, and then making a power series expansion in e (Baer and Tribbia, 1977; Tribbia, 1984). Despite these seemingly very different defining conditions, these two manifolds are often equivalent to one another, as shown by Kasahara (1982) and Daley (1991, p. 289). The reason is that (32) is implicitly an assumption that the desired, seiche-free solution is varying slowly in time so that its jth derivative at the origin will be small in comparison with unity. Our first comment is that the range of this equivalence of the two manifolds is still very confused. It is easy to prove that (32) must, like the multiple scales series, give (21) for the linear problem (1). However, for various nonlinear problems, Lorenz (1986), Warn and Menard (1986) and Curry et al. (1994) have found that

J.P. Boyd / Dynamics of Atmospheres and Oceans 22 (1995) 49-75

67

the superbalanced manifolds either converge or diverge slowly. In contrast, the multiple scales series diverges as j! where j is the order, as confirmed in Lorenz (1986). The reason for the dramatic differences between these two (usually) closely related methods is not known. Our second comment is that the multiple scales series defines a 'superasymptotic' approximation in the sense of Berry and Howls (1990) and Berry (1991). That is to say, if the multiple scales series is truncated, for a given fixed e, with the smallest term in (21) (the minimum of the appropriate curve in Fig. 5, for example), then the error is not a power of e, but rather an exponential function of 1/E. As illustrated by Boyd (1995a,b), an optimally truncated perturbation series converges to everything but the exponentially small part of the solution--in this case, to everything except the pseudoseiches. For the purposes of a wiggle-free initialization, a superasymptotic multiple scales approximation is as optimal as any of the alternatives. Our third comment is that it is possible to go beyond the superasymptotic approximation to a 'hyperasymptotic' approximation that explicitly calculates the exponentially small terms, too. Although this is useless for practical initialization-we want to kill the pseudoseiches, not compute them--hyperasymptotics are reassuring. Boyd (1995a) has catalogued a large number of hyperasymptotic schemes including the double Pad6 method (Reinhardt, 1982), complex-plane matched asymptotics (Segur and Kruskal, 1987), and the transform-and-restart perturbative method (Boyd, 1995b). Hyperasymptotics reiterates the point that there is a slow manifold, albeit a wiggly one.

9.5. Differential-algebraic manifolds The mth-order differential-algebraic manifold, which we shall abbreviate as DA(m), is obtained by replacing the differential equations for the time-evolution of the gravitational modes by algebraic equations, such as those obtained by truncating the Baer-Tribbia multiple scales series at mth order. Quasi-geostrophy, which filters gravity waves entirely, is a kind of DA(0); balanced systems are a species of DA(1); and the superbalanced manifolds computed by Lorenz (1980) and Warn and Menard (1986) are DA(m) of large m. For the linear problem (3), the first few DA(m) are u(O) = 0 u-

u=

f(t) iw

DA(O) DA(1)

f(t) df/dt iw + 0.)2

(33)

DA(2)

Daley (1980) made a thorough comparison between DA(0), DA(1), and the usual semi-implicit time-marching for the integration of a baroclinic forecast model, initialized with real data. Our comment is that there is a close analogy between such schemes and the

J.P. Boyd/ Dynamics of Atmospheres and Oceans 22 (1995) 49-75

68

so-called 'nonlinear Galerkin' method of Foias et al. (1988), Marion and Temam (1990) and Jauberteau et al. (1990). The key idea of the nonlinear Galerkin method is the observation that in many dissipative systems, the flow rapidly evolves to a narrow region in phase space, the so-called 'inertial manifold'. In the words of Jauberteau et al.: 'An approximate inertial manifold is a smooth finite-dimensional manifold which contains the attractor in a thin neighborhood and attracts all orbits in this neighborhood'. In chaotic flow, the dimension of the attractor is nonintegral and the flow is never on but only near the attractor for finite times. Computing the inertial manifold (whose dimension is an integer) is much easier than calculating the fractal attractor, and just as useful for reducing the dimensionality of the system for large times. The goal of the nonlinear Galerkin family of algorithms is to first approximate the inertial manifold and then approximate the flow itself by motion confined to the inertial manifold. The lowest-order nonlinear Galerkin method expands the flow in a series of standard basis functions, such as the terms of a Fourier series, and then truncates to N = dim(Im), the dimension of the inertial manifold. This tactic, which is analogous to a DA(0) approach to geophysical flows, is subject to Machenhauer's (1977) criticism: the nonlinear interaction of the modes included in the DA(0) model in fact forces activity in the modes outside the truncation. In other words, the dance of the slow modes casts a shadow on the fast modes, and accuracy is much improved by incorporating this shadow into the answer. The most popular nonlinear Galerkin algorithm calculates this 'shadow' in the same way as Daley's (1980) Scheme B: the DA(0) approximation

u(t) = 0

(34)

is replaced by the DA(1) approximation

u( t ) = f ( t) / ( ioJ) (35) When there are many 'fast' modes so that u(t) is a vector, the division is to be interpreted as the element-by-element division of the jth element of the vector f by the frequency of the jth mode. The major difference between differential-algebraic systems for weather forecasting and previous applications of the nonlinear Galerkin method is that in the latter, 'fast' and 'slow' refer to rates of dissipation rather than to frequencies of oscillation. The divisor in (33) is not i times frequency of oscillation of the mode but rather is a real number which is the mode's damping rate, the eigenvalue of the Laplacian operator (or a similar dissipative operator). However, the algorithms are the same. The nonlinear Galerkin method was independently developed, almost a decade after the work of Daley, Machenhauer and other meteorologists, by Foias et al. (1988) and Marion and Temam (1989). Debussche and Temam (1991, 1994) have asserted the close relationship between the slow manifold and an inertial manifold, but this relationship has to be stated carefully. The inertial manifold existence proofs always rely on dissipation to pull the flow onto the manifold. In contrast, the proposition that the flow of the atmosphere lies on a manifold containing only low frequencies is a statement about relative oscillation periods, not about relative damping time scales.

J.P. Boyd / Dynamics of Atmospheres and Oceans 22 (1995)49-75

69

Nevertheless, the damping of gravity waves, rapid compared with Rossby modes, is what keeps the atmosphere in a state which is mostly-slow-but-a-little-fast (Jacobs, 1994). It seems probable that there is an inertial manifold for the atmosphere, and that the slow manifold includes this. However, the slow manifold is a stronger concept. Both the slow manifold and inertial manifolds in general assume that the energy in strongly damped m o d e s - sufficiently far into the 'dissipation range' of turbulence t h e o r y - - i s negligible. The meteorologist's slow manifold also ignores those modes whose damping time scales are large compared to the oscillation time scale and are therefore not in the dissipation range. A proof that the gravitational modes are less than some small fraction 6 of the rotational modes (in energy) is still lacking.

9. 6. The optimal control manifold As the goal of initialization is to minimize something, optimal control methods are the most natural formulation of algorithms for adjusting the flow onto the slow manifold. Such methods have been little used, however, partly because of cost and partly because optimal control became an important branch of applied mathematics (as opposed to engineering) only recently. This manifold is defined by choosing the initial conditions on the gravitational modes to minimize a chosen 'cost' function over the forecast interval. This is a more natural definition than the bounded derivative condition, which is expressed entirely in terms of behavior at t = 0. However, the cost function is a mathematical proxy to control the gravity waves. The manifold is only useful to the degree that the cost condition is high for seiches and low for Rossby waves and pseudoseiches. To illustrate the idea, we used the LK Quintet without forcing or damping and defined the cost function to be

C( xo,Zo) - fo'FX2t + z2t dt

(36)

where x 0 and z 0 are the initial values for the gravitational modes and x(t) and z ( t ) are defined by integrating the five coupled differential equations for the chosen values of x 0 and z0 (plus the three initial conditions for the Rossby waves, which are fixed during the variational calculation). Optimal control seems not to have been previously applied to initialization. However, Courtier and Talagrand (1990) added a term of the form of (36) to control gravity wave oscillations during variational assimilation. They noted, rightly, that initialization is in some sense merely a special case of assimilation, and therefore the cost function (36) should be useful in controlling gravity waves both with and without the additional data-fitting terms that are added to the cost function for assimilation. Lorenc (1988) has used a similar strategy. Courtier and Talagrand (1990) found that the cost function (36) is considerably better than no initialization at all. We have confirmed this for the LK Quintet. Unfortunately, Fig. 7 shows that the oscillations in the gravitational modes are considerably larger than for a low-order multiple scales initialization. The problem

70

J.P. Boyd/ Dynamics of Atmospheres and Oceans 22 (1995) 49-75

.001 8e-4 6~4 4e-4 28-4

-2e-,

C~r~rison ofz

i

-4e-4

-6~4

Fig. 7. A comparison of the gravitymode z(t) for the LK Quintet (b' = b, a = F = 0) and two different sets of initial conditions. U(0)= e, V(0)= 0, and W(0)= - • on the interval t e [0,120] with • = 1/6. Continuous curve (large amplitude): optimal control manifold, x(0)= 6.23E-5, z(0)=-0.00463. Dashed curve: first-order multiple scales initialization, x(0)= 0, z(0)=-0.00527. To compare the high-frequency oscillations more easily, only a subinterval in I is illustrated. is that, like the algorithmic goals of the other flavors of the slow manifold, the cost function (36) is only an indirect measure of high-frequency activity. What we really want might be coarsely summarized as ' n o wiggles'. Unfortunately, this verbal definition must always be replaced by a mathematical proxy. The minimization of (36) is obviously a rather poor proxy, though better than nothing. A much more effective cost function would be frequency selective. We could, for example, compute the temporal spectra of a forecast and minimize the power spectra in some small neighborhood of each of the gravity mode frequencies. We have not put this improved cost function to the test because of the complexity and expense. Nevertheless, in spite of their cost, optimization methods are becoming more and more important in meteorology and oceanography. As noted by Courtier and Talagrand, the minimization of (36) is closely related to four-dimensional data assimilation (Daley, 1991). Perhaps, with a more skillfully chosen cost function, the optimal control manifold will be useful. 9.7. The fuzzy manifold The exponentially fuzzy slow manifold is an unpublished idea of Warn's. He observed a decade ago that, as the Baer-Tribbia series is only asymptotic, we cannot use its sum to specify a unique slow manifold. However, because the error in truncating the asymptotic series at optimum order is an exponential function of the reciprocal of the Rossby number, it follows that we can use the standard initialization procedures to compute flows in which the gravity wave oscillations

J.P. Boyd / Dynamics of Atmospheres and Oceans 22 (1995) 49-75

71

have been reduced to O[exp(- constant/e)]. (In other words, he advocated using a superasymptotic approximation as an optimum initialization, though Berry coined the term 'superasymptotic' years later.) The true slow manifold, if it exists, is replaced by a thickened manifold of exponentially small thickness. For operational forecasting, initializing anywhere within this fuzzy manifold is satisfactory. Indeed, current models normally use schemes equivalent to only the first- or second-order Baer-Tribbia series. Diabatic effects and sources of gravitational noise (e.g. mesoscale phenomena such as squall lines and thunderstorms) are more serious than the error in a low-order multiple scales series (Daley, 1991; Errico, 1984, 1989a,b). Because it appears that pseudoseiches are as inevitable as death and the tides, the fuzzy manifold is a useful concept. Any initialization scheme which adjusts the flow to within this thickened manifold is, in a practical sense, optimal, and as good as any alternative algorithm that also lands within the fuzzy manifold. However, the algorithms we have described above are precise and definite; once we have specified the scheme and the order, all ambiguity is removed. It is the goal of (practical) initialization which is fuzzy, not the algorithms.

10. Open problems Although rigorous proof is lacking, it seems clear that: (1) a pseudoseiche-free slow manifold, that is, one which is completely free of oscillations of frequency w ~ w where w is a gravity mode frequency, does not exist. (2) A seiche-free slow manifold, that is, one which has minimal (but non-zero) oscillations, does exist. (3) The amplitude of the irreducible minimum of gravity wave activity is O[exp(- constant/e)] where E << 1 is the Rossby number. (4) Existing initialization schemes are generally divergent. However, because the minimum error of the divergent, asymptotic series is exponentially small in 1/E, good results can be obtained by applying the scheme at finite order. (5) In a practical sense, neither the inescapable minimum of gravity wave activity nor the resulting divergence of the initialization schemes is a problem. Rather, the important difficulties in initialization are other topics such as cumulus convection and other diabatic heating effects, tropical initialization where the separation between fast and slow modes is less extreme, mesoscale generation of gravity waves, and so on (Errico, 1984, 1989a,b; Ko et al., 1989a,b; Daley, 1991). Nevertheless, a number of important mathematical questions are unsolved. First, is it possible to prove rigorously (1) that a pseudoseiche-free slow manifold does not exist and (2) an oscillation-minimizing slow-with-a-little-fast manifold does exist? Closely related problems have been solved in the theory of nonlocal solitary waves. Vuillermot (1987, 1988) and Kichenassamy (1991) have proved nonexistence-of-wiggle-free solutions for the '04 field theory' equation and a model of dendrite formation in crystal growth. Hunter and Scheurle (1988), Beale (1991) and Sun (1991) all have proved the existence of generalized solitary waves which

72

J.P. Boyd / Dynamics of Atmospheres and Oceans 22 (1995) 49-75

take the form of a classical soliton combined with a small-amplitude oscillation which fills all of space. It seems likely that their method could be extended to the slow manifold. Second, what is the relationship of the gravity mode servitude manifold to the other definitions of slow manifold? Does the analyticity condition reproduce all of say, the bounded derivative manifold, or perhaps only a portion, as Lorenz's (1992) numerical calculations suggest? Does the analyticity condition minimize highfrequency oscillations or do other manifolds have significantly smaller oscillations than the gravity mode servitude manifold? Third, how rapidly does the slow manifold attract perturbed flows with nonzero seiches? Jacobs (1990) has shown that in a simple damped model, the highfrequency free gravity waves decay rapidly, leaving a flow dominated by slow modes. The existence of some sort of mechanism that evolves the general solution to the special family of solutions which is the slow manifold is strongly suggested by observations. However, this relaxation mechanism has been little studied. Fourth, is it possible to prove, as strongly suggested by numerical evidence and the arguments given above and by Vautard and Legras (1986) and T. Warn (unpublished data, 1983), that the multiple scales initialization series is divergent? As explained above, divergence of the series is equivalent to the statement that the nonlinear forcing of the gravity mode, f ( t ) , is not an entire function. We offer two conjectures which, if proven, would remove this uncertainty: (1) chaotic motion cannot be represented by entire functions: (2) modal amplitudes in atmospheric models have fractal natural boundaries in the complex t-plane. Fifth, when is the bounded derivative initialization equivalent to the multiple scales method and when and why are they different? As noted above, l_x)renz (1986), Warn and Menard (1986) and Curry et al. (1994) have all found problems for which the bounded derivative (superbalance) manifolds somehow escape the expected divergence. Yet one can prove that for the forced, linear problem (3), the multiple scales and bounded derivative manifolds are identical, order-by-order (Daley, 1991, p. 289). The argument is a simple one, based on repeatedly differentiating (3) to express high derivatives, as needed for (32), in terms of low-order derivatives. The words of Ko (1987) are still true: 'the relationships among these.., procedures is not clear'. Lastly, the frequency spectra of the slow modes for simple flows, such as the analytic solution to the inviscid LK Quintet, decay exponentially fast with frequency w. This in turn implies that the irreducible minimum of gravity wave activity is exponentially small, too, as already asserted as point (3) above. However, atmospheric flows are believed to be chaotic, and there is no proof that the Fourier transform of a chaotic atmospheric solution must decay exponentially with w. It is possible that the decay is only algebraic, as is known to be true of the spatial spectrum, in which case the minimal gravity wave activity is probably proportional to E-a for some constant d > 0. In the atmosphere, how rapidly does the frequency spectrum actually decay with frequency? Clearly, the slow manifold is still full of challenges.

J.P. Boyd / Dynamics of Atmospheres and Oceans 22 (1995) 49-75

73

Acknowledgments I t h a n k J o s e p h T r i b b i a f o r h e l p f u l c o m m e n t s , e s p e c i a l l y for p o i n t i n g t h e a n a l o g y b e t w e e n t h e m a t h e m a t i c s o f h y s t e r e s i s in t h e L o r e n z - K r i s h n a m u r t h y model and t h e m a t h e m a t i c s o f w e a k l y n o n l o c a l s o l i t a r y w a v e s . I a m a p p r e c i a t i v e to T h o m a s W a r n f o r e x p l a i n i n g his e a r l y u n p u b l i s h e d w o r k a n d also f o r r e m i n d i n g m e o f t h e e s s e n t i a l r o l e o f d i s s i p a t i o n in i n e r t i a l m a n i f o l d t h e o r y . I also t h a n k S t a n J a c o b s for s u p p l y i n g m e w i t h his w o r k in m e d i a s res. T h i s w o r k was s u p p o r t e d by t h e N S F t h r o u g h G r a n t s O C E 8 8 1 2 3 0 0 , DMS8716766, E C S 9 0 1 2 2 6 3 a n d O C E 9 1 1 9 4 5 9 , a n d by t h e D e p a r t m e n t o f E n e r g y t h r o u g h G r a n t K C 0 7 0 1 0 1 .

References Baer, F., 1977. Adjustment of initial conditions required to suppress gravity oscillations in non-linear flows. Contrib. Atmos. Phys., 50: 350-366. Baer, F. and Tribbia, J.J., 1977. On complete filtering of gravity modes through nonlinear initialization. Mon. Weather Rev., 105: 1536-1539. Beale, J.T., 1991. Exact solitary water waves with capillary ripples at infinity. Commun. Pure Appl. Math., 44: 211-257. Berry, M.V., 1991. Asymptotics, superasymptotics, hyperasymptotics. In: H. Segur, S. Tanveer and H. Levine (Editors), Asymptotics Beyond All Orders. Plenum, New York, pp. 1-14. Berry, M.V. and Howls, C.J., 1990. Hyperasymptotics. Proc. R. Soc. London, Set. A, 430: 653-668. Boyd, J.P., 1987. Generalized solitary and cnoidal waves. In: J.J. Tribbia and R.C. Madden (Editors), Proceedings of the NCAR Colloquium on Low Frequency Variability in the Atmosphere. National Center for Atmospheric Research, Boulder, CO, pp. 715-722. Boyd, J.P., 1989. New directions in solitons and nonlinear periodic waves: polycnoidal waves, imbricated solitons, weakly non-local solitary waves and numerical boundary value algorithms. In: T.-Y. Wu and J.W. Hutchinson (Editors), Advances in Applied Mechanics, 27. Academic Press, New York, pp. 1-82. Boyd, J.P., 1991a. A comparison of analytical and numerical methods for the reduced wave equation with multiple spatial scales. Appl. Numer. Math., 7: 453-479. Boyd, J.P., 1991b. Weakly non-local solitary waves. In: A.R. Osborne and L. Bergamasco (Editors), Topics in Nonlinear Ocean Physics: Proceedings of the Fermi School, Course 109. North-Holland, Amsterdam, pp. 527-556. Boyd, J.P., 1994. The slow manifold of a five-mode model. J. Atmos. Sci., 51: 1057-1064. Boyd, J.P., 1995a. Weakly nonlocal solitary waves and other exponentially small phenomena. (Monograph.) in preparation. Boyd, J.P., 1995b. A hyperasymptotic perturbative method for computing the radiation coefficient for weakly nonlocal solitary waves. J. Comput. Phys., in press. Browning, G., Kasahara, A. and Kreiss, H., 1980. Initialization of the primitive equations by the bounded derivative method. J. Atmos. Sci., 37: 1424-1436. Chang, Y.F., Greene, J.M., Tabor, M. and Weiss, J., 1983. The analytic structure of dynamical systems and self-similar natural boundaries. Physica, DS: 183-201. Courtier, P., 1987. Application du contr61e optimal ~ la pr6vision num6rique en M6t6orologie. Th~se de doctorat de l'Universit6 Pierre-et-Marie-Curie, Paris. Courtier, P. and Talagrand, O., 1990. Variational assimilation of meteorological observations with the direct and adjoint shallow-water equations. Tellus, 42A: 531-549.

74

J.P. Boyd / Dynamics of Atmospheres and Oceans 22 (1995) 49-75

Curry, J.H., Haupt, S.E. and Limber, M.N., 1994. Low order models, initialization, and the slow manifold. Tellus, in press. Daley, R., 1980. The development of efficient time integration schemes using model normal modes. Mon. Weather Rev., 108: 100-110. Daley, R., 1991. Atmospheric Data Analysis. Cambridge University Press, Cambridge, 457 pp. Debussche, A. and Temam, R., 1991. Inertial manifolds and slow manifolds. Appl. Math. Lett., 4: 73-76. Debussche, A. and Temam, R., 1994. Inertial manifolds and slow manifolds in meteorology. Differential and Integral Equations, in press. Errico, R.M., 1984. The dynamic balance of a general circulation model. Mon. Weather Rev., 112: 2439-2454. Errico, R.M., 1989a. The degree of Machenhauer balance in a climate model. Mon. Weather Rev., 112: 2723-2733. Errico, R.M., 1989b. Forcing of gravitational normal modes by condensational heating. Mon. Weather Rev., 112: 2734-2752. Foias, C., Jolly, M.S., Kevrekidis, I.G., Sell, G.R. and Titi, E.S., 1988. On the computation of inertial manifolds. Phys. Lett. A, 131: 433-436. Hunter, J.K. and Scheurle, J., 1988. Existence of perturbed solitary wave solutions to a model equation for water waves. Physica, 32D: 253-268. Jacobs, S.J., 1991. On the existence of a slow manifold in a model system of equations. J. Atmos. Sci., 48: 893-901. Jacobs, S.J., 1994. The decay rate of gravity waves. To appear. Jauberteau, F., Rosier, C. and Temam, R., 1990. The nonlinear Galerkin method in computational fluid dynamics. Appl. Numer. Math., 6: 361-370. Kasahara, A., 1977. Numerical integration of the global barotropic primitive equations with Hough harmonic expansions. J. Atmos. Sci., 34: 687-701. Kasahara, A., 1978. Further studies on a spectral model of the global barotropic primitive equations with Hough harmonic expansions. J. Atmos. Sci., 35: 2043-2051. Kasahara, A., 1982. Nonlinear normal mode initialization and the bounded derivative method. Rev. Geophys., 19: 450-468. Kichenassamy, S., 1991. Breather solutions of the nonlinear wave equation. Commun. Pure Appl. Math., 44: 789-818. Ko, S.-D., 1987. Vertical modes and energetics of gravitational and rotational modes in a multilevel global spectral model. Ph.D. Thesis, University of Michigan, p. 5. Ko, S.-D., Tribbia, J.J. and Boyd, J.P., 1989a. Energetics analysis of a multilevel global spectral model. Part I: Balanced energy and transient energy. Mon. Weather Rev., 117: 1941-1953. Ko, S.-D., Tribbia, J.J. and Boyd, J.P., 1989b. Energetics analysis of a multilevel global spectral model. Part II: Zonal and meridional gravitational energy. Mon. Weather Rev., 117: 1954-1964. Kopell, N., 1985. Invariant manifolds and the initialization problem for some atmospheric equations. Physica, D14: 203-215. Kreiss, H., 1980. Problems with different time scales for partial differential equations. Commun. Pure Appl. Math., 33: 399-440. Leith, C.E., 1980. Nonlinear normal mode initialization and quasigeostrophic theory. J. Atmos. Sci., 37: 958-968. Lorenc, A.C., 1988. Optimal nonlinear objective analysis. Q. J. R. Meteorol. Soc., 112: 205-240. Lorenz, E.N., 1980. Attractor sets and quasi-geostrophic equilibrium. J. Atmos. Sci., 37: 1685-1699. Lorenz, E.N., 1986. On the existence of a slow manifold. J. Atmos. Sci., 43: 1547-1557. Lorenz, E.N., 1992. The slow manifold--what is it? J. Atmos. Sci., 49: 2449-2451. Lorenz, E.N. and Krishnamurthy, V., 1987. On the nonexistence of a slow manifold. J. Atmos. Sci., 44: 2940-2950. Lynch, P., 1985. Initialization using Laplace transforms. Q. J. R. Meteorol. Soc., 111: 243-258. Lynch, P., 1992. Richardson's barotropic forecast: a reappraisal. Bull. Am. Meteorol. Soc., 73(1): 35-47. Lynch, P. and Huang, X., 1992. Initialization of the HIRILAM model using a digital filter. Mon. Weather Rev., 120: 1019-1034.

J.P. Boyd / Dynamics of Atmospheres and Oceans 22 (1995) 49-75

75

Machenhauer, B., 1977. On the dynamics of gravity oscillations in a shallow water model, with application to normal mode initialization. Beitr. Phys. Atmos., 50: 253-271. Marion, M. and Temam, R., 1989. Nonlinear Galerkin methods. SIAM J. Numer. Anal., 26: 1139-1157. Marion, M. and Temam, R., 1990. Nonlinear Galerkin methods: the finite elements case. Numer. Math., 57: 205-226. Miyakoda, K. and Moyer, R., 1968. A method of initialization for dynamical weather forecasting. Tellus, 20: 115-128. Nitta, T. and Hovermale, J., 1969. A technique of objective analysis and initialization for the primitive forecast equations. Mon. Weather Rev., 97: 652-658. Platzman, G.W., 1969. A retrospective view of Richardson's book on numerical weather prediction. Bull. Am. Meteorol. Soc., 48: 514-550. Reinhardt, W.P., 1982. Pad6 summation for the real and imaginary parts of atomic Stark eigenvalues. Int. J. Quantum Chem., 21: 133-146. Richardson, L.F., 1922. Weather Prediction by Numerical Processes. Cambridge University Press. Reprinted in 1965, Dover, New York, 236 pp. Segur, H. and Kruskal, M.D., 1987. Nonexistence of small amplitude breather solutions in t~ 4 theory. Phys. Rev. Lett., 58: 747-750. Sun, S.M., 1991. Existence of a generalized solitary wave solution for water waves with positive Bond number less than 1/3. J. Math. Anal. Appl., 156: 471-504. Takaoka, M., 1989. Pole distribution and steady pulse solution of the fifth order Korteweg-deVries equation. J. Phys. Soc. Jpn., 58: 73-81. (Addendum, 58: 3028.) Tribbia, J.J., 1984. A simple scheme for high-order nonlinear normal mode initialization. Mon. Weather Rev., 112: 278-284. Vautard, R. and Legras, B., 1986. Invariant manifolds, quasi-geostrophy and initialization. J. Atmos. Sci., 43: 565-584. Vuillermot, P.A., 1987. Non existence of spatially localized free vibrations for a class of nonlinear wave equations. Comment. Math. Helvetici, 64: 573-586. Vuillermot, P., 1988. Varieties lisses a certains systems dynamiques et solitons quasiperiodiques pour les equations de Klein-Gordon non-lineaires sur R e. C.R. Acad. Sci. Paris, 307, Ser. I, 639-642. Warn, T. and Menard, R., 1986. Nonlinear balance and gravity-inertial wave saturation in a simple atmospheric model. Tellus, 38A: 285-294.