THE UNIVERSITY OF MICHIGAN COLLEGE OF ENGINEERING Department of Meteorology and Oceanography Technical Report No, 4 AN ASYMPTOTIC SOLUTION OF THE TIDAL EQUATIONS S. Jo Jacobs ORA Project 06372 Supported by: NATIONAL SCIENCE F'OUNDATION GRANT NO. GP-2561 WASHINGTON, D.C. administered through: OFFICE OF RESEARCH ADMINISTRATION ANN ARBOR April 1967

TABLE OF CONTENTS Page LIST OF FIGURES v ABSTRACT 1.1 INTRODUCT I ON 2 2o FORMULATION 3 3. THE RAY APPROXIMATION 6 4,o ASYMPTTIC SOLUTIONS NEAR THE CAUSTICS 19 5o DISCUSSION 24 APPENDIX A 26 APPENDIX B 29 REFERENCES 33 iii

LIST OF FIGURES Figure Page 1o Ray paths in (x,y) plane for Rossby mode, with arrow denoting direction of group velocity. 34 2. Ray paths in (x,y) plane for + gravity mode. 35 3. Caustics in (y,t) plane for + gravity mode for k = 1. 36 4. Pre-image of curve y = 1 for k = 1. 37 5. Ray paths near a caustic~ 38 v

ABSTRACT The Cauchy problem for the P-plane form of the tidal equations is solved for oscillatory initial data. The radius of deformation is assumed to be much less than the radius of the earth, and in accord with this assumption a ray approximation is employed. It is shown that, due to the rapid rate of propagation of inertio-gravity waves, the motion in its initial development tends toward geostrophic balance. However, the solution given by the ray approximation is singular on certain surfaces in space and time, the envelopes of the rays. A local boundary layer theory is employed to correct this deficiency. The existence of these caustics implies that the process of geostrophic adjustment is more complicated than hitherto imagined.

1. Introduction. One possible explanation for the well-known fact that gravity waves play an unimportant role in the large scale motion of the atmosphere is that these waves, because of their high rate of propagation, spread out rapidly from local sources, leaving the slower Rossby waves behind. This process is called geostrophic adjustment. In the standard treatment (Obukhov, 1949) a flat earth idealization of the tidal equations is solved as an initial value problem, and the motion indeed tends toward geostrophic balance as time - oo. Among the defects of this model is the neglect of curvature effects and the consequent omission of Rossby waves and of refraction of the gravity waves. Accordingly, there is some interest in a treatment in which the effect of the earth's curvature is not totally neglected. In this paper we consider the:-plane model and thus allow for the above mentioned omissions of Obukhov's study. A parameter N, the ratio of the radius of the earth to the radius of deformation, is assumed to be large, and an approximate solution is found on this basis. The approximation is similar to the ray method used by Keller (1958) for the treatment of diffraction problems, though here it is used to solve an initial value problem.

3 As in diffraction theory the ray approximation is invalid on envelopes of the rays. A local boundary layer theory is needed to correct this deficiency and is supplied in part. Near and on the envelopes the amplitudes of the dependent variables are large. This indicates that the process of geostrophic adjustment is considerably more complicated than would appear from a study of the flat earth model. 2. Formulation. Let ~ and 0 measure longitude and latitude on an earth of radius R rotating about a polar axis with angular velocity n, and let g be the gravitational constant. It is assumed that 02R/g is small so that ellipticity of geopotential surfaces can be neglected. We consider here a homogeneous fluid of uniform depth H with a free surface, but bear in mind that the theory applies equally well to inhomogeneous fluids provided H is replaced by an appropriate scale height. Let t be the surface elevation, u and v velocity components to the east and north, and let x and y be Mercator coordinates defined by x =, y = log ( s n ) (1) cos e

4 We ignore the tide-generating forces. Then, in accord with the usual approximations of tidal theory (Lamb, 1932, Chapt. 8) the equations of motion are - 2Qfv + =0, (2) + 2fu gm = 0, (3) a~t R ay +t m2[ (u/m) + (v/m)] = (4) at R x where f = sin e = tanh y, m = sec 0 = cosh y. (5) In the B-plane approximation used here we replace m by unity and f by y. Introducing this approximation, scaling the variables through r = hi*, (u,v) = h(g/H)/ (u*,v*), (6) (x,y) = (x*,y*), -1/2 t = R(gH) - t*, where h is a characteristic amplitude of the surface elevation, letting N = 2CR/(gH), (7) and omitting the asterisks, we obtain t - Nyv + = 0 (8) at + Nyu + = 0. (9) at ay~~~~~~~~~~~~~~9

5 and a u av + $- + 0. (10) at + ax+ = o These are to be solved subject to initial conditions and to the requirements that the dependent variables be bounded at Iyl = c and periodic in x with period 2r. It proves convenient to cast the problem in terms of v alone. Elimination of C between (8) and (10) and between (9) and (10) and of u between the resulting equations yields V2 3v 6v 22v t sv + N N-2y2 0. (11) Similarly, the initial data becomes v = V9 v_ = - Nyu - r, (12) v = Ny(X - N) + u + v YY where a tilda over a variable denotes its value at t = 0. We will suppose that N is large and exploit this in seeking an asymptotic solution. Also, we take as initial data u = v = 0, iNkx (13) = A(x,y)e where k is a constant of order unity and where A(x+27r.y) = e2Nk A(xy) (14) in agreement with the sentence following (10). Thus the

6 initial disturbance is in the form of a wave with variable amplitude and whose wavelength is short compared to the radius of the earth. 3. The ray approximation. We assume that the asymptotic expansion of v is of the form iNpR (xy,T) v = w (x,y,T,N)e + ~R ~ iNcp+(x,y,t) + w (x,y,t,N)e + iNc_ (x,y, t) + w (x,y,t,N)e, (15) where T = Nt, (16) where the first term of the sum is to satisfy a Y " " aS +,, a iNcpR (N2V2 + N4 - N 4y2)wR e =0. (17) ^T ^T3Sx J R and where each of the last two terms satisfies iNcp+ (V2 - + N - N2y2 w+ e 0. (18) at at3 x at + The motivation for this assumption is due in equal measure to known facts about long wave motions in the atmosphere and oceans and to what is now classical singular perturbation theory. We note first that there are two classes of long waves which are often observed in geophysical contexts, the inertio-gravity and the Rossby waves. In elementary treat

7 ments their dispersion relations are derived by assuming that the Coriolis parameter may be taken to be constant once a scalar equation for a single unknown has been obtained. In the present case this amounts to taking y to be constant in (11), whence the equation admits plane wave solutions of the form = e iN(px + qy - ct) where W3 - (p2 + q2 + y2)w3 - Nlp = 0. With an error of order N-1 the cubic has approximate solutions - r r(p2 + q2 + y2) the dispersion relation for Rossby waves, and D = ~(p2 + q2 + y2)Y the dispersion relation for inertio-gravity waves. This calculation implies that v should be written as v = v (x,y,r,N) + v (x,y,t,N) + v (x,y,t,oN) R + where the notation is obvious. Now the location of the large parameter N in (11) or in the amended form of (11) with T as time variable indicates that this is a singular perturbation problem, with either the time or spatial derivatives of v large, of order N. Hence we are led naturally to assuming the asymptotic expansion of v to be of the form

8 v+ = v+ (x,y,t,Ncp+(x,y,t),N), (19.a) VR = vR(x,yTNcpR(x,y,T),N) (19.b) *z R where I+ and ~R are additional unknowns. This prescription is due to Mahoney (1962) and is designed to insure that the asymptotic expansion be uniformly valid in space and time. The next step would be to substitute (19.a) into (11) and (19.b) into the amended form of (11), and to solve by expanding v+ and vR in a perturbation series in powers of hN1. It is reasonably obvious that if this procedure is carried out the lowest order equation for either v in (19) involves derivatives only with respect to Ncp and is translationally invariant with respect to this variable. Mahoney's advice is to pick p so that, if the solution involves a function of NFcp, where F is a functional of c F is a constant. These last two sentences lead immediately to equation (15), and we point out that assuming this form of solution is the first step in the geometrical optics approach to wave propagation problems. In both (17) and (18) we impose the condition that the coefficient of the highest power of N be equated to zero. This leads to cpx c = (cp2 + cp2 + y2)' (20) ~z (~ y

9 which is obeyed by cpR and -t = +( 2 + cp + Y2), (21) -~ t - ~x y which is obeyed by c+. These will be referred to as dispersion relations but, in contrast to the case of wave propagation in a constant medium, the phase is not a linear function of its arguments. To satisfy the initial conditions we require that eR = + = ^ = kx. (22) It can then be shown that (p+t = ~(k2 + y2)/, _ -1,~-, (23) ((+)tt = +y(k2 + y2) at t = 0, and these equations together with (12) and (13) imply that at t = 0 w + w + w = 0, (24) R +h - (wR) + i(CR) R + (w + w t = -A, (25) R r R R + - t y and N- ( (w) + iNF[2(p R) (WR) + (pR) WR]- (CR) 2 + R T R R tT R Rt R + (w + w_) + iN[2(k2 + y2)/(w - w ) + + - tt- + t -/ (26) + y(k2 + y2) (w - w ) ] - N (k2 + y2) ( + W ) = y(iN2kA + NA). Thus far no approximations have been made. Now, however, we expand each w in an ordinary perturbation series of

10 the form w = w() + N1 w(l) + N- w2 +..., (27) and substitute into (17), (18) and the initial conditions. It follows that w() obeys (p2 + cp2 + y2)w + (2Xcp - 1)w + 2c y w + x y T xT x yT y (28) + [cTV2p + 2V cp VT ]W = 0, where cp is cp and that both of w(o)satisfy R + 2(Vcp.Vw - tt + (V - - t) - ic /Pt )w = 0, (29) where p is p+. These will be referred to as the transport equations. Substitution of (27) into the initial conditions yields ~(o) w(o) (o) w + w +w =0, R + - " (0) (0) = (3) (o) ~(o) ( kyo w+ + w - i 2+Y2 A(xY), 0W -% = 0, (3=0) +w We consider first the Rossby mode, and as before omit subscripts and superscripts whenever possible. Our first task is to solve the dispersion relation (20). To this end, let ~(o) _ iky R ~x, -), (31. a) r(o) = ~(o) 1 iky w _ =F2 -2A(x, y). (31lb) We consider first the Rossby mode, and as before omit subscripts and superscripts whenever possible. Our first task is to solve the dispersion relation (20). To this end, let P = x q = cpy a = cp' (32) so that

11 O(p2 + q2 + y2) _ p = 0. (33) The associated characteristic system of differential equations is (Courant and Hilbert, 1962, Chapt. 2) x = 2p - 1, yV =2qo, r = p2 + q2 + y2, p =0, (34) q= -2 y, a~ =0, v = 2o(p2+q2), where v is a parameter. Solving (34) subject to X =?, y = p, T = 0, p = k, (35) q = 0, p = kx, at v = 0, where? and pJ are also parameters, we obtain the solution of (20) in parametric form, cP = = k\v - 2 sin(14kv), (36.a) ~R 4 k~ - R k2+PI2 4 k2 where k2 -12 2kv = + k2+2, Y = COS(k2+ 2) T = (k2+2)v. (36.b) We designate the family of curves generated by assigning constant values of A and P in (36.b) as rays. These are not orthogonal to surfaces of constant phase, instead, as may be verified, they are everywhere tangent to the local vector group velocity of the waves. We turn now to the problem of solving the transport equation (28). Now from (34), - =x + y + - _v v ax v ay v aT (37) = (2cp P-1)- + (2pcp ) (cp +cp2 + y2 x'yxT y Y x y t

12 so (28) is actually an ordinary differential equation, the differentiation in (37) being (p2 + p2 + y2) times x y differentiation following points which move with the group velocity. Furthermore, letting j = (x,y,T) (38) R a (3,8Lv) be the Jacobian of the transformation (36.b), we find after a short calculation that (JR)v = 2JR[ T2(p + 2Vcp (39) so that (28) becomes 2J w + (JR)V = 0 (40) with solution w = W(?,) [JR(, L)/JR] (41) Evaluating the Jacobian and using (31.a), we obtain (= wO) ikj. 2kv + w = w V — A (~, os R )= k2+2 A(7,) [cos(k2+) + + 8 ki2 sin( 1 (42) (k2+i2) 2 k (2+ ) 42) and the first term in the ray solution for the Rossby mode is (o) (o) iNNR VR R Since the treatment of the + and - gravity modes is identical, we exhibit the calculations only for the former. To solve the dispersion relation, let P = Tcpx q = py C (43) y T

13 and we write the dispersion relation as 1 (44) -(P2 + q2 + y2 - D2) = 0. (44) the characteristic system of differential equations is x = p y = q t = C, p = 0, (45) q = -y, = = - = = (y2 where s is a parameter, and these are to be solved subject to x = ~, y = 7, t = 0, p = k, 1/ (46) q = 0, = (k2+y2), cp = kx, at s = 0, where 4 and 77 are also parameters. The solution is C = = k - 277s - 772 sin 2s (47oa) CP = CP + 2 4 where x = i + ks, y = 7? cos s, t = (k2+72)2s, (47 b) and as before the rays, generated by assigning constant values to ~ and T7 in (47.b), are everywhere tangent to the group velocity. Since cp- x + cP y -cp (48) as = *x x + y ay - Tt at' the transport equation (29) is an ordinary differential equation. Also, letting =_ (x,y,t) (49) + a ( 7 s ) S be the Jacobian of the transformation (47.b), we find that

14 (J)s = J (V2P - _ ),(50) so that (29) becomes 2J w + [(J ) + iJ +k(k2+2) /]w = 0 (51) with solution w = w(,7) [J(1,/J ]2 e-ik(kn (52) Evaluating the Jacobian and using (31.b), we obtain w = W - - (k2+27)A(,7) [cos s + + k2+T2 s sin s]-/ e- (53) and (o) (o) - Np+ v = w e The solution for v () is obtained by changing (k2+72) 2 to -(k2+72) in the above formulae. The parametric solutions given above are complicated and merit some discussion. Before undertaking this, we give the solutions for u and C. For the Rossby mode, let the column vector < uR, vR, CR > be written in the form iNcpR < uR v R > R e (54) where 4R is another column vector. Substituting (54) into (8), (9), and (10), and recalling the definition of T, we find that

15 0 iy cp -iy 0 P ~)= y Ko) = (55) cp (P 0 Dx *y where cp is pR. The matrix is rank 2 and one can solve for the first and third components of () in terms of the second. This yields (o) Cy (o) U = _ v(~ R px R (o) = _ (~) (56) R Cp R Similarly, for the gravity waves, we let iNcp+ < U+ V+, + > = t+ e - (57) and obtain iy x <?t C x -iy P cPy (0) o 0 (58) X t y y i+ cPx My t where p is c+. By virtue of the dispersion relation (21), the matrix in (58) is singular and of rank 2, and u( and (o) are found to be given by (o) _x t - yy v(o) u y v+ i (59)'y'Pt Ycx

16 y2 2 (o)' (o) = v. (59) + py pt + iycpx We now turn our attention to the behavior of the solution in the initial stage of the motion. Now for small T, the transformation (36.b) becomes k2 -_2 x= + k2+L2 y = -, (60) T = (k2+pt2)v with an error of order T2, and the Jacobian becomes k2+,2, with a similar error. Thus for small T (o0) iky y2-_k2 ) R k2y2 A(x + (Yk2)2 Nlt, y) R k ~~~~+y'( ~~~y ~(61) 3 exp(iN[kx + k2+y2) hNt]) + 0(T2). Similarly, (47.b) becomes x = ~ + ks, y = 7, (62) t = (k2+n2) s, the Jacobian is (k2+72)/, and (o) 1 iky kt v+ - - k A(x y) + 2 k2+y2 k2+y2 J^~ ~~ k ~+y2~ ~(63) 1 k * exp(i[N(kx ~+ k2+y2 t) - 2 k2+y t]} + 0(t2). Equations (61) and (63) state that the energy of each mode travels with the group velocity of that mode. Since

17 the group velocity of the gravity modes is much greater than that of the Rossby mode, the former disperse more rapidly from local sources. This much could have been anticipated, especially in view of recent work on the concept of group velocity (c.f. Landau and Lifshitz, 1959, section 66). What is surprising is the formation of envelopes of the rays on which the amplitude of the dependent variables is infinite, at least according to the ray approximation. Consider, for example, the solution w = w( )[J () /J, ] e- ik(k2+n2) (52 bis) CD = k. -/s - 1/42 sin 2s, (47.a bis) where x = ~ + ks, y = 7 cos s, (47.b bis) t = (k2+772) s for the + gravity mode. The Jacobian, given by J = [(k2+72) cos s + 72 s sin s] ( +................... (64) V k2+72 vanishes for those values of 7 and s such that n 2 k2 cos s ( s sin s + cos s and consequently for the curves described parametrically by y2 = _k2 cos s t2 = k s sin s (66) cos s + s sin s s sin s + cos s On these curves, shown in Fig. 3, the ray solution for

18 v( is infinite. Also, the curves are envelopes of the rays, or caustics, and we must anticipate that at least in their neighborhood the transformation is not one-one. Actually, the situation is far worse. Consider first the cusps of the caustics, y = 0, t = k(j - Y2)wT, where j is any positive integer. The preimages are found by solving 0 = 77 cos s, k(j - 2)7 = (k2+72)12 (67) One solution is 7 = 0, s = (j - 2)T. The others are s= (m - ), 2 = k2 (j /2 T (mn /2)2 (68) (m4/2)2 where m is any integer which is less than j. Consequently, the first cusp has one pre-image, the second three preimages, the third five pre-images, and so forth. For considering ordinary points on the caustics and points off the caustics it is useful to plot surfaces of constant y in the 77-s plane. Using (47.b) and the derived relation t) - = (y t) - (y. ) / (ys ) = J Os St (69) sy (y,s) a(77.s) (77s JSCOS / s we find that every ordinary point on the jth caustic has 2j pre-images, while every point between the jth and (j+l)th caustic has (2j+l) pre-images. The only region for which the transformation is one-one is for those values of t between zero and the first caustic. Elsewhere the solutions are either multi-valued or infinite.

19 A similar analysis can be made for the other gravity mode and for the Rossby mode, and similar results are obtained. Consequently, further discussion is necessary. This is supplied in part in the next section and in Appendix A. 4. Asymptotic solutions near the caustics The remarks in Appendix A indicate that in a region in which solutions are multi-valued the correct solution is obtained by adding the different values. Thus, for example, at any point between the jth and (j+l)th caustic for the gravity mode the ray solution has (2j+l) values corresponding to 2j+l pre-image points. The correct solution is obtained by summing over all the pre-image points. In order to carry out this procedure one must assign the correct branch to the square root of J/J. This is known as finding a phase shift rule. Both the phase shift rule and the correct asymptotic expansion near caustics are obtained by making a. local boundary layer theory. This theory is supplied here, in part, for the gravity modes Consider the equation t = s(k2 + cos2 s) F(sy) (72)

20 which is obtained by eliminating 77 between the last two equations of (47.b). Let overbars denote values on a caustic, and let t = t + a, y = y + b, s = s + a, 7 = 77 + p, (73) where a, b, a, and p are small. We note from (69) that F = 0. Consequently, F is the slope of the caustic in a (yt) plane, and z = a - F b measures distance from the caustic. Now from (72) and (73), a = Fb+(F b2 + 2F b + F b2) +... (74) y ss sy yy and with the aid of Newton s diagram we find that (74) considered as an eauation for a, has two solutions which tend to zero as a and b - 0. The lowest order approximation for these is = - I+ -- ss and we similarly find that the corresponding values of p are given by p = (I tan s) C. (75.b) Consider the first caustic. Since F = (J) /Cos S ss + s is positive and z is positive for t > t and negative for t < t, we see that there are real solutions for a and p for values of y and t on only one side of the caustic. In terms of rays, shown in Fig. 5, this means that through

21 a point r near the caustic in region A there pass two rays which touch the caustic and which coalesce if the point is on the caustic. There is also a'non-singular' ray through r which does not touch the caustic in the neighborhood of r and which need not concern us. Our starting point in the boundary layer analysis is equation (18) which, when the dispersion relation is taken into account, becomes ( + iNcp) (iN[2(Vp-Vw - twt) + (V2 - tt )w]+ + (V2w -w t) + Nw + iN2 w = 0. (76) This equation was solved in section 3 by expanding w in a perturbation series in powers of NS1 under the assumption that all derivatives are 0 (1). This assumption is not valid near the caustic; instead we assume that the derivative in the direction normal to the caustic is much larger than the tangential derivatives. The ensuing analysis is straightforward but very tedious, and we give only the results. Let 2/3 (a b) N 2/3 Z = N (a - F b (77) Q = + (2/F (78) - S and 6 = 2Q(k2+n,2)/ It is found after a standard boundary layer scaling and much algebra that the boundary layer version of (76) is

22 wzz + i6(zz + -) = 0, (80) and that near the caustic Ncp = P + 1 6 z3/2 (81) where P is a linear function of x, y, and t. The solutions for w and p corresponding to the positive value of Q are associated with a ray which has already touched the caustic, while the negative value of Q denotes a ray which has not yet touched the caustic. Letting subscript 1 label the former case and subscript 2 the latter, we find, upon solving (80), that -1/3 idZJ 2l) 2/3 w 2 = e [(A12 Ai[-( ) Z] + 16rsl 2/3 + B1,2 Bil-(- -) Z]) (82) where A and BI are independent of Z and where the 1, 2 1,2 standard notation for the Airy functions is used. Now from (64) and (75), we find that near the caustic (J+) 12 is given by (J +)1 ~ -z (J+)2 +z, as z - 0, with the same proportionality constant. Therefore w2 as given by the ray approximation is given by w = S z2 as z - 0, where S is independent of z, and we pick A2 and B2 so that w2 as given by (82) matches this as Z -+ + a.

23 Also, we require that iNcp1 iNcp2 v = e w +e 2 0 as Z -> -o, since there are no singular rays through points on that side of the caustic with z < 0. This determines A and B and hence the solution, which is _1/3i6Z3/2 wl = - C e13i (Bi(q) + i Ai(q)}, (83.a) 1/3i6Z 3/2 w2 = C e 1/3Z Bi(q) + i Ai(q) ), (83.b) where q is the argument of the Airy functions in (82) and where C = S 2 ( )1/6 e /4 (84) This calculation serves two purposes. First, it shows that near the caustic the amplitudes are large, of order N /6 Secondly, since in the limit Z -> o (83.a) becomes =e —i7/2S -l / it provides a phase shift rule, namely that the correct branch 1/f Y2 is e-iv/2 branch of J2 when J < 0 is e/ |J A similar calculation for the - gravity mode gives instead a phase shift of + T/2, but otherwise the results are identical. This analysis is incomplete in that it applies only to ordinary points on the caustic. It is conjectured that

24 near the cusps, which lie on the equator, the amplitude is even larger. However, no treatment of this case has been made. Also, the corresponding analysis for the Rossby mode has not been developed, although at least for the ordinary points of the caustics the treatment is probably not difficult. 5. Discussion. The work presented here emphasizes how important it is to take account of the earth's curvature, since it is the curvature which causes the waves to refract and to form envelopes. If the gravity waves are considered to be ~meteorological noise,," it must be said that at each latitude there are periodically recurring bursts of noise which have greater amplitude than the "signal," the Rossby waves. It is no longer obvious, at least to the author, that geostrophic adjustment through horizontal dispersion is a feasible explanation for the dominance of Rossby waves in the large scale motion of the atmosphere. Another effect of the curvature is to produce the phenomenon of equatorial trapping, in which the rays are always confined in a latitude belt around the equator. Bretherton (1964) and Longuet-Higgens (1965) have discussed this effect at length for time periodic waves.

25 An important assumption of the present work is that N = 2QR</(gH) the ratio of the radius of the earth to the radius of deformation, is large. The commonly accepted scale height of 8 km. for the atmosphere makes N = 3.3, which is not particularly large. However, this scale height applied to only one of an infinite number of vertical modes, and for the other modes the scale height is smaller and N is larger. The present work is therefore probably qualitatively inaccurate for that mode corresponding to a scale height of 8 km., but accurate for the smaller scale modeso

APPENDIX A Despite important work by Lax (1957) and Lewis (1964), the ray method has not been used extensively for the solution of initial value problems. In this Appendix we discuss another motivation for the method and indicate how to treat difficulties arising from the existence of caustics. We begin by considering a function O~ t) = ( )liN[ x-u( )t+f(~)]d~ (A 1) e(x,t) = (2-) / afE('e (A.1) -00 which may be thought of as the solution of an initial value problem. When N is large, 0 = q|'F(k) I / a(k) exp{i[N(kx-o(k) t + f(k)] + - 7/4 sgn F(k)) + 0(-1) (A.2) by the method of stationary phase, where F(k) = tWg" (k) - f" (k) (A.3) and where the sum is over all values of k such that x - c& (k)t + f' (k) = 0. (A.4) Equation (A.4) defines a one-parameter family of straight lines in the (x,t) plane, the slope of each line being the group velocity cu (k). We define these straight lines as rays and note that 26

27 for any x and t the sum (A.2) is over those rays which pass through (x,t). Alternately, we could introduce a new parameter s such that x = C' (k)s - f' (k) t = s (A.5) so that (A,5) defines a mapping from the (k,s) plane to the (x,t) plane. Then for any point (x,t) the sum (A.2) is over all the pre-image points. On the locus of points (x,t) such that 3x = F = 0 (A.6) (k,s) the method of stationary phase is invalid. This locus is an envelope or caustic of the family of rays. Near ordinary points of caustic (A.2) is replaced by the wellknown uniformly valid expression of Chester, Friedmann, and Ursell (1957). Presumably it is possible also to find a uniformly valid expression which is applicable to cusps of caustics, but this has not yet been done. These observations are useful for interpreting the ray method of solution. The expression (A.2) implies that, to find an asymptotic solution of the problem which has (A.1) as its exact solution, one should assume 0 to be of the form e(x.t) = A(xtN)e (t) (A.7) and then solve by expanding A in powers of N1. When the ray method is applied to problems which can be solved exactly by Fourier analysis the solution should

28 be and is identical to that obtained by approximating the exact solution through an asymptotic integration. Difficulties arise when the rays have envelopes, for on such envelopes the ray solution is singular and a local boundary layer theory is needed. Calculations by Buchal and Keller (1960) and by the author for some sample problems indicate that the boundary layer theory yields the same solution as is given by the method of Chester, Friedmann, and Ursell. A second difficulty associated with the envelopes is that more than one ray can pass through a point, and at this point the solution is multi-valued. This is resolved by summing the different values in accord with the remarks following (A.4) and using the boundary layer theory to provide the phase shift rule. A third difficulty is that the initial data may not be of the form (A.7). Lewis (bp. cit.) indicates how this case may be treated and provides a sample calculation. A great advantage of the method is that it can be used to solve equations with non-constant coefficients provided that the scale over which the coefficients vary significantly is much greater than the scale of the dependent variables. Thus many apparently intractible problems can be treated, An alternate procedure is to use the closely related method of multiple scales (c.f. Carrier, 1966). However, the ray method is more systematic and lends itself

29 more readily to obtaining higher approximations. APPENDIX B. The calculations here have been carried out under the 3-plane approximation, and it is natural to inquire what changes in the previous results occur when this approximation is not made. The dispersion relations, at least, are easily obtained, and are similar to those obtained under the P-plane approximation. The method used is due to Lax (op.cit.). Let u and v be the non-dimensional velocity components, and let u = mU, v = mV, (Bol) and recall that f = sin c = tanh y, m = sec p = cosh y. (5 bis.) Then the non-dimensional eauations of motion in Mercator coordinates are at -NfV+ = 0. (B.2) -+ NfU+ at= 0+ (B.3) -+ NfU + -0, +t m2 (a + =0 (B. 4) s b r we a shat soutios ay be it a s before, we assume tha t ions may be spli t into a Rossby mode with time variable z = N't and two gravity

30 modes. For the Rossby mode, let iNpR < UR VR' R > R= e(B5) where as before the vectors are column vectors, and substitute into (B.2) through (B.4). This yields i A R + N D R + N ('R) = 0 (B.6) where 0 if c x A = -if O c m m2cp 0 X y m x Y! -_ (B.7) cp 0 = and p i. T y (0> 0_o 0 m a m2 a and p is R.' Expanding ~R in powers of Nr1, we obtain AI O, (A ) D 0, 3.8) A?R i A D R 0 (B.8 and higher order equations. The matrix A is of rank 2 and has right and left nullvectors such that eA = Ar = 0, (B.9)

31 where e is a row vector and r a column vector. These are uniquely determined except for a multiplicative constant, and a convenient choice is = (m2p, -m2cX, -if), r = < p, -cp if >. (B.10) y ^. x y From the first equation in (B.8) () is a scalar multiple of r, (0) = () r; (B. R R and if we substitute this into the second equation in (B.8) and multiply by Q, we obtain (o) Q D e( r = 0, (B.12) R which after the necessary algebra is fcpT[m2(q2 + cp) + f2] - xp} (o =( (B 13) Hence the dispersion relation obeyed by cR is T =[m2(cp2 + c2+ f2] (B.14) x y The dispersion relation for the gravity modes is obtained much more simply. Let iNcp+ < U+, V+ = + e (B.15) substitute into the equations of motion and expand A+ in powers of N1. The lowest order equation is

32 (Pt if (Px -if (0) -if cp t |Y) = _ (B.16) m2Tx m2Py cp where cp is cp+, and for (o) not to be zero the matrix must be singular. Its determinant is p[cp2 - f2 m2(p2 + p2)] tx y and the case Tt = 0 has already been incorporated in the Rossby mode, Consequently, c obeys - = ~[m2(p + + ) + f2]1 The ray paths for the different modes can be obtained by quadrature, and are qualitatively similar to those obtained through use the 3-plane approximation. Acknowledgments This work was supported by the National Science Foundation under Contract GP-2561.

REFERENCES Buchal, R.N. and Keller, J.B. 1960 Boundary layer problems in diffraction theory. Comm. Pure Appl. Math. 13, 85. Bretherton, F.P. 1964 Low frequency oscillations trapped near the equator. Tellus, 16, 81. Carrier, G.F. 1966 Gravity waves on water of variable depth. J. Fluid Mech. 24, 641. Chester, C., Friedman, B. and Ursell, F. 1957 An extension of the method of steepest descents. Proc. Camb. Phil. Soc. 58, 599. Courant, R. and Hilbert, D. 1962 Methods of Mathematical Physics, II, Interscience. Keller, J.B. 1958 A geometrical theory of diffraction. Calculus of variations and its applications, Amer. Math. Soc., L.M. Graves, Ed., McGraw-Hill. Lamb, 1932 Hydrodynamics, Dover. Landau, L.D. and Lifshitz, E.M. 1952 Fluid Mechanics, Pergan oni Press. Lax, P.D. 1957 Asymptotic solutions of oscillatory initial value problems. Duke Math. J. 24, 627. Lewis, R.M. 1964 Asymptotic methods for the solution of dispersive hyperbolic equations. Asymptotic solutions of differential equations and their applications, C.H. Wilcox, Edo, John Wiley and Sons. Longuet-Higgens, M.S. 1965 Planetary waves on a rotating sphere. II. Proc. Roy. Soc. A, 284, 40. Mahoney, J.J. 1962 An expansion method for singular perturbation problems. J. Aust. Math. Soc. 2, 440. Obukhov, A. 1949 K voprosu o geostroficheskom vetra. Isv. Akado Nauk SSSR, Ser. Geograf. Geofiz., 13, no. 4. 33

34.rd — LO 11 0i II.OH 0 / i H - ~ 0 10 11 ci 0 Hon o1 D / lo 0 1 \rc HC H O b F IIr o 11 0 ~ Inb.O In \^^"~~~~~~ r-t.r-i

35 -d11 H II Cd rI OD I" b0 -, U) // // a O ^, I" (1C\i 11 /0 I -)0 C^/ ~om

36 II / k o 0 bD ~\/ ~ ~a+ o 4-; CO IfI I~~~~~~~~~~~~~~~~. N, N

37 crr Ue 0 11 Ho II Q) C) m I 0 bQ. ^ ^I~~~~~ p__ I ____ ^\~~~~~~~~~~~~~~~~~' lnVQ 0 ^ T 33 / ^a Ic) ~ 0 Ic8)

58 C / Y A /T - /. /2 NS I Figure 5. Ray paths near a caustic. (1) and (2) are rays which touch the caustic (c), and NS is a non-singular ray.

UNIVERSITY OF MICHIGAN 3 90111115 03025 4497II 3 9015 03025 4497