THE UNIVERS I TY OF MI C H I GAN COLLEGE OF ENGINEERING Department of Meteorology and Oceanography Technical Report ON BAROCLINIC AND BAROTROPIC AGEOSTROPHIC STABILITY Joseph Sela Stanley J. Jacobs ORA Project 07344 under contract with: OFFICE OF NAVAL RESEARCH DEPARTMENT OF THE NAVY CONTRACT NO. Nonr-1224(55) WASHINGTON, D.C. administered through: OFFICE OF RESEARCH ADMINISTRATION ANN ARBOR September 1968

ACKNOWLEDGMENTS The author wishes to express gratitude to the Graduate School for the Rackham scholarship and to the Office of Naval Research for their support. Sincere appreciation is expressed to Dr. Stanley Jacobs for his valuable discussions and suggestions and to Dr. A. Wiin-Nielson for availing his experience in interpreting the results of this study. Thanks are due to Mr. F. Brock for his guidance in the analog computer phase and to Mrs. K. Tingle for her assistance in using the department's computer. Thanks are also due to Mr. E. Downing of I.B.M. for providing encouragement during the time sharing experiment. For their exercise in hieroglyphics thanks are extended to Mrs. B. Paddock and to the efficient staff of the University's Reports Office. ii

TABLE OF CONTENTS Page LIST OF TABLES v LIST OF FIGURES vi ABSTRACT vii I, INTRODUCTION 1 1.1 Background and Purpose of the Baroclinic Study 1 1.2 Background and Purpose of the Barotropic Study6 1.3 The Stability of the Gulf Stream 8 II. THE BAROCLINIC MODEL 12 2.1 Formulation 12 2.2 Analytical Considerations 19 2.3 The Singular Regions in the c Plane as Neutral Boundaries 21 2.4 The Nonsingular Regions of the c Axis as Neutral Boundaries 26 2.5 A Necessary Condition on the Eigenvalues 28 2.6 The Numerical Study of the Eigenvalue Problem 30 2.7 Analog Computer Solution of the Eigenvalue Problem 32 2.7.1 Method of solution 33 2,7.2 The digital solution 41 208 Results of the Baroclinic Study 42 2.8.1 Growthrates and propagation speeds 42 III. THE BAROTROPIC GULF STREAM MODEL 62 3.1 Formulation 62 5.2 The Numerical Treatment of the Gulf Model 67 5.4 Results of the Barotropic Study 73 3.4.1 Quasigeostrophic unstable modes 73 354.2 Ageostrophic unstable modes 78 IV. CONCLUSION 84 4.1 Summary of Results 84 4.1.1 The baroclinic problem 84 4.1.2 The barotropic problem 86 4.2 Concluding Remarks and Suggestions for Future Work 88 4.2.1 The baroclinic problem 88 4.2.2 The barotropic problem iii

TABLE OF CONTENTS (Concluded) Page APPENDIX 90 A.l Transformation of Equation (2.9) Into Canonical Form 90 A.2 Determination of the Exponents From the Indicial Equation 91 A.3 Quasigeostrophic Solution of the Baroclinic Problem 92 A.4 Calculation of Nongeostrophic Pressure and Temperature 94 A.5 Quasigeostrophic Pressure and Temperature 95 BIBLIOGRAPHY 97 iv

LIST OF TABLES Table Page I. Table of Eigenvalues and Growth Rates for I = 0.5 44 II, Table of Eigenvalues and Growth Rates for I = 1.0 45 III. Table of Eigenvalues and Growth Rates I = t 46 IV. Unstable Modes for the Gulf Stream Model with Countercurrent (A = - 7/6) 7 = 3, y = 2.5 74 V. Unstable Modes for the Gulf Stream Model with No Countercurrent (A = - 5/6) 7 = 2.5, 7 = 2 75 VI. Unstable Modes for the Gulf Stream Model with Countercurrent (A = - 7/6) y = 2, 7 = 0.5 76 v

LIST OF FIGURES Figure Page 1. Separation of the c-plane by the integral condition. 22 2. Analog computer schematic for eigenvalue search. 36 3. Percent deviation of ageostrophic growthrates for ~ = 1, B = 0.5. 48 4. Percent deviation of ageostrophic growthrates for I = 1, B = 1.0. 49 5. Percent deviation of ageostrophic growthrates for ~ = 1, B = 1.5. 50 6. Percent deviation of ageostrophic growthrates for I = i, B = 0.5. 51 7. Phase lag between temperature and pressure fields. 55 8. Vertical cross sections of vertical velocity perturbations. 57 9. Vertical cross sections of vertical temperature perturbations. 58 10. Vertical cross sections of vertical pressure perturbations. 59 11. Tilt of vertical velocity cells as a function of k, ~a R, and B. 60 12. Basic state velocity profiles of the barotropic model. 68 13. Unstable quasigeostrophic modes of the barotropic Gulf Stream model. 77 14. Unstable nongeostrophic modes of the barotropic Gulf Stream model 7 = 0.5, 7 = 3. 79 15. Unstable nongeostrophic modes of the barotropic Gulf Stream model = 2. 80 16. Unstable nongeostrophic modes of the barotropic Gulf Stream model 7 = 2.5. 81 vi

ABSTRACT Ageostrophic influences on the hydrodynamic stability of two geophysical phenomena are investigated. The first part of the study treats the baroclinic stability properties of an atmospheric model. Linear theory is used to study an adiabatic frictionless atmosphere regarding the Coriolis parameter as constant with the additional simplifications of the hydrostatic and Boussinesq approximations. The resulting eigenvalue problem formulated in terms of the vertical velocity is studied both analytically and numerically. The object of the study is to examine the influence of the transversal coordinate. The motivation for this examination is found in the fact that the north-south variation and the first order Rossby terms are very strongly coupled in the w equation. In quasigeostrophic theory a zero Rossby number removes the influence of lateral variations leaving only a quantitative distortion in growth rates through a modified wave number. It is found in this study that when the ageostrophic terms are not neglected the lateral wave number has a selective influence. For long lateral waves the short zonal waves are stabilized whereas the long ones are more unstable than their quasigeostrophic counterpart. For short lateral waves all zonal modes are found to be destabilized. The method of solution consists of an analytic study removing large portions of the c plane from the domain of possible eigenvalues, while the real c axis is studied as a neutral stability boundary. It is found that if the product of the Rossby and zonal, wave numbers is less than unity the neutral stability boundary is confined to the range of variation in the basic zonal speed. The numerical, solution employs a shooting method for searching the c plane for eigenvalues. The search is carried out for various combinations of wave numbers, static stability, and Rossby number. For long transversal waves the ageostrophic effects are found to be greatest for a preferred static stability. It is also found that an increase in the Rossby number is such as to always increase its effect. Vertical distributions of vertical velocities, pressure, and temperature perturbations are computed and plotted. The second part of the study treats the barotropic stability problem of the Gulf Stream. A two-layer model is considered with no motion in the lower deep layer and a basic state zonal jet in the upper one. The purpose of the study is to examine the validity of hydrodynamic stability investigations of the Gulf Stream in view of the overestimated growth rates resulting from quasigeostrophic theory. It is found that the ageostrophic influence depends on the basic state velocity profile and on the radius of deformation y7 In the case vii

of basic velocity profile without a countercurrent an increase in the Rossby number results in a more stable flow. When a countercurrent is included the results depend on the value of the radius of deformation and the wave length. It is found that for a realistic value of y the short zonal waves are destabilized while the long ones are stabilized. It is also found that an increase in y stabilizes the flow for both profiles studied, in agreement with analytical results obtained by Jacobs. viii

CHAPTER I INTRODUCTION Two geophysical phenomena to which stability theory can be applied are the large-scale atmospheric traveling waves associated with the westerly current in midlatitudes and the Gulf Stream meanders near the North American continent. The present study will concern itself with an investigation of nongeostrophic effects on the stability properties of the two flows. Chapter II deals with an atmospheric baroclinic model which except for its derivation is identical to Eady's 1949 formulation. Chapter III studies a twolevel barotropic model of the Gulf Stream. Since the studies of the two types of flow have developed rather independently we shall proceed to describe the problems separately. 1.1 BACKGROUND AND PURPOSE OF THE BAROCLINIC STUDY The first attempt to establish a relationship between the equations of hydrodynamics and the observed upper air waves in extra-tropical and midlatitude regions was made by Rossby (1939). The importance of the local variation of the vertical component of the earth's rotation was thus first recognized as the main mechanism determining the speed of such perturbations coupled to the basic zonal westerly flow. The two main shortcomings of that model were the retrogression of the very long waves and their neutral character. Since the synoptic map indicates a continuous process of birth and death of individual disturbances, models that can support energetically active perturbations had to 1

2 be constructed. The energy conversion properties of atmospheric models came into even sharper focus in the light of observational evidence as demonstrated by studies of momentum transport in relation to the maintenance of the general circulation Starr (1954). It is now generally accepted that large-scale disturbances are an energy conversion agent participating in transporting momentum and heat. Regarding these disturbances as mature stages of instability, built in energy sources become essential. The source of energy responsible for the motion of the atmosphere is the external solar radiation. Since we shall be dealing with the initial stages of developing perturbations and shall limit our discussion to time intervals that are short compared to "time constants" of differential radiative processes it is reasonable to treat the flow as adiabatic. Once the hydrodynamic system under consideration is assumed isolated, energy sources from within must be included. The forms of energy on which an atmospheric disturbance can grow are either the kinetic energy of an existing flow or the available potential energy residing in a given temperature density distribution. In a hydrostatic atmosphere the thermal energy is proportional to the potential energy and the two can be combined. In large-scale motion of the type envisaged the hydrostatic approximation can and will be adopted. It is pointed out that the existence of these forms of kinetic or available potential energy and the state of the atmosphere required to bring about their existence is assumed; they are specified according to the type of effects one wishes to study. In specifying the basic state of the flow a zonal steady current is prescribed. If the current is assumed constant with height and latitude neutral

35 Rossby waves are obtained. The variation of the zonal current with height, which must be associated with a horizontal temperature gradient (through the thermal wind relation) incorporates a source of available potential energy derived from the latitudinal differential heating. On the other hand, a'horizontal shear in the basic state wind can convert the kinetic energy of the zonal flow through horizontal Reynolds stresses. The large scale vertical stresses are an order of magnitude smaller since they involve the vertical velocity fluctuations and are therefore too weak to account for observed instabilities. In forinulating a detailed stability study it has been customary to specify either a horizontally sheared zonal current, the barotropic case, or a vertically varying flow in a baroclinic basic state. When both variations are included the results within a quasigeostrophic formulation appear to be more general; they are conditions under which instability may develop, and bounds (semicircle theorem) on the phase speeds, e.g., Pedlosky (1964).'The first study of a baroclinic model was performed by Charney (1947). His model includes a beta effect but the perturbation quantities are assumed independent of latitude. Stability criteria are derived, and it is shown that instability increases with vertical wind shear, lapse rate and latitude, and decreases with wavelength. A comparison with average seasonal values of the relevent parameters suggests that the Westerlies of midlatitudes are a source of constant instability. The unstable solutions of Charney's model possess similarities to observed atmospheric disturbances; e.g., the motion is generally to the east, the speed is approximately equal to the surface zonal current, and

4 the vertical structure has a westward tilt with height. An attempt to include latitudinal variation in the perturbation variables was made by Eady (1949) at the expense of dropping the variation of the Coriolis parameter with latitude. In addition, the terms involving the Rossby number (and its square) were neglected. Thus Eady's model does not include any ageostrophic effects and should be considered as the zeroth order approximation in a framework of series expansions in the Rossby number. Once again broad agreement with the observed atmospheric flow was found. The interval 1949-1964 saw the emergence of a large number of studies of various versions of the baroclinic stability problem. Common to all investigations was the simplification of the eigenvalue problem governing the vertical velocity. These simplifications are usually traced to the argument of small Rossby number (and of course its square). In a short critical paper by Phillips (1964), titled "An Overlooked Aspect of the Baroclinic Stability Problem," the question of approximating the exact characteristic equations arising from model atmospheres is treated. The choice of example fell on Eady's model, probably because it was the first to include north-south variation in the perturbation quantities. It is there shown that the effect of neglecting the first order Rossby terms (but only the first order) are interchangeable with dropping the north-south dependence in the linearized w equation, insofar as the simplified resulting equations are identical. The 2 remaining term in R may still manifest itself in some nongeostrophic modification. Thus the inclusion of latitudinal variation in the linearized quantities together with the quasigeostrophic assumption has the effect of modifying the

5 computed growth rates but quantitatively only (through the modified wave num2 2 1/2 ber (k + i ) / ) and in the direction of stabilizing the flow when the northsouth wave number I is increased. On the other hand if only the second order Rossby terms are neglected the y variation becomes much more significant and may affect the nature of the influence. The paper reaches the uncomfortable conclusion that since the first order Rossby terms are so strongly coupled to the y variation of the perturbation vertical velocity, greater care must be taken in solving nongeostrophic stability problems. The difficulties arising in such a treatment are related to the fact that under these more general conditions the characteristic equation is no longer separable in y and z in such a way as to satisfy vanishing transversal velocities at one or two straight vertical walls. The equations are, however, separable if such boundary conditions on v are not imposed. A variant on "nonvertical walls" was given by McIntyre (1965). In his treatment vanishing boundary conditions can be satisfied by the north-south velocity on two slanting walls, the slope being a function of the geometry of the channel and the Rossby number. The equation then becomes separable in y and z to the second order in R. The eigenvalue problem is formulated but not o solved and the agestrophic effects on speeds and growth rates of unstable disturbances are not apparent. At this stage of the discussion two possible methods of approach are: (1) Treating the partial differential equation describing the vertical motion by the method of expansions in terms of a small parameter

6 (R ) and including boundary conditions on the transversal velocity. (2) Ignoring the vanishing requirement of v on walls of a channel and separating the equation, thus reducing it to an ordinary differential equation. The ageostrophic influences will be retained and the solutions will indicate whether the y coordinate can introduce more than quantitative distortions. The first method was adopted by Dolph* and Derome at The University of Michigan, whereas the separation technique was applied in this study. Before proceeding to the formal presentation it should be mentioned that the further approximations incorporated are the neglect of friction, which is justified away from the boundary layer and for short-time intervals, and the IBoussinesq approximation which is adequate under low Mach number conditions as are found in the free atmosphere away from weather phenomena and when the vertical extent of the fluid is not too large, Phillips (1967). 1.2. BACKGROUND AND PURPOSE OF THE BAROTROPIC STUDY The stability study of atmospheric flow in a barotropic environment where the basic zonal wind varies with latitude was first undertaken by Kuo (1949). As in Charney's model the introduction of perturbations whose amplitude is a function of one coordinate only (in this case the north-south coordinate) allowed the retention of the beta term. Implicit in the derivation are simplifications due to the assumption of small Rossby number through the two-dimensional nondivergent character of the basic flow. The main result of the study is a der*Three-Dimensional Disturbances in the Non-Geostrophic Eady Model of the Atmosphere. Currently submitted for publication in J. Atmos. Science.

7 ivation of a necessary condition for instability to occur. It is shown that the condition for the presence of neutral and amplifying waves, with a phase velocity whose value is between the maximum and minimum wind velocity in the belt of a strong westerly jet, is the existence of critical points where the absolute vorticity has an extremum. The wave moving with the velocity of the current at a critical point is neutral, while waves with speeds intermediate to the minimum speed and critical point speed are unstable. The study is concluded by a successful comparison with a synoptic example. The problem arising in allowing a two-dimensional variation in the amplitudes of the perturbation quantities has been bypassed in barotropic stability problems through the introduction of multilayer models. An exception is found in the work of Jacobs and Wiin-Nielsen (1966) where divergence effects in a quasigeostrophic atmospheric model are found to give rise to several internal unstable modes. Wiin-Nielsen (1961) applied a vorticity equation to the midtroposphere. The only nondivergent term was the term involving the advection of vorticity. It was found that the inclusion of divergence effects reduces the rate of growth of unstable disturbances and confines the instability to a narrower band of wavelength. Lipps (1962) considered a two-layer fluid in which the basic velocity profile was given by a (hyperbolic) westerly jet. The model was applied to both atmospheric and oceanic stability. The main features of the model were a beta plane approximation on an infinitely wide plane as opposed to a finite channel in the treatment by Wiin-Nielsen, and no motion was allowed in one layer, the

8 lower layer in the oceanic problem. The results of the study when applied to Gulf Stream instability give rise to growth rates that are too large when compared with observations, See Stommel 1965, p. 196; also, Gulf Stream Summaries). The object of the present study is to determine the nongeostrophic influences on the growth rates and speeds of small perturbations superimposed on a zonal jet resembling the Gulf Stream. A comparison of stability is also made between the cases of a simple jet and one including a courter-current. The study is conducted with the radius of deformation y and the Rossby number R as parameters. 1.3 THE STABILITY OF THE GULF STREAM The Gulf Stream meanders can be explained by the effects of topography when the zonal current is lowered in leaving the coast while conserving the potential vorticity of its columns. Their time variation is attributed to hydrodynamic instability. When instability techniques are applied to Gulf Stream theory it is convenient to consider a two-layer model thus circumventing from the outset the problem of separating variables in the resulting partial differential equation determining the transversal velocity v. The resulting problem can be characterized as an ordinary two point boundary value problem for the complex amplitude function of the transversal velocity perturbations in the upper layer, subject to vanishing boundary conditions at two latitudes. The motion in the lower deep layer is neglected and as a result no potential energy converting processes can take place.

9 The two-layer formulation is equivalent to single layer flow under the influence of a modified (reduced) gravity if the motion is assumed hydrostatic and the lower layer is much thicker than the upper one. The hydrostatic approximation is valid since the vertical scales are small compared to the horizontal scale of the model. In addition the assumption R < 1 and F < 1 must be valid for the one layer equivalence to hold. R is the Rossby number and F is the o Froude number. The situation F> 1 does not affect however the applicability of the two-layer model for comparison with observations, Stern (1961). As in atmospheric stability problems, the ordinary differential equation governing the unknown on which boundary conditions can be imposed from physical considerations, has solutions which cannot be easily found. Once again the starting point is the quasigeostrophic approximation with the sharp contrast that the only analytical solutions that have been found in the barotropic case, corresponding to the basic state profiles that have been studied, give rise to neutral modes only. In order to determine possible unstable modes the "perturbation off the neutral boundary" technique can be employed. This method is applicable since it can be shown that certain neutral modes that were found by analytical methods are indeed on a stability boundary, in the sense that they are limiting cases of unstable modes. The main limitation of the perturbation method is the necessity to stay close to the neutral boundary without the ability to cover longer ranges of variation in the physical parameters, wave number, Rossby number, and Froude number~ The problem of varying two parameters simultaneously is clearly much more cumbersome. The jet studied by Stern (1961) gives rise to an analytical solution for a

10 zero wave number thus corresponding to very long waves only. In the present study the two profiles treated are trigonometric jets with and without a counter-current; in both cases the analytical quasigeostrophic neutral modes correspond to nonzero wave numbers. The first phase of the study is to determine the dependence of the complex wave speed on wave numbers considerably far from neutral stability. It was mentioned before that the decisive problem arising in Gulf theory is that of the ageostrophic effects on the stability of the flow. Since quasigeostrophic theory overestimates the rate of growth of unstable perturbations, compared with observational values based on Stommel (1965), the direction of nongeostrophic influences is decisive in judging the validity of quasigeostrophic formulations as a first approximation. The recent work by Jacobs (personal communication, and to which this section of the study is a natural continuation) shows that complex phase speeds which are obtained as the Rossby number approaches zero are located inside a 1 1 semicircle of radius 2(u -b) and center at C = -(u+b). b = u if u < o and 2 2 2 2 1 1 b = 0 otherwise; u1 and u2 are the minimum and maximum of the zonal speeds in the basic state profile. Since quasigeostrophic phase speeds can be shown to lie in this semicircle too, it may be concluded that quasigeostrophic theory is avalid first approximation. On the other hand a calculation by Jacobs using a variational method perturbing off the neutral boundary to investigate the sense of nongeostrophic effects near neutral stability indicates that for a trigonometric jet resembling the Gulf Stream and including a countercurrent, the Rossby number effect is destabilizing in the sense that a marginally stable

11 configuration for zero Rossby number is unstable for a finite nonzero value. This is a surprising and uncomfortable result inasmuch as it is necessary for nongeostrophic effects to stabilize the flow in order to conform with observed rates of growth. Further calculations based on different profiles of the basic flow have shown that the Rossby terms in some cases stabilize the motion, in others destabilize it. A general conclusion can therefore not be drawn. To investigate the problem of nongeostrophic effects away from the neutral boundary a numerical approach is suitable and can provide the answer to the question of the validity of the quasigeostrophic theory as a zero order approximation. It is interesting to note that while analytical methods are effective near the neutral boundary, and hence close to singularities of the differential equation, numerical methods fail there. On the other hand further from the singularities, analytical methods fall short whereas numerical schemes are reliable and successful. The second part of the study will concern the unstable nongeostrophic modes for finite Rossby number. The investigation will be carried out for a large spectrum of wave numbers and for various values of Froude number. The question of behavior very near the neutral boundary must be answered by the analytical treatment.

CHAPTER II THE BAROCLINIC MODEL 2.1 FORMULATION We consider a hydrostatic incompressible energetically isolated atmosphere in the usual meteorological cartesian frame of reference. Let (x,y,z) be the position of a point relative to some origin on the surface of the earth and let t be time. The coordinates x,y,z increase to the east, north, and vertically up, respectively. Denote the velocity components by (u,v,w) and let g, p, p, represent the acceleration of gravity, the density, and pressure. Finally denote by f the Coriolis parameter (considered constant 0 in this study). The system of equations describing the flow are the equations of conservation of momentum,mass and energy modified to incorporate the physical assumptions of the model. These are: au u 8u u 1 2-1 an+u- + V -+ w =- 1 - +fv (2.1) at ax Wy az p ax o 0 av aV v 8v av 1 5p Y+ u- + v 4-+w —==- - - fu (2.2) a-t ax + y az po y o ap - a + g ~ o = +g 0o (2.3) p dz ~u av aw + u 5+v -+w - (2.3) - +U +V +W =0 ( 12

13 where the further assumption of constant density p in the horizontal pressure force terms is introduced. The Boussinesq approximation is justified if the flow is of slow meteorological character and the depth of the atmosphere is not too large. It is convenient for both analytical and numerical purposes to scale the variables to nondimensional quantities. The scaling adopted here is designed to confine the variation of most variables to the interval zero to one, while forcing a vanishing vertical velocity for zero Rossby number. The pressure and density are represented by some realistic standard distributions depending on height only and hydrostatically related, plus deviations from the standard. The deviations are so scaled as to bring about maximum simplification in the scaled equations. Let L be a characteristic horizontal dimension and V a typical horizontal speed. We shall consider an atmosphere of finite vertical extent and denote its depth by H, the upper boundary being a rigid plane. Let primed quantities represent the nondimensional scaled variables. The scaling is defined by: (x,y) = L(x'y') z = H-z t = L tt V (u,v) = V(u,v') w = R - o L

14 V wiwere R = — is the Rossby number, 0 f L o VfoLPo P P s(zl) gH g = (z') + VfoLpo wvhere the subscript s refers to the standard atmospheric distribution under consideration. The system (2.1)...(2.5) transforms into: R(t + u + v + R w )u = - + v (2.2a) o at dx 5y o z y R ( u aX+ v + R w )v = - U (2.2a) p = 0 (2.3a).az bu av bw - + - + - - - 0 +x +Y o dz (2.4a) ( + u x+ v + Rw )p = Bw (2.5a) dt x v y o z where gH2 Ps 22 dz fi L p o o is a measure of the static stability and the primes have been dropped. The o: dimensional system (2.la)..,(2.5a) is now linearized about a basic geostrophic state with linear vertical shear given by u =z v w = o p = -yz, P = y,

15 and the perturbations are assumed of the form ik(x-ct) t = i*(y,z)ek( t where k is the nondimensional zonal wave number and c = c + ic is the complex r i phase speed. For k>0 and c. >0 the motion will be unstable through the growth factor multiplying all the dependent variables. The growth factor is eki, kci being the growthrate referred to later in this study. If ci i 0, the solutions of the linearized system are complex and either the real or imaginary part of any variable represent physical solutions. Substituting the prescribed perturbation forms into (2.la)...(2.5a) the linearized system becomes: R [ik(z-c)u+R w] = - ikp+v (2.lb) 0 0 R ik(z-c)v = -u (2.2b) o I + p= 0 (2.3b) az av 5w iku + - + R a =0 (2.4b) 6y o 6z ik(z-c)p + v - Bw = 0 (2.5b) where the asterisks have been dropped. The natural procedure at this point is to eliminate as many variables from the system (20lb)...(2.5b) as possible thus reducing the number of equations at the expense of raising their order. When this course is followed p, p and u can be eliminated resulting in the system:

16 v - (z-c)v - R (z-c)v - R (z-c) v = Bw R (z-c)w (2.6) z o y o yz o zz (z-c)(v -k v)= w +Rw - R (z-c)w (2.7) yy z o y (2.7) where subscripts denote partial differentiation. Equations (2.6) and (2.7) together with boundary conditions w(y,z=O) = w(y,z=l) = 0 v(y=O,z) = v(y=l,z) = 0 constitute an eigenvalue problem from which eigenvalues c characterizing the speed and stability of the flow are to be determined. The problem can be further simplified by eliminating v from (2.6) and (2.7) thus leading to a formulation in terms of a single dependent variable. In terms of w the eigenvalue problem becomes: 2 22 15 (z-c)[w - k Bw + Bw ] - w + 2R [(z-c)w -w k R (z-c) w = 0 (2.8) zz yy z o yz y o zz and the boundary conditions are as before. A possible method for solving equation (2.8) is suggested by the fact that the coefficients of the equation are independent of the lateral coordinate. A Fourier Synthesis method would reduce the equation to an ordinary one, significantly simplifying the problem. The obstacle in persuing this path is the difficulty in obtaining expressions for w and its derivatives at y=o. Considering the method of separation of variables it is observed that the equation is not separable as pointed out by Phillips (1964), if the vanishing boundary con

17 ditions on v are to be imposedo It was suggested in the introduction that the ageostrophic effects resulting from the terms involving the Rossby number and its powers can be studied independently of the lateral boundary conditions' effect on the possible modes of flow. That this is plausible can be inferred from the fact that if the basic state is replaced by a uniform current the only effect of the lateral boundary conditions is to introduce the additional class of neutral Kelvin waves to the possible modes of flow. It is not claimed that this is necessarily the case in baroclinic flow; however, it is felt that this approach can lead to profitable results shedding some light on the possible quantitative effects of the ageostrophic terms. To this end the separation of w according to: A iiy w= w(z)e will result in an ordinary differential equation when iequation when introuced ino euation (2.8). The result is: (z-c)[l-k R (z-c)2] + 2[iR (z-c)- - 2 B(k 2)(z-c) 0 (2.9) 0 ZZ 0 Z 0 with the boundary conditions: A A w(z=o) = w(z=l) = 0. It is noted that for consistency all the other variables must be separated in a similar manner. It is also observed that while the transversal velocity v does not vanish on straight vertical walls, it does vanish periodically in y ithrough the factor e through the factor e for,

18 A/ i y 1A A v = v(z)e = v(z J[cos(argv+y) + isin(argv+iy)] and both the real and imaginary parts oscillate in the lateral direction. The possibility of unrealistic growing transversal velocities as y increases therefore does not exist. In the present formulation I is a parameter and must be specified. Since I plays the role of a lateral wave number it is reasonable to assign values to it of the order of the zonal wave number k in the meteorological range of interest. It will be shown in the subsequent development that the value of i has more than a quantitative effect on the growthrates. In fact the matter of stabilization or destabilization of the motion with respect to its quasigeostrophic counterpart as a function of the Rossby number will prove to be highly dependent on the particular value of I at hand. The observation by Phillips (1964) that the dropping of the ageostrophic terms from the w equation removes with it the effect of the y variation is amply verified. Equation (2.9) subject to its boundary conditions will be investigated along two lines. 1. A study of the complex c plane to determine the possible location of eigenvalues. It will be shown that under realistic atmospheric conditions fast moving unstable perturbations, essentially for Ic > 1, are not present in the r model. This does not preclude the existence of fast moving neutral disturbances. It indicates, however, that most of the real c axis is not a stability boundary in the sense of being marginal to unstable modes. In addition a necessary condition on the eigenvalues will be derived. This condition affects a partition

19 ing of the c plane according to possible versus prohibited regions for eigenvalues to occur. In this respect the condition plays the role of the so called semicircle theorem developed in quasigeostrophic theory. The integral condition will prove most useful in the subsequent numerical study. 2. A numerical investigation leading to actual eigenvalues and growthrates, as well as vertical distributions of the physical variables, as a function of various combinations of static stability, Rossby number and zonal and lateral wave numbers. The results of this part of the study will be graphically displayed thus permitting a visual comparison of the effect of each of the four parameters B,k,~, and R o 2.2 ANALYTICAL CONSIDERATIONS In the first part of this section the real c axis will be investigated as a possible neutral stability limit of slightly unstable modes. The latter part will consist of a derivation of a necessary integral condition on the eigenvalues by which the complex c plane will be separated into allowed and forbidden regions for c. The first aspect of the study is related to the filtering properties of the model. It will be shown that if R k < 1, which covers the meteorologically interesting range of wave numbers even for a moderately large Rossby number, fast moving slightly unstable modes, essentially for c outside the interval 0 < c < 1, are not present in the model. Fast moving neutral modes, however, are not ruled out and may be present. The entire c axis with the exception O < c < 1 is thus not a neutral stability boundary. A line or surface is - r -

20 said to be a neutral stability boundary if it contains eigenvalues which are the limiting case, as c. approaches zero from above, of eigenvalues with c >0. The existance of fast moving highly unstable modes can not be dealt with directly and is partially treated by the integral necessary condition. The integral condition may be considered a relative to what has come to be known as a semicircle theorem. It will be shown that the numerical value of a certain definite integral, as a function of the parameters k,.,R,B and the assumed value of the eigenvalue c, determines whether the tested c may be an eigenvalue. It is not sufficient for the condition to be satisfied to render c an eigenvalue. A sufficient condition can be formulated only in terms of at least one analytical solution of the characteristic equation. The results obtained in this chapter will prove most valuable in the numerical phase of this study. The search procedure employed to determine eigevalues numerically will be considerably simplified in view of the large excluded regions of the c plane. Before proceeding to the analysis it is advantageous for both analytical and numerical purposes to recast the eigenvalue problem in canonical form. Let w(z) in equation (2.9) be written as: w(z) = g(z)w(z). Substituting in (2.9) and requiring the coefficient of w to vanish leads to z (see appendix) + S(z)w = 0, w(O) = w(l) = 0 (2.11) zz

21 where P [1-kR (z-c)] 1 o 1 iu 1 iR g(z) = (z-c), p = _ - + p =- + - [l+kR (z-c)]2 1 2 2k' 2 2 2k and B(k2+ 2)+2k2R2 R (k2 2) S( z ) + S(2) 22, - _, -2-,.__.22 2 (z-c)2 -k2R2(z-c)2 [1-k R2(z-c) 2 0 0 Equation (2o11) is observed to have three regular singularities. For c = 0 these are at: 1 1 Z = C Z c + -- = c - -- r r Rk' r R k o o The existance of regular singularities on the real c axis suggests a possibility of utilizing power series expansions for formal representation of solutions. In the subsequent discussion the real c axis will be separated into its singular and regular domains. In the singular intervals the power series are applicable while the regular domains will be studied through a variational technique. 2.3 THE SINGULAR REGIONS IN THE c PLANE AS NEUTRAL BOUNDARIES The singularities of equation (211) are on the real c axis. They are defined by the three cuts 1 1 1 1 I. - - < c < 1 -- II. 0 < c < 1 IIT. - < c < 1 +Rk - r- Rk - r - R k- r - R k 0 0 0 0 For R k < 1 the three cuts are mutually exclusive as is shown by the heavy segments of the c axis in Figure l;considering the indicial equation (c.f. Ince r

22 INTEGRAL CONDITION TEST C. k= I.5- -- B= 1.5 R=0.5 0.4.3.2 Singularities for Ro k < I / Domain of Possible Eigenvalues 1+ Rok Rok C Rok Rok Cr Figure 1. Separation of the c-plane by the integral condition.

23 p. 160) corresponding to equation (2.11) the exponents relative to the three singular points are found to be (see appendix): relative to I: 1 I 1! ~H n= (li ), n =- (1-i ) + 2 k 2 k' relative to II. n -1, n = 2 1 2 relative to III as in I. In regions I and III power series expansions are sufficient for formal representation of solutions. In region II, where the exponents differ by an integer, a logarithmic term must be incorporated. If the three singular intervals of c are mutually exclusive and ci = 0, equation (2.11) has one singular point as z transverses its interval of definition 0 < z < 1.'he power series expansion will then be valid for the entire range of z thas enabling us to apply the boundary conditions at the end points of z. (If there is more than one singular point a power series expansion around one singu lar point is valid only as far as the next singularity). Consider the case R k < 1 and c in region III. For I o the general o r solution of equation (2.11) can be expressed as: w =A [z - (c - I- n f + A Z - (c -— ) f1 (212) i r' L n (2.12) w 1 Z`Ir Rk' + 2 r Rk -' o 0 1 where f and f are regular power series in T = z - (c - — ), n and n are + r R k + the exponents relative to region II, and A and A are constants to be determined from the boundary conditions. In orde-r to use the boundary conditions branches of the multivalued functions must be specified. Select

24 T = T for > 0 -ix T = I Te for r < 0 The boundary conditions can now be applied leading to the homogeneous algebraic set: n A T n+f +() + A T -f () = 0 at z = 0 T l+f (1) + A rn- f (1) = 0 at z = 1. 1i1 + 2 - The condition for nontrivial A1, A is 1' 2 nn n n+ Tor T f (O)f ()f1) - -T f (O)f () = O. rearranging: n + IT- f (O)f+(l) 1 _ _ + n_ n+ f (O)f (1) T T + - T0 1 Simplifying the LH side by using the explicit form of n+, n and the branch under consideration we find: r-i 1 8i ~. 1 -i n+ n_ 1 - (c IT n+ T T n —n+ r ___ k _ 1 k T LoH - 1 -i 1 e k n. n+ T c - -- )e T Tr Rk - O 1 o o If we assume that c is in the interior of interval III we can further write r TtR - iLnb L.H. = e * e where

25 s =. — _- -1 > 0 C -- r Rk 0 Consider now the RH side. The conjugacy of n+, n carries with it the result that f and f are also conjugates of each other. i.e., f = f + _ + Miles (1961). Therefore: f (o)f (l) f(f ( f (() f (O)f (l) i R.H. = f()f ) = () (1) = f ()f ) = e for some real The LH side can have this form with modulus unity if and only if i = 0. But for the form of the solution (2.12) to be valid I cannot vanish (the exponents differ then by an integer). The contradiction thus shown can be traced to the assumption that c is r in interval III. We therefore, conclude that since no solution in Cut III exists as c +- 0, this interval of the c axis is not a limit of an unstable 1 r regime. A similar development can be applied in interval I. excluding that cut from a neutral boundary as well. The end points of intervals I and III cannot be analyzed by the present method and require special consideration. It will be shown in the next section that the intervals of the c axis devoid of singularities are not marginal to r unstable regions. Consider an end point of interval I or III. Approaching the end point from above and from below on the c axis, that is from regions that are not neutral stability limits the end point cannot be a neutral stability point. For if it were it would be possible to approach an arbitrarily close interior point from an unstable region contrary to what was proved. The above

26 argument leaves the interval 0 < c < 1 as the only possible neutral stability limit. 2 4 THE NONSINGULAR REGIONS OF THE c AXIS AS NEUTRAL BOUNDARIES r In this section the nonsingular regions will be studied. They are: 1 1 1 c < - 1< c c + -- r R k' r R k r R k o o o In these regions c is a continuous function of the parameters k,,R,B, henceforth denoted by p.o Consider the variation of c with respect to any of the dc parameters. If - is real as ci - O+- and some c, that point (c, 0+) is not pi r r adjacent to modes with c. f 0 and therefore does not belong to the neutral stability boundary. If - has a complex part at (e, 0+) that point is a limit 1dp. ~~r of unstable modes and belongs to the neutral surface. We are thus led to comac pute p- in the regions where c is continuous. To this end a variational method Pi can be profitably used. Equation (2.11) can be written compactly in the form 2 D E w + L- - + - w = 0 (213) zz m ( -m) where 2 1 2 2 22 2 ~ 2 (z-c)2, m = 2 D = mLB(k + ) + 2R2k2, E = m R (k + Ro o 0 Multiplying equation (2.13) by w, integrating over the range of z while invoking the boundary conditions we obtain:

27 2 D1 E2 12 - w dz - 2 f - dz + f(-_ + )w dz 0 (2.14) 0 z o o n-m2 (Tm) which can be rewritten as n 2= - = F D where F is a functional and 1 D E'2 l2 w I = f(- + 2)w dz - f w dz, I - dz. n o'n-m 2 z' D o T (n-m) If w is a solution of (2.11) satisfying the boundary conditions the functional F is stationary with respect to small variations in w. Therefore D F aF 5c ~op. +c dp. and the differentiations are performed only when c and p. appear explicitly. 1 We are now in a position to isolate the desired derivatives ~F ~i ac Pi bc Let the limiting point (ci-O+,c ) be a regular point. If c = O0 a real solution of (2.11) can be found. Since w is continuous in c any solution for ci approaching zero must approach the real solution in the limit. Consider now lim (F )F (2.15) c ~0+ a7' ) ac The parameters p. are real. In the limit ci+0+ w is also real, up to a multiplicative complex constant, since a real solution for w exists. Since c is r

28 a nonsingular point there is no source for complex contributions in any of the bci derivatives of (2.15). -- is therefore real and the conclusion that the nonpi singular regions on the c axis are not stability boundaries is proved. 2.5 A NECESSARY CONDITION ON THE EIGENVALUES Consider equation (2.11) rewritten in the form: w + S (z)w = -S (z)w (2.16) where S(z) = So(z) + SL(z) (2.17) and the splitting of S(z) is such that S (z) is not identically zero. Equation (2.16) can be transformed into an integral equation through the use of a Green's function. The integral equation is: 1 w(z) = -f G (z,)S1(I)w(S)d (2.18) where G (z,) is the Green's function corresponding to 0 Go + S (z)G = b(z-S),G (o) = G (1) = o Taking absolute value of (2.18), squaring, and applying the Cauchy Schwarz inequality: 01 0 1 1 o il(z)2 = C ( Z)s ( )w(e ) dl < (z,o Ils(e)z|> 1l dJ21s ( i)2|de. Integrating over the range of z we have: 1 1 1 (z) wd ~~~~~oS Wz wi ~ b (~L1 5()2dz

29 Transposing and factoring the constant w(r) 2dr For nontrivial w J|w(r)| dr>O. leading to the condition: 0o If - fIGO (z ) 21 S ()12d dz<l (2.19a) the value of c in G and S1 is not a possible eigenvalue. At this point S (z) and S (z) must be specified. A variety of combinations has been tried, all leading to rather prohibitive expressions. The simplest workable pair is clearly S (z) = 0 S (z) = S(z) The Green's function corresponding to this choice is readily found to be: G(z,) =((Z-1) < z) carrying out the integration of (2.19a) with respect to z: 1 2 22 1 2 2 2 JIG(z,) dz = fz2dz + f(z-1) 2dz = 1 2(1- )2 and we have the workable condition on c: 1 2 2 2d f d2(- )2 S(5)|2d 3 (2.19b) for the existence of a solution satisfying the boundary conditions. The analytical evaluation of the integral is cumbersome~ The alternatives are either an approximation or a numerical integration. With an eye to needs imposed by the differential analyzer a numerical integration is appropriate

30 since its accuracy will certainly exceed the best of approximations. At this stage of the study the c plane can be divided into possible regions for eigenvalues and regions where eigenvalues cannot be found. A typical result for given ki,B, and R is shown in Figure 1. The integration of (2.19b) was performed by using Simpson's formula. 2.6 THE NUMERICAL STUDY OF THE EIGENVALUE PROBLEM Two point boundary value problems are of sufficient interest and complexity to warrant a detailed exposition in this phase of the study. In solving the problem numerically there are essentially two methods of attack. (1) A search process based on an iterative algorithm whereby a crude eigenvalue is continuously changed until the boundary conditions are satisfied. This approach is delicate to program and requires considerable preliminary knowledge about the eigenvalues' confinement and dispersion~ The iterative process can be started by a first guess which is continuously improved upon, or by a scanning procedure through the domains known to contain eigenvalues. (2) The other possibility, which is more commonly found in the meteorological literature, consists of a finite difference scheme applied to the differential equation at points in the interval of definition. The finite difference scheme of a homogeneous differential equation leads to an algebraic system of homogeneous linear equations whose coefficients depend upon the required eigenvalues. The condition for the parameter c to be an eigenvalue is the condition for the existance of nontrivial solutions, namely the vanishing of the coefficients' determinant.

31 Results obtained by the first method, while often difficult to achieve are trustworthy and can be tested in a straightforward manner Results based on the determinantal condition may raise problems of interpretation. To illustrate the problem consider the finite difference version of equation (2.11) applied at the (N-l) internal points z = I.D W(I+i} + w{I-1.- 2w(I) )+w(-l) - (I)+ S(Ic)w(I) = 0 I = 1...N-l (2.20) D where N + 1 is the number of grid points including the boundaries and D = N is the grid size. (recall 0 < z < 1) The N-l equations (2.20) give rise to an (N-l)x(N-l) determinant which when expanded results in a polynomial in c. The degree of the polynomial and, therefore, the number of its roots depend on N and the form in which c enters into S(z,c). Since N is arbitrary, there arises a problem of the relationshp between the true eigenvalues and those computed. To illustrate the possible difficulty just mentioned, we refer to a case study made by Arnason (1967). It is shown in the study that near singular points of the differential equation the finite difference approximation method results in true eivenvalues when the singularity is approached from one side but not from the other. In view of the shortcomings of the second method coupled with the rather cumbersome expression for S(z,c), the first method was adopted. The results of sections 2.3 and 2.4 can now be applied and the general area of the eigenvalues confinement can be considered as known. There remains however the possibility of the existance of more than one eigenvalue for a given set of parameters k, Q, R, B. If more than one mode does exist, and this aspect is dif

32 ficult to investigate, there is a possibility that during the process of modifying an eigenvalue an oscillation may occur between two neighboring modes and the numerical process will not converge. To overcome such a possibility the changes in c can be programmed to be very small (inefficient) but the question will remain unanswered. At this point it becomes clear that if analytical methods are not available to study the existence of more than one mode and the computing time is to be kept to a reasonable amount, resort should be made to technical methods which fortunately are available and most suitable for this type of an investigation. The tool in mind is the differential analyzer (analog computer). Its merits in the study of eigenvalue problems are striking) the greatest advantage being the immediate contact and control. the programmer has with the problem, program and tool. The advantage of using the analog computer is two fold. Not only can the possibility of the existence of a spectrum be settled by inspection, actual. eigenvalues can be read and determined directly and where accuracy is not too demanding this constitutes the final solution. Unfortunately the problem under consideration is, for some combinations of k, I, B, R, quite sensitive to changes in c which the analog computer cannot resolve. The results of the analog solution can, however, be trusted to within a few percent of the true solution and thus constitute excellent initial guesses. 2. 7 ANALOG COMPUTER SOLUTION OF THE EIGENVALUE PROBLEM The most suitable technique for invesigati the most suitable technique for investigating the solution of equation (2.11) as a function of the parameter c is to employ the computer in the rep

33 etitive operation mode. In this mode the differential equation is integrated many times per second and the solution is displayed on a cathode ray tube in graphical form. Changes in c are instantly observed and can be affected until the boundary conditions are satisfied. Before proceeding to the analog program it is valuable to reconsider equations (2.9) and (2.11) from the point of view of nonlinear equipment and ease of coefficient generation. The advantage of equation (2.11) from the standpoint of minimizing the number of multipliers is quite clear. A closer examination however, indicates that the first order pole at z = c involved in equation (2.9), as compared with the second order pole of S(z,c), outweighs the economy consideration, particularly in regions close to neutral stability. This conclusion was also verified in practice, the first order pole resulting in a much smoother, and therefore more accurate, curve than the second order one. It was therefore decided to return to equation (2.9) for the analog computer search program. 2.7.1 Method of Solution 1. Equation (2.9) will be separated into its real and imaginary components thereby resulting in two second order coupled equations with real variable coefficients. 2. The coefficients of the system, which are functions of c, ci will be generated in such a way as to enable c and ci to be varied independently by a single potentiometer each. 3. A computer program for solving the system of equations will be utilized

54 in the repetitive operation mode, c and ci will be varied by their respective A potentiometers until the boundary conditions w (1) = wi(l) = 0 are satisfied. A A0), (O) The initial conditions are chosen as w (0) = w.(O) = O0 w (0) = A, w'(O) = B r i r i 2 2 with A + B2 O A and B arbitrary convenient constants. The arbitrariness of ABis a consequence of the homogeneous boundary conditions. When the boundary conditions are met c and ci are true eigenvalues, within the accuracy of the r i program, and can be read off their respective potentiometers. It can be anticipated that the program schematic for integrating the system of equations will be of considerable size. It is also observed that the sequence of steps in the parameter study is rather lengthy. In view of the large number of coefficient potentiometers required in the program and the number of combinations of the data parameters k, ~, B, R, a hand computation of scale factors, potentiometer settings and static voltages checks is prohibitive. The observations mentioned above lead to the following convenient sequence of steps: 1. The variables that are prescaled are z9 c, c, only, the latter two through the application of the integral condition thus determining their maxima for each set of parameters. 2. An unscaled schematic is written leaving all other scale factors unspecified. 3. A digital program is utilized to compute the following: (a) maxima of unsealed output variables of all amplifiers as functions of ci, c and z. (b) scale factors based on computed maxima (c) potentiometer settings based on computed scale factors (d) static voltages at key points of the schematic based on all previous computations.

35 The method just described proved very efficient and very reliable. With the schematic shown in Figure 2 eigenvalues can be obtained quite readily. In addition the existance of more than one mode for given k, ~, B, R has been 0 ruled out for a number of data sets. The Analog Computer Program. Equation (2.9) can be written (dropping hats) w + Pw + Qw = O zz z where 2i~R 2iRo 2 P = P + i P 1-k2R2(z-c (z-c)(l-k2R2z-c )2 r 0 0, B(k2+ ) 2R = + i Q. l-k 2R(z-c) (z-c)(l-k2R2(z-c)2 r 1 0 0 Let primes represent differentiation with respect to z, then w" + iw' + (P + iP.)(w' + iw) + (Q + iQ.)(w + iw.) = 0. r i r 1 r 1 r 1 r 1 Upon separating real and imaginary parts: = - P w' + P.w' - Q w + Q w r rr 1 1 rr r 11 (2.21) = - P.w' - P w' - Q.w - Q w i 1 r ri 1 r ri The initial conditions for integrators are: wr () = w.(0) = 0. r1

0. 1 2C. 2 +100V - 1.0 1. 0 I OO1.0 01 +11- | O1 — l/ \ p^q1.0 1. 0 1. 01, GP 100Z TR 1A1 - — D L^OOD,.1.0 _1 0 Cr+ -I 1 B5 1 0 (100z from Amp. 1B2, GP 1B8) G 10OG G -G200(z - Cr) +G200(z - Cr) +100V 0 200(z - Cr)1B7 1B4 RK RK 2 2 G200(z - Cr) G - 100 M) -G200(z - Cr) 100 NI V>2 -0 V -1000G -100G 3 G200(z - Cr) Observe Gain RK G2 R2DK2 R -20 * 100(D + R2K D) GlOOD +G1OOD -100G L IOO)G I 4 -G200(z - Cr) RK 6 19 Amplifiers 9 Pots Figure 2. Analog computer schematic for eigenvalue search.

37 -G1OOD -G100D+ G100D_ G100D G ~ 100 Mli G - 100 NV8 M+GM 100I D -GM -IO +GM ~ 100 DG G M GM0100M GM 100i10 O 1 1. 0 OR 181 * (_ Ci) G (D D ) NI 1.0 1. 0 +GM 100 Ml V GM.1 00 l r^ -TR 1B1R G M 100 N ~- — O. 10+GM 100 ( Z — Cr.) -=IGM ~ 100 - +GM ~ 100D -GM ( + O. 1 G1000 GIOODOOD+ V9 Vlo 200OC. 2200Ci 62G02 ci 62 62 200' 200 62 C.LO0 1..O G2 1.0Ci G D 00 10 TR lC7 TR 1C9 6 C.1o 0+ -1. C. C 1.0 -G2 Ci 1.0 -200 A200C 0.1 0.1 TR 1C8 TR IC10O'^J 20 -- G __\D G2 C_ LO -G2 C 1Ci LO- - -200- -2006 D0 6 0, -200D -200D 2000 200D TR 1A9, [ V12 TR 1A10 0.1 0.-CI) 1 Ci 60GRG.100 (z-C 6-16 C) 100 20 Amplifiers 6 Pots 6 Multipliers Figure 2. (Continued).

38 P( 2(z - Cr) -RK Di GM 100 -l - PD ) O GM.100( —^(l J —-^- KG KR -zR \ (V - Kp P D r P D r GR 50 -G2 C( Ci c+) K G I 2 200G2 Q = (Z - Cr IKR2M1 N1+ KR i Qi D D_ D+ - + D+ GM.100( ) (M N10 ). \ KU i /. 21R 100- GM \ V Q = KQ 0 (z -~ C r) C —- ~ -K-Q / 100GR / -G2 C( i C ) _F / G 0^*^ K-U.G-KR y 200G2 2C i M1 N1 2K 2 C Ci P=: —-D- + RI D ~ -, + + )+ 4~2 ji + G ('D D+0 KG.R2K2 \ 100 1 G2 Vp - Kp P. ID Ci_ -K. * P. GI-~ 1002 P+ 100 GI / -GM 50 1 + N / 50 ~GM Q, 21R 2Ml + N1 2_ D D3 ) I Q G2 \ V KQr (QR -Qr) 100 Gl / 50 ~ GM Figure 2. (Continued).

-K" WI,0.1 +Kw'iWI K Pi -K W 0.1 r "W.' i w — -Kw''WQ W'r i KQi Q -KW Kw' ri\',.,mK.C 100 \ KW,,.-w,'W -100V -KQ L= KW KPi Kw. \ OK r 01, +K Wi ~~~~~100-K Q P' -KWr -K W 0. PW1 Wr + -W P r -PrW r K; 0K1 -KW fr RT Wr Qr K W - 1 Q Q I WI l KW W. Wr O K.! KQr' KWr K Qr W -K W K W —---- 100 rr Q-K - I rr Wr ~~ iu r 2, r C d \-K,__.Fg KWr' 2. K W'0.1 W, K Pr P r/ i +Q~^ r - 100 * 100Q - KK K wi W i' -'' PWl'rKVO. W,__ -. Kr rWi, K W KrI r 2 (Co 10 0cl W. -KP'r W'io KWir 100 i K,W KWW K K —----- RT KW'. W'r 0. 1 -KWM i t l e o.) i l 100 QiWi -KQ Qi KWr'W Figure 2 (Concluded).

40 Since the boundary conditions are homogeneous, initial conditions on the integrators of w' and w' are arbitrary, i.e., r 1 w'() =A, w'(O) = B, A + B 0 r 1 Using the following notation Ml = 1 - kR (z-c o r N1 = 1 + kR (z-c ) o r D = (z-c )2 + c r i I) = 1 - 2kR (z-c ) kR2D O r o D = 1 + 2kR (z-c ) + kR2D + o r o 2 B(k2+ 2) 2 P 9 P.i Q Q. are found to be: r r 1 2(z-c ) M C r Ml i N1 P = + KR [- - + R + - o R ] r D o D o D D oD _ _ + + 2c c c i Ml 2 i N1 2 i. = - + -- + k R - - + k R - i D o D oD D oD _ + iQ = MIR i i C. M1. N Q 2 - ~ i- (+ --- 2 R2( - + - + o P o P oP D - + - + The schematic in Figure 2 shows the generation of the equations' coefficients and the analog circuit for integrating the system (2.21).

41 The analog computer investigation has led to initial good guesses of eigenvalues and the following findings: 1. To a given set of data, k, ~, R, B there corresponds one and only one mode for c's satisfying the integral condition. 2. If cO is an eigenvalue corresponding to k, X, B and R = 0, and c corresponds to k, i, B and R i 0 the imaginary part of the difference c1 - c may vary by as much as two orders of magnitude for different combinations k, ~, B, while the real part of the difference is always very small, within the order of the error introduced by the circuit's limitations. 2.7.2 The Digital Solution The information and first approximations of eigenvalues obtained in the analog study provide us with the necessary guide lines in programming the digital part. The absence of more than one mode per data set and the change (c1-c0) as a function of k, i, R, B led to the following sequence of computations. 1. Select k, ~, B. coordinate be of the order Im(c.-c ) where c. is the approximate eigenvalue 1 0 1 for Rol i 0, and Re (c -c ) be the grid size in the c direction. 4. With R = 0 integrate equation (2.11) forward in z up to z = 1. The theoretical value of w (1), w.(l) should be zero. Truncation error, however, will cause a residue at z = 1, the truncation error for a centered difference

42 1 -2 -c scheme being O(- w"). This residue, which for d = 10 was of the order 10 d is next used to determine fulfillment of the boundary condition with R i O. 5. Set R = R 1 corresponding to c. Integrate equation (2.11) at each 1 of the nine points (- + ic ) and the eight surrounding neighbors. Out of the nine sets of w (l), wi(l) choose the set with the smallest absolute values. Repeat (5) with the c corresponding to best end values. 6. When no reduction in the values at z = 1 can be obtained halve the grid sizes and continue. 7. When the end values are less then or equal to the residue for R =0, o the corresponding c is an eigenvalue. 8. Use last eigenvalue with initial grid sizes and repeat step (3) to (7). With R = R repeat (5) to (7), etc. The above sequence proved very efficient and the computing time for about 1.5 eigenvalues was only 2 minutes on an IBM 7090. The results are displayed and discussed in the next section. 2.8 RESULTS OF THE BAROCLINIC STUDY 2.8.1 Growthrates and Propagation Speeds The major unknowns of the stability problem are the growthrates of sinusoidal disturbances, assumed to be initiated by an unspecified mechanism, and their speed of propagation relative to the embedding current. The reality of the model, in relation to the physical assumptions introduced in its formulation can be studied and tested against the observational gross features of the atmosphere. To this end typical vertical velocities, pressure and temperature

43 perturbation fields can be computed and plotted and their mutual geometrical relation compared with representative observational conditions. We begin with a discussion of the growthrates and propagation speeds resulting from the numerical study. The first and most striking result is the finding that all waves, regardless of their lateral extent and for all static stabilities and Rossby numbers investigated, travel with the same speed, c = 0.5, from west to east. It is observed that this is also the propagation speed in a quasigeostrophic formulation and is equal to the mean speed of the basic current. Since true atmospheric wave motions are certainly dispersive and since the linear theory applied here is certainly not indicative of quasisteady state conditions observed, it is impossible to conclude whether the missing mechanism which would bring about dispersion is the beta effect or the nonlinear nature of the primitive equations. In this respect the ageostrophic effects do not improve on quasigeostrophic theory. The ageostrophic effects on the growthrates kci were studied as a function of the Rossby number (considered to be the main parameter). Beginning with a pair of wave numbers k, i and a given static stability parameter B, the quasigeostrophic complex part of the phase speed is computed from C = [-(1/4 + - - coth)] /2 52 = B(k2+ 2) (2.22) oi 02 6 For each fixed set of k, X, B, the Rossby number R was varied in increments of 0.1 and the new ageostrophic growthrates determined numerically. The results are given in Tables I through III.

44 TABLE I TABLE OF EIGENVALUES AND GROWTH RATES FOR I = 0.5 Ro B k 01 0.2 0. Eigenvalues 0.5 0.279118 0.2791.57 0.279147 0.279167 0.5 1.0 0.264886 0.264761 0.263980 0.262691 1.5 o.241857 0.241233 0.259539 0.236136 0.5 0.269678 0.269717 0o269716 0 269755 1.0 1.0 0.261858 0.24164 0,240940 0 2 39788 1.5 0.195869 0. 1955313 0 1,93614 0o1.90133 0o5 0.260335 0 260374 0 260584 0 260423 1.5 1.0 0.218939 0o 218743 0o218099 0.217026 1.5 0.147192.o 146621 0.1.449053 0o141980 Growth Rates 0.5 0.139559 0o 139578 0Oo 159575 0.139583 0.5 1.0 0.264995 0. 761 2676,65980 0o262691 1.5 05362787 o.561849 0.559008 0 354204 0o5 0o134839 o 134858 0. 134858 o 154877 1.0 1.0 0.241858 0 241645 0o240940 0 239788 1o5 0.293803 0.292969 0.290421 0o286099 0.5 0.130167 0ol,0187 Oo150192 0130211 1o5 1o0 0.218939 0o21.8745 0o218099 0217026 1.5 0.220788 0 219951 0 o 21735 0 212970

45 TABLE II TABLE OF EIGENVALUES AND GROWTH RATES FOR I = 1.0 Ro B k _B k 0 0.1 0.2 0.3 0.4 0.5 Eigenvalues 0.50 0.264995 0.265295 0.266145 0.267589 0.269545 0.272095 0.75 0.259173 0.259360 0.259891 0.260766 0.261954 0.263423 1.00 0.251068 0.251118 0.251218 0.251368 0.251518 0.251568 0.5 1.25 0.240710 0.240585 0.240179 0.293429 0.258304 0.236679 1.50 0.228102 0.227790 0.226759 0.225009 0.222415 0.218821 1.75 0.213201.212654 0.210967 0.208061 0.203795 0.197936 2.00 0.195869 0.195082 0.192632 0.188419 0.182198 0.173316 0.50 0.241858 0.242158 0.243008 0.244458 0.246433 0.248958 0.75 0.230393 0.230595 0.231158 0.232096 0.233377 0.234955 1.00 0.214350 0.214425 0.214625 0.214912 0.215912 0.215637 1.0 1.25 0.193538 0.193495 0.193209 0.192725 0.191975 0.190850 1.50 0.167317 0.167005 0.166243 0.164829 0.162712 0.159747 1.75 0.133868 0.133356 0.131791 0.129064 0.124986 0.119261 2.00 o.o86883 0.085858 0.082708 0.77023 0.067992 0.053570 0.50 0.218939 0.219239 0.220114 0.221564 0.223564 0.226126 0.75 0.201674 0.201893 0.202502 0.203502 0.204877 0.206580 1.5 1.00 0.176998 0.177098 0.177367 0.177798 0.178348 0.178973 1.25 0.143268 0.143218 0.145067 0.142750 0.142238 0.141433. 50 0.094140 0.093843 0.092937 0.081351. 088929. 085445 Growth Rates 0.50 0.132497 0.132647 0.133072 0.133794 0.134772 0. 136047 0.75 0.194380 0.194520 0.194918 0.195584 0.196465 0.197567 1.00 0.251068 0.251118 0.251218 0 8 5 168 0.251518.251568 0.5 1.25 0.300887 0.300731 0.300224 0.299286 0.297880 0.295849 1.50 0.342153 0.341685 0.340138 0.337513 0.333622 0.328231 1.75 0.373102 0.372144 0.369192 0.364107 0.356641 0.346388 2.00 0.391738 0.390164 0.385264 0.376838 0.364396 0.347238 0.50 0.120929 0.121079 0.121504 0.122229 0.123216 0.124479 0.75 0.172795 0.172946 0.173368 0.174072 0.175033 0.176216 1.00 0.214350 0.214425 0.214625 0.214912 0.215262 0.215637 1.0 1.25 0.241922 0.241869 0.241511 0.240906 0.239969 0.258562 1.50 0.250975 0.250507 0.249364 0.247243 0.244068 0.239620 1.75 0.234269 0.233373 0.230634 0.225862 0.218725 0.208707 2.00 0.173766 0.171716 O.1654i6 0.154o46 0.135984 0.107140 0.50 0.109469 0.106619 0.110057 0.110782 0.111782 0.113063 0.75 0.151255 0.151420 0.151876 0.152626 0.153658 0.154955 1.5 1.00 0.176998 0.177098 0.177367 0.177798 0.178348 0.178973 1.25 0.179085 0.179022 0.178821 0.178437 0.177797 0.176791 1.50 0.141210 0.140764 0.139405 0.137026 0.133393 0.128167

46 TABLE III TABLE OF EIGENVALUES AND GROWTH RATES I = Ro _B k 0.1 0.2 o.3 Eigenvalues 0.50 0.083229 0.0890089 0.105065 0.128268 0.75 0.072968 0.079325 0.096239 0.120107 1.00 0.055947 0.063603 0.082720 0.108032 1.15 0.039351 0.049351 0.071382 0.098570 Growth Rates 0.50 o0.41614 0.044544 0.052532 o.o64134 0.75 0.054726 0.059494 0.072179 0.0900ooo80 0.5 1.00 0.055947 0.063603 0.082720 0.108032 1.15 0.045254 0o056754 o.082089 0.113355

47 The ageostrophic deviations from the quasigeostrophic values of growthrates are displayed in graphical form in Figures 3 through 6. The ordinates represent percentage deviations obtained from c - c Ageostrophic percentage deviation = 100 x c oi where c o is computed from (2.22) and ci corresponds to the numerical value computed for a nonzero Rossby number. Figures 3, 4, 5 indicate that for large lateral wave configurations, namely i1. the ageostrophic effects are selective, stabilizing the shorter zonal waves, k > 1, and destabilizing the longer ones, the separating value of k being slightly greater than unity. The above observation holds for the three values of the static stability investigated9 B = 0.5 1) 1.5, with the maximum effects occuring for B = 1. This is not to say that the actual maximum is at B = 1., but somewhere between 0.5 and 1,5J It is further noticed that the percentage deviations are, for these values of I and B, of the order of a few percent even for large values of the Rossby number. The behavior of the shorter lateral disturbances, I = Tc is shown in Figure 6. The maximum zonal. wave number k that could be used with a static stability B = 0.5 was k = 1..2, the limitation resulting from equation (2.22). It is noted that the zonal short wave cutoff for a model incorporating lateral variations in the perturbation quantities, occurs at a longer zonal wave length than is found in the absence of y dependence. The ageostrophic deviations in the present case are found to be much more sensitive to the Rossby number9 reaching more than l0o variation for moderate R and short zonal waves. The o

5 K =0.75 ~~~~/ K\~~~o.~ _=1.75 ~fo^ \I - O 2 - -1.00 415 0.2-3 —F?040.5 Figu~re 3- Percent deviation of ageostrophigrothrates for 1, 0.5

K= 0. 75 5/F —~~~ =.50:: X-.K=.75 LU 155 K2:O0 0.0 8.1 0.2 o 0.4 r nt dof ageostrophic growthrates for

K= 0:;'7 0.50 — ~~~~ ~~~ — 5~ 1 -~~~~~~~~~~1. 2~~~~~.5 LU 0. 0.2 0.3 0.5 Figure pr 0 0.5 c cet deviatio, Of eot 8geostrophic grO-thrat

51 150 k =1. 15 100 k =1 k =.75 / z k -.50 B=.5 50 0.1.2.3 R Figure 6. Percent deviation of ageostrophic growthrates for ~ = x, B = 0.5.

52 effect is not selective and all zonal disturbances are destabilized. It is generally accepted that a characteristic Rossby number for mid and high latitudes is of order 0.1. In lower tropical regions where the circulation is less geostrophic a larger Rossby number is more representative. Recognizing the limitations of this model, in particular the missing beta effect which is largest at equatorial regions, it is of interest to note the general behavior of the ageostrophic influences. Since middle and high latitude flow corresponds to small Rossbynumber and low tropical flow to a higher value, a partial explanation for the birth of violent instabilities in low latitudes as opposed to their absence in mid and high latitudes is borne out by the model. 2.82 Cross Sections of the Perturbation Fields Once the eigenvalues of equation (2.9) are known the determination of the various perturbation fields becomes an algebraic problem. The boundary conditions on the vertical velocity component are sufficient to determine u, v, p, T (or p) up to the arbitrary multiplicative constant incorporated in w. Since w is considered known the system of equations can be written (from (2.lb)... (2,5b)) ikR (z-c) -1 ik 0 u -Rw 0 0 1 ikR (z-c) i~ 0 v 0 0 = ik ii 0 0 p -R w 0Z 0 1 0 ik(z-c) p BW the solution of the system is readily found to be

53 p i([l-k R (z-c) w + [i~R + k R (z-c)w) /. k(k + )(z-c) (2.24) 0 Z 0 0 p = -i[[l-i~R (z-c)]w + [i~R + B(k +A )(z-c)]w) 7/. k(k + )(z-c) (2.25) 0 z 0 and it can be verified that the hydrostatic equation is satisfied. We are also interested in the temperature field. To this end, recalling the constancy of density in the horizontal momentum equations, the Boussinesq approximation p = p (l-a(T-T )) (2.26) will be used. Tn order to obtain vertical cross sections of w, p, and T perturbations, the real or imaginary parts of rw(zn i[~y +k(x-ct )] p(z) e T(z) are plotted, where y and t are fixed. o o As for horizontal distributions of the mentioned variables the plotted fields are given by the real or imaginary parts of: w(zo i[~y+k(x-ct )] p(z ) e Since the coefficient of the exponential is a constant for a given height all fields will be similar except for a phase displacement. It is thus sufficient to compute the phase relations. In geostrophic theory the lateral heat and momentum transports can be obtained from the pressure and temperature fields

54 alone, since the meridional velocity is related to horizontal pressure field in a simple way. In this model the heat fluxes cannot be inferred directly, nonetheless some phase relations will be computed and will serve, in addition to giving some more insight into the details of the model, to check the numerical computations. In Figure 7 some representative numerical results are shown. As a check the case R = 0 was computed both from the analytical solutions: 0 w = (l-yc)((z-c) - 1)e7 + (l+c)((z-c) + -e 7 7 p = 17i [(l-yc)eY - (l+yc)eZ] ( k(k2+i2) +227). 2 T = ----- [(l-yc)eY + (l+yc)e Z] k(k +) and from the numerical scheme. The results obtained by the two methods were found to be in complete agreement. The phase difference between the temperatureand pressure fields was defined by -1 Im(T) 1 I P) phase lag = tan - tan R (T) R (p) and was plotted in Figure 7 as a function of height. In the present formulation of the model the separation technique results in an arbitrary multiplicative factor in all the variables. Consequently a quantitative comparison of actual magnitudes as a function of the parameters, in particular the Rossby number9 is not meaningful. To circumvent this arbi

I" \ \ \ ^'' Q=f77- Ro.0 \ \ \ \ R -.I Ro =.3 -=I Ro -0 0.5 \ \ \ \ \ \ \ \ 0 O 71T 2 PHASE LAG Figure 7. Phase lag between temperature and pressure fields.

56 trariness one must resort to a comparison of normalized fields. In plotting vertical cross sections of w, p, and T each of the variables was therefore normalized by its maximum absolute value. Figures 8 through 10 show some typical vertical cross sections of vertical velocity pressure and temperature perturbations. Each figure corresponds to two values of the Rossby number for comparison purposes. The isolines of the plotted variables correspond to a normalized intensity scale designed to portray the zero lines and the areas of maximum intensity. The vertical velocity wave is seen to consist of two cells one of rising the other of sinking motion. The temperature field agrees with this vertical motion inasmuch as the warm air cell ascends and the cold cell descends. The pressure field has associated with each cell two sub-cells, or a high and low, located in the upper and lower levels of the basic current. A typical result is a low level cold high lagging a sinking cell while a cold low aloft leads the sinking air. This mode of operation is in broad agreement with typical synoptic conditions in midlatitudes. The slope of the vertical systems is a strong function of the lateral wave number and the Rossby number. While for long lateral waves, I = 1, the slope is always from east to west the occurrence of the opposite tilt is found for short lateral waves, I = -r. The effect of the Rossby number on the slope is somewhat unexpected; an increase of R tending to make the systems more vertical. A collection of some typical results is shown in Figure 11. The slope referred to is the axis of symmetry of the vertical velocity cells, or the zero line of vertical velocity. In quadrant (a) B and ~ are fixed and k, R are varied it is observed that the short zonal waves, k = 1o5, are associated with a more

- >< ~ \ \ wlo \ l w' O \ \ \ I U\, // \ z L,\ I I I \\ wI= O U.J 0 7r 27r x -= I k= B=.5 Ro=.I —- Ro=.3 Figure 8. Vertical cross sections of vertical velocity perturbations.

I I I /I I0~~~~~~~~~~I 0 2 - - I I II ~ T I T'< 1 T'=O T'O I I I I =.5 =.I R =.3 I II I III I I I I I Ij I Ii I i 2-1 k- I B-.5 R. —- R o-.

H V\ \L \\ p'>o \\ \\I \ /' p'<o' \ \ 11/ N 0 pO p 2O x \ =1 k=l B=.5 R=. I —I R = 3 Figure 10. Vertical cross sections of vertical pressure perturbations. 00 0 Figure 10. Vertical cross sections of vertical pressure perturbations.

60 ki\ \\.5 R0:.3 \ f k=-.15 Ro.3 k\ \\ 1.5 Ro 0 \ k.15 Ro 0 k\\ = I Ro.3 k I Rb.3 k\ I R,0 0 k\ I R' o0 k *5Ros=. \3 5 Ro.3 k.5 Ro= I k.5 R B=.5 B=.5.-.5 I=.5 B=1.5 k=.5 B \ B.5 k=.5 Roe. \ Ro\ o.3'B=1.5 k=1.5 \B8.5 k=.5 RoO \ \ \ RO = I kB 1.5 \ \ k= Ro a.3 \ Ro =.3 8Bl1 k1.5 - B=.5 kl Ro =0 \ Ro=O B=.5 k=1.5 \ \.5 k=1.5 Ro:.3 3 \ \ Ro.3 -B=.5 kl.5 \ \ \ -.5 k1.5 R1O Ro=O c d Figure 11. Tilt of vertical velocity cells as a function of k, I, R, and B.

61 vertical structure than the long waves. It is also noted that for increasing Rossby number the tilt of the vertical velocity and pressure systems becomes less westerly, reaching an easterly tilt for large Rossby numer. The effect of the static stability and zonal wave number k is similar as is evidenced by quadrants c and d. As for the phase relation between the temperature and pressure field Figure 7 indicates an antisymmetry about the mid stream. The phase difference is - at z = 1/2 and increases towards the ground while diminishing upward. This 2 behavior is found in the lower stratosphere and near the ground and is in good agreement with observations. The effect of increasing the Rossby number for a particular set k, ~, B, is to increase the range of variation of the phase difference. Since our model is too crude to recognize a tropopause and does not have a jet stream incorporated, the matter of heat and momentum fluxes is not further pursued.

CHAPTER III THE BAROTROPIC GULF STREAM MODEL 3.1 FORMULATION* We consider a two-layer fluid on the P-plane, the density p of each layer being constant. Let x, y, and z measure distance to the east, the north, and the vertical, and let subscripts 1 and 2 refer to the upper and lower layers, respectively. The layers are of finite depth, with depths D (a = 1, 2) in the absence of motion. We take the upper boundary of the fluid to be the free surface z = + D D + (Zp/p )r), where Ap p - p, the interface between the lay1 2 21' 2 1 ers to be z = D2 + 2 (1/P2 )l and the lower boundary to be the rigid plane z = 0. Also, we assume shallow water theory to be valid and denote by q the horizontal velocity and by D/Dt the material derivative, t- ( + q a V) dt 5t ct Since the motion is assumed hydrostatic we have P1 P2 -+ g Z - C (x,y), p + g Z C2(xy) p1~ 1p2 2' The conditions p = at the free surface and p p at the interface imply CL(x,y) = g[DL + DP + P C (x,y) = g[D1 + p D2 + 2 The pressure force can therefore be written as xExcept for minor changes the formulation is a verbatim copy of an unpublished paper by Jacobs. 62

65 - VP = g'V n. p a a where g' = gipp/p is the reduced gravity. The assumption of constant density in each layer and no mixing between the layers leads to the continuity equation in terms of thte cpVhs. The equations describing the flow are then Dq + f k x q + g'V =0 (5.1) Dt a a a (no sum convention) Dh + h V = 0 (35.2) Dt a a a where the depths h and h are given by 1 2 h D +r h =D +' i 1 l 2' 2 2 2 p2 l () and where f = f + Py is the Coriolis parameter. 0 Our aim is to study the stability of small perturbations to the flow = (VU(y/L),O), q = 0, V being a characteristic velocity and L a characteristic length. Two approximations will be made at the outset. First, we will neglect entirely motions in the lower layer and thus eliminate potential energy conversion as a source of instability. This is justinied in the quasigeostrophic case, provided that a source of kinetic energy is present and that D < D (Pedlosky, 1965, Sec. 4). We assume that it is true also for the ageostrophic case. Secondly, we

64 restrict our attention to horizontal length scales so small that L2 < V/P. It is then permissible to neglect Py next to f provided we simulate the trapping effect due to the earth's curvature (Jacobs, 1967). This will be achieved here by supposing the fluid to be confined between walls at y = A, y = LB, B>A. It is convenient at this point to scale the variables. Omitting subscripts, which are now superfluous since q = 0 by assumption, we define nondimensional variables by (x,y) = L(x*,y*), t = (L/V)t*, q = V q*, n = (f VL/g')q*, h = Dh*. (.4) The nondimensional equations obtained through use of this scaling are, with asterisks omitted, Dq +A R + k x q + V = (5) o Dt and Dh -- + hVq = 0, (3.6) Dt where h = 1 + R 7 2 (.7) The nondimensional parameters are the Rossby number, R = V/(foL), (.8a) and a nondimensional radius of deformation y, given by 7 = foL/(g'D)1/2 (3.8b)

65 If we now take q = (U(y) + u, v), n = ~ + cp, (3.9) where V' = -U, and neglect products of the perturbation quantities3 we obtain the linear equations F3xd (9a10)cp o (t +U ax )(u,v) + (-(VU) + ( ax = 3 -10 ) and 2 ( - u (Hu) 3(Hv) (,o0 (a + 6) + += (3t 1) where H, the unperturbed depth of the layer3 and 5, the unperturbed total vorticity, are given by H = 1 + R 7, = 1 - R U (U.12) 0 0 As is usual in stability problems (see3 however, Case3 1960)9 we bypass the initial value problem and instead separate variables according to F ik(x-ct) k i =c F(x,y,t) = F(y)e, k 0 c = c + ic (5.1) where F is any of u, v, or cp, Equations (5.10) and (5.11) then become a homogeneous system of ordinary differential equations with the homogeneous boundary conditions v(A) = v(B) = 0^ and c plays the role of an eigenvalueo The values of c for which B' 21 12 22 2, fl (u| + iv] + y Ipj )iy o will be called the spectrum for this system3 and if ci > 0 the motion is un

66 stable. Invoking (3.13) and eliminating u, we obtain two equations for the two remaining unknowns. With the definitions w = U(y) - c, v = vH/ik, q = H( - R2k2)-, = (H - R1 w ), (5.14) 0 0 and with the tildes omitted, we have q(R wp' - p) + = 0, (3.15) Q(R w*' + At) - p = 0, (5.16) where primes denote differentiation with respect to y. These equations are to be solved subject to r(A) = V(B) = O. Alternate formulations in terms of a single unknown are 2 2 (Wqr) t -_ (k + 7 Qt-. i )* = =0 (3.17) H H Rw 0 or (qcp')' -(kq + + q' /Row)cp = 0, (3.18) the latter equation being subject to cp = Rwcp' at y = A and y = B. 0 We note that if I c = 0(1) as R + 0, then, in the limit, p =, and (3.17) o becomes the quasigeostrophic equation w," - (kew + U" - yc)r = 0, (3.19)

67.52 THE NUMERICAL TREATMENT OF THE GULF MODEL The obvious difficulty in applying quasigeostrophic theory to the Gulf Stream problem is that the Rossby number, based on a relative vorticity of 0.4 -~[0 sec (St sommei c loc, cit.) is not particularly small. Not so obvious but equally importLalt is the fact that the slope of the interface must also be smrnall if the fluid is to be nondivergent in the lowest approximation. In the 2 present paper this slope is the quantity R 7y and the conditions for quasigeostrophic theory to be valid are R << 1, R << 1. The first condition is 0 O sa-tisfied marginally, the second not at all, for based on the data in Stommel's bo.ck R Y7 if the 10~C isotherm is taken to be the interface. This can be o seen either by taking as the characteristic depth D the depth of the 10~C isotherm at the midpoint of the stream or else by directly computing the slope of thie isotherm. Substituting representative values in the definitionof R and - o we use K- 2 1-l tak2 2 Zp/p 2C g' g 010 - m/sec and taking V' = 2 m/sec 5 I) -- 800 m R = 0.2 7 - 2.5

68 5 6 y Quasigeostrophic Solution |\\ ^ without Countercurrent 1 I \ /I / I / I/ - Basic State / I/ Velocity Profile / / —---, Quasigeostrophic Solution - / with Countercurrent 6 / / / / 7 6 Figure 12. Basic state velocity profiles of the barotropic model.

69 Tl the numerical study the value of R will be varied in small increments and the values of y will be chosen as y = 2, 2.5, 3, and the extreme case y = 0.5. The choice as to which of the alternate formulations in terms of a single eq:at i:n is settled by the boundary conditions. Since equation (3.17) is subjected t( vanishii-tg boundary conditions r(A) = -V(B) = 0 while equation (5.18) requires cp = R wq)' at y-A and y = B the selection of (5.17) is obvious.'Ve now specify the basic state zonal current. Two cases will be studied. Case i. u(y) = cosiry + A 7/6, B = 5/6 Thlis velocity distribution is plotted in Figure 12. A countercurrent is included. The marginally stable mode corresponding to this profile with R 0 has the solut:i: n cos [2 (y +)j (7.20) o 2 6 wV.Ith 2 0o 2 as can be verified by direct substitution into equation (3.19). n tas e In this case the "south wall" has been moved northward from A = -ed/c-: two y = - 5/6 thereby eliminatingg the countercurrent and reducing the width of channel. It is readily verified that the quasigeostrophic solution corresponcing to this basic state is given by:

70 53 65 5 K = cos 5 Xy, A -, B 6 with k = 0.8 I C = /(t + ) + io+ In programming a second order differential equation for a "shooting method" the most convenient form is the one in which the first derivative term is absent. Since the boundary conditions on equation (3.17) are zero values at the end points this form can be readily obtained by the substitution:,(y) = g(y) @(y) Substituting in (35.17) and equating the coefficient of the first derivative to zero we obtain: g(y) - [H - R2 y )21/2 0 and the canonical form is: "n + S(y) = 0,;(A) = (B) - (3.21) where * S(y): = - P 1 2 a1 y (5.22) Pi Q' 2 QR Ro w (5R25) The appearance of R in the two denominators in the definition of p2 is not as alarming as may seem, as the expression for S(y) can be written as: s(y) = -k2+ ( [R2 k(u-c)2 - l](u-c) - u + 2 (')/(-) Oi O

71 - 3 R2 Q [u + 2R (u-c)u'] + R y Q23 u' - 3R u' - R (u-c)u"] 4 o o o 2 o o and in the limit R - 0 equation (3.19) is recovered. o The numerical scheme consisted of a centered difference approximation for the second derivative in equation (3.21) and a forward integration treating the problem as an initial value one while continuously modifying the eigenvalue c until the boundary condition at y = B is satisfied, subject to the error inherent in integration. This error was determined by the residue at y = B when the known analytical case was integrated numerically. In addition the step size for the integration was determined by minimizing the error in the numerical solution as compared with the known analytical solution (3.20). It was found that the minimum error occurs for 80 points. In order to allow for changes due to the ageostrophic terms, which are not present in the testing procedure, one hundred points were actually used. The problem of finding eigenvalues efficiently from the standpoint of comuting time was discussed in the baroclinic model. The writer was fortunate in this part of the study to participate in the experimental phase of the "time sharing system" conducted by The University of Michigan Computing Center and the IBM staff, whereby computing time was abundantly available. The extra effort in finding first guesses was thereby diminished and the eight point method described in the digital part of the baroclinic study could be utilized from the outset. The digital program was constructed as follows: (1) Selecting a value of 7, the corresponding c was evaluated from ro

72 J3 2/ 2 2 ro 2 (2) The value of k was taken as k = - ror k = 0.8 according to o 2 o whether Case 1 or Case 2 were investigated (with or without countercurrent). (3) A decrement Ak was set to 0.2 and an initial guess for the unstable quasigeostrophic mode c = (C + c + A + i A ci was chosen. Ac and ro - r Aci were varied until fast convergence of the shooting method was obtained. These values were retained in subsequent changes in k, k = k - n'zAk reducing k as far as the method permitted. 0o (4) When the residue of 4 was within the error described above the corresponding c was considered an eigenvalue, still a quasigeostrophic one. (5) At this point the Rossby number was set to a nonzero small value and the search in a square around the quasigeostrophic value of c just found was repeated until a nongeostrophic mode was found. (6) Using the nongeostrophic value last found, the Rossby number was further increased and the search repeated. (7) The sequence of steps was now reiterated by starting with the last quasigeostrophic eigenvalue and changing the wave number to search for the next quasigeostrophic unstable mode. The nongeostrophic modes were then found by calling on the ageostrophic subroutine that varied the Rossby number. The results of the search are shown in Figures 13-16 where the complex phase speed in the nongeostrophic cases were plotted splitting the real and

73 imaginary parts for clarity. Some representative numerical results are shown in Tables IV through VI for reference. 3.4 RESULTS OF THE BAROTROPIC STUDY 3.4.1 Quasigeostrophic Unstable Modes The numerical technique using the shooting method, described in the previous section, was first applied to equation (3.19) to obtain unstable quasigeostrophic modes. The process consisted of selecting a value of the radius of deformation 7 and varying the value of the wave number k in decrement of 0.2 from its initial critical value at neutral stability to as small a value as the numerical integration permitted. The process was repeated for the values of 7 = 3, 2.5, 2 and the extreme case 7 = 0.5. In order to study the effects of the countercurrent the cases 7 = 2, 2.5 were repeated with the southern wall 7 3 at A = 7 displaced northward to the position A = To facilitate a com- 66' ofci parison the results of both computations are plotted in Figure 13. The solid graphs correspond to the case with countercurrent basic state, the dashed ones to the simple basic state jet. Figure 13 reveals that the existance of a countercurrent destabilizes the flow considerably. The difference in the widths of the channel in the two cases is not believed to be of major consequence. Since the Gulf Stream possesses a countercurrent this case will be considered in more detail. It is observed that for each value of y there is a wave number which is most unstable. Considering the realistic case of 7 = 2.5, the most unstable wave number is found near k = 1.8, which corresponds to a wave length of approximately 350 km. and

TABLE IV UNSTABLE MODES FOR THE GULF STREAM MODEL WITH COUNTERCURRENT (A = - 7/6) (a) = 3.o R = R =0.05 o o k r 1 r i 2.52 0.415218 0.026874 0.445218 0.036875 2.32 0.380217 0.049375 0.410217 0.054375 2.12 0.346216 0.063124 0.376216 0.063125 1.92 0.311965 o.o6g 069999 41966 0.067500 1.72 0.276965 0.069999 0.301965 0.062500 1.52 0.240964 0.063749 0.258464 0.056250 1.32 0.202463 0.052499 0.214964 0.045000 1.12 0.161463 0.037499 0.171463 0.032500 (b) 7 = 2.5 R 0 R 0.1 R = 0.15 k 0 0 0 o o o c c c C. r i r 1 2.72 0.530244 0 2.52 o.488944 0.034999 0.506131 -------- 0.574100 o.040937 2.32 o.448994 0.0612499 0.505244 0.067812 o.540400 0.065625 2.12 0.411494 0.079999 o.469306 0.082187 0.503681 0.082969 1.92 0.372744 0.090624 0.428994 0.086250 0.466494 0.090625 1.72 0.372742 0.094375 0.382743 0.085625 0.423368 0.087500 1.52 0.292743 0 089375 0.334930 0.078437 0.373993 0.073125 1.32 0.248992 0.078125 0.280243 0.069375 0.308368 0.055000 1.12 0.203993 0.060624 0.224305 o. o54608 -- - 0.92 0.153993 0.040624 0.169618 0.040624 -------- -

75 TABLE V UNSTABLE MODES FOR THE GULF STREAM MODEL WITH NO COUNTERCURRENT (a) = 2.5 R = 0 R= 0.1 o 0 k C C C C r 1i r 1i 2.51 0.530244 0 -------- -------- 2.26 o.468916 0.020516 -------- -------- 2.01 0.410322 0,305547 0,418759 0,013555 1o76 0,346259 0.034359 0.340009 0.018422 1.51 0.282197 0.033922 0.279697 0.024547 1.26 0.216572 0.026547 0.21,3446 0,021859 (b) ^ = 2 R =0 R = 0.1 o o C C C C r _ 1 r i 2.31 0.561512 0.022499 0.569824 0.013375 2.11 0.512012 0.035600 0.516387 0.024125 1.91 0.458012 0.045750 0.460824 0.035124 1.7.40 1 O.401761 O 403325 0.042875 1,51 o.345261 0.051374 0,344511 0.046375 1.31 0.284760 O.046875 0.285072 0.043750 1.11. 0.225509 0.038249 0.225510. 036375 0.91 0.165883 0.027562 O.165884 0.026793 0.71 o. 111884 0.013319 0.1n1571 0.013312

76 TABLE VI UNSTABLE MODES FOR THE GULF STREAM MODEL WITH COUNTERCURRENT (a) y = 2.0 R = R = 0.1 k 0 _____ c c c c r i r i 2.52 0.570013 0.o04875 0.605013 0.056875 2.32 0.528763 0.074374 0.564763 0.086375 2.12 0.486262 0.099575 0.522262 0.111375 1.92 0.443761 0.116875 0.482762 0.126875 1.72 0.398762 0.124374 0.434762 0 132375 1.52 0.353762 0.124374 0.387512 0.129375 1.32 0.303762 0.114375 0.333762 0.118375 1.12 0.251262 0.096875 0.275262 0.100875 0.92 0.198762 0.073124 0.213762 0.077124 0.12 0.141262 0.045624 0.15326 0.049624 (b) 7 = 0.5 R O R =0.5 o o. c c. c c 2.52 0.794623 0.058750 0.814630 0.053750 2.32 0.747130 0.112500 0.767130 0.107500 2.12 0.699630 0.158750 0.719650 0.153750 1.92 0.652130 0.197450 0.672130 0.197500 1.72 0.597189 0.227500 0.617129 0.229500 1.52 0.537129 0.247499 0.557129 0.252500 1.32 0.469627 0.257499 0.489626 0.267500 1.12 0.394629 0.252500 0.414629 0.262500 0.92 0.314629 0.233500 0.324629 0.248499 0.72 0.234629 0.197499 0.239626 0.212499

.28.24.20 \ k 0 2.5 Velocity Profile with Countercurrent 0 2.3.1 6 ----- Velocity Profile without Countercurrent A 2. 1 0 1.9 ci x. 7.12 -.08 0 2.04 2.5y =2. Y.2.4.6.8 1.0 0.2.4.6.8 1.0 Cr Figure 13. Unstable quasigeostrophic modes of the barotropic Gulf Stream model.

78 an e folding time of just over three days. The quasigeotrophic results obtained by Lipps correspond to a wavelength of 180 km with e folding time of four days. In his case, however, the horizontal scale was taken as L = 31 km compared with L = 100 km in the present study. The speed of the most unstable mode is slightly less than 0.2 m/sec. The speed of all unstable modes decreases with increasing wavelength. The effect of increases in the radius of deformation is to stabilize the motion for both profiles treated. This result is in complete agreement with analytical results obtained by Jacobs (personal communication) where quasigeostrophic divergent flow is investigated. 3.4o2 Ageostrophic Unstable Modes The quasigeostrophic unstable modes discussed in the previous section served as initial guesses in the search for ageostrophic eigenvalues. Selecting a triplet 7, k, and the eigenvalue c corresponding to R = 0, the Rossby number was set to some small value and a search was conducted. When the boundary conditions were satisfied the Rossby number was incremented and the sequence repeated. When the process no longer converged the triplet 7, k, c was changed. The results of the present computation are shown in Figures 14, 15 16 and Tables IV through VI. Figure 14 displays the two extreme cases 7 = 0.5 and 7 = 35 It is noted that when the radius of deformation is very small,; = 0.5) the ageostrophic effects are negligible even for Rossby vales as large as R = 0.5. The other case, 7 = 3, is very sensitive to the Rossby number as can be seen by the lower curves for R = 0 and R = 0.05. The propagation speeds and growth rates for o o small. are much greater than those for large 7.

79 0.8 y-0.5 Ro=0.5 - -y=0.5 Ro0O Cr 0.5 -y3. Ro-0.05 O —.O5 - — y3-. R O0.2 O 0.2 0.1;01111 -yz3. Ro=0.05 \ =r3. Roso I 2 3 k Figure 14. Unstable nongeostrophic modes of the barotropic Gulf Stream model 7 = 0.5, = 35.

~* = A Tapotu wReazgS jTn9 oTdoqIxoBq aq% jo sapour oTqdoq.soaSuou a[Iqe.sufn [ s jnST, ~ Z I \ I\ \ \ \" A" - \ \ \tI ~// \ \ O.o o \A9/-o=VO~___-, I'o. 9/-=V -= \ 0= o _. O, -0 =~0 —-—.'0 9/-= V \9 9/.-=V 0=0 I 0=0 9/-= V0 o9

*^3z = A Tapour ua;cq.S JTn3D ordojIo.Osq aGqU jo sapoul T-qdojq~soDSuou aTqs4sun 91 aGx[2TJ 9~.. I I I I 9/S-~v^ so'o =~y -\ - =0 ~ - \, ~0'0 = o~._.._ ~'~ 1-0o __ I'O 9/9-:V 0=~ — J0=o0d'*0 9i~~~3 /.- = TO. 18r-'OZ

82 Figures 15 and 16 show the results for the cases y = 2 and y = 2.5. Figure 15 shows that for r = 2 and R = 0.1 the ageostrophic influence is to destabilize the flow and speed up the motion of unstable perturbations if the basic state velocity profile has a countercurrent. Figure 16, corresponding to the realistic case y = 2.5 shows a selective ageostrophic effect; the short zonal waves are destabilized while the intermediate and long waves are stabilized, the effect becoming more pronounced as the wave number decreases below the point of cross over. This ~behavior is found for a basic state'with countercurrent as opposed to the stabilization of all wavelengths in the simple jet basic state case. It is pointed out that the no countercurrent case, with appropriate rescaling and inverting of layers, can be made to correspond to an atmospheric problem. The results pertaining to this problem will be quantitatively different but the ageostrophic stabilizing influence will be retained. The effect of the radius of deformation is found to destabilize the flow for decreasing 7. With the exception of the extreme case 7 = 0.5 the ageostrophic modes near the neutral boundary are more unstable than their corresponding quasigeostrophic ones. Considering the propagation speeds of unstable perturbations it is observed that the ageostrophic modes are faster than their quasigeostrophic parents in the case with countercurrent, The propagation speeds in the simple jet case are rather insensitive to the Rossby number. The main purpose of the study, as outlined in Chapter I, was to examine the direction and extent of the ageostrophic effects on the rates of growth of unstable modes. The finding that the short waves are destabilized is perhaps

83 an indication that viscosity effects are missing. The rather strong stabilization of the long waves is in accord with the observation by Stommel that some eddies associated with the current can be identified for long periods of time. In this respect the results of this study justify the stability approach to the Gulf Stream problem.

CHAPTER IV CONCLUSION 4.1 SUMMARY OF RESULTS 4.1.1 The Baroclinic Problem A linearization technique was applied to the equations of hydrodynamics subject to the physical assumptions of a frictionless adiabatic atmosphere in which the flow is incompressible and hydrostatic. The resulting system of scaled equations was linearized about a baroclinic basic state while utilizing the Boussinesq approximation and prescribing a zonal wave character to the perturbation quantities. The elimination of pressure, density, and the horizontal components of the velocity from the linearized system led to an eigenvalue problem formulated in terms of a partial differential equation in the vertical velocity w, subject to vanishing boundary condition at the top and bottom of a finitely deep atmosphere. The stability of the flow is determined by the imaginary part of the complex eigenvalue, which represents the phase speed of the perturbation quantities. In order to study the ageostrophic effects resulting from terms involving the Rossby number of the flow, the partial differential equation governing w was treated by a separation of variables technique resulting in an ordinary w equation. The consequences of this separation are a simplified mathematical problem at the expense of the ability to satisfy vanishing boundary condition on the transversal velocity v. 84

85 The deductions from the ordinary boundary value problem have shown that, if the product of the zonal wave number and the Rossby number are less than unity, the stability boundary consists of the segment 0 < c < I on the real c axis. This result was obtained by treating the singularities of the w equation by the method of Frobenius, while the regular part of the c axis was investigated r by a variational method. A general condition on eigenvalues dependent on the parameters of the problem was also derived. The condition, formulated in terms of an integral inequality, separates the complex c plane into possible domains of eigenvalues versus areas in which eigenvalues cannot be located. Using the information provided by the analytical considerations a differential analyzer was employed to give crude solutions and eigenvalues. These first guesses were introduced into a digital computer program for refinement. The numerical results of the eigenvalue problem show that ageostrophic effects, unlike quasigeostrophic theory with latitudinal variation in the perturbation quantities, are selective. The short meridional waves are associated with destabilized zcnal waves for all zonal wavelengths, while the long meridional disturbances are associated with destabilized long zonal waves and stabilized short ones. The references for stabilization (destabilization) are the quasigeostrophic growth rates. A detailed numerical study of the physical variables has led to zonal vertical cross sections for the vertical velocity, pressure, and temperature fields. The gross features of these fields are in good agreement with general meteorological conditions including the westward tilt of pressure systems and the tem

86 perature pressure phase relations. The unexpected result that all perturbations, regardless of size, are traveling from west to east with the mean speed of the basic current cannot be considered an improvement on quasigeostrophic theory and is probably traced to the infinitely wide plane or the missing beta effect, or both. 4.1o2 The Barotropic Problem A two-layer model with constant density in each layer and no motion in the deep lower one was chosen to represent Gulf Stream conditionso Two basic state zonal currents were studied, one with, the other without a countercurrent. A linearization technique resulted in a characteristic value problem in terms of a single dependent variable. The parameters influencing the flow modes are the radius of deformation 7, the zonal wave number k, the Rossby number and the basic state current distribution. The numerical treatment consisted of determining the modes for various combinations of parameters and currents. The findings indicate that an increase in 7 results in a more stable flow under all conditions. This result was found for both velocity profiles9 in both the quasigeostrophic and nongeostrophic cases, and for all wavelengths considered. In the quasigeostrophic case the stability dependence on the zonal, wavelength is such as to increase the growth rates when the wavelength, is increased above its value on the neutral stability boundary. At a certain wavelength, depending on the value of 7- and R and on the particular basic current chosen, a further increase in wavelength is associated with decreasing growth rates. This result was found for all values of 7, for both

87 basic currents. The presence of a countercurrent was found to destabilize the flow considerably, in both the quasigeostrophic and nongeostrophic cases. As for the propagation speeds they decrease with increasing wavelength in all the cases studied. The ageostrophic effects were found to stabilize the simple jet case. This observation should therefore hold for the corresponding atmospheric problem. For a basic flow resembling that of the Gulf Stream the short waves are found to be more unstable for nonzero Rossby numbers. For a representative value of 7 this effect is reversed below an intermediate value of the wavelength and considerable stabilization takes over. This last remark validates the application of stability techniques to the Gulf Stream problem. 4.2 CONCLUDING REMARKS AND SUGGESTIONS FOR FUTURE WORK 4.2.1 The Baroclinic Problem The present study has successfully shown that the incorporation of latitudinal variation in the perturbation quantities of a linearized baroclinic model atmosphere is strongly coupled with the ageostrophic character of the problem. The observation by Phillips (1964) that quasigeostrophic theory loses its selective qualities with respect to the meridional wave number through the loss of the meridional derivative is confirmed for a flow in an infinitely wide plane. The opposite effects on different transversal wavelengths indicate that the north-south variation in the flow has more than the distortive influence on the stability of the motion as predicted by quasigeostrophic theory.

88 The unresolved problem in the present treatment is the possible influence of lateral walls at two finite latitudes. The natural next step would then be to study the partial differential equation governing w. To this end a numerical search in a rectangle would be appropriate. The eigenvalues obtained in this study could be tried as first guesses for the rectangular region scheme. It is felt that this suggested study would bring about the desired variation in the perturbation speeds as a -function of the physical parameters. Once the effects of the vanishing north-south velocity at two walls are established, the effects of the variation of the Coriolis parameter with latitude should be examined. This problem will result in a partial differential equation with coefficients which are functions of both y and z. The separation technique employed in this study is not expected to apply but a numerical study on a rectangular domain should provide insight into the mechanism of the beta effect. Once again the eigenvalues corresponding to the previous finite channel case could provide good initial guesses for the present suggested problem. 4.2,2 The Barotropic Problem The major result of the study is the evidence that ageostrophic effects operate in such a way as to bring some computed unstable growth rates into closer agreement with observation. The conclusion that stability methods are indeed valid for Gulf Stream theory should initiate a more general formulation of the problems in particular a stratified fluid model. In addition the results obtained so far can be related, under the appropriate interpretation,

89 to atmospheric conditions. To this end the relative positions of the layers has to be reversed, the upper motionless layer representing an approximation to the stratosphere. The qualitative results are expected to be of the same nature as in the Gulf model, namely, a stabilizing ageostrophic modification. The actual values of the pertinent quantities would obviously be altered through. the different scaling involved. The problem arising from the destabilization of the short waves should be investigated by incorporating viscosity in the two layer model. Since numerical methods are likely to be applied the complexity of this new problem should be comparable to the present problem.

APPENDIX A.1 TRANSFORMATION OF EQUATION (2.9) INTO CANONICAL FORM Equation (2.9) can be written as: 2ilR (z-c)-2 -2iQR - B(k 2+ )(z-c) A 0 A 0 A w + -----— 2 —) w + 2 w 2 o (A.1) (z-c)[l-k2R2( -c) z (z-c)[1-k R (z-c) ] O 0 with boundary conditions w(o) = w(l) = o, Let w(z) = g(z) w(z) (A.2) Substituting (A.2) in (A.1) and equating the coefficient of w to zero we have: z 2g + Pg = o (A.3) z where 2iQR (z-c)-2 0 ~P =~22 N (z-c ) lk2R2 (z_ c )2] The equation for w is: gzz gz w -+ - P- + Q)w = o, w(o) = w() = o (A.4) zz g g where 2iQR + B(k2+ 2)(z-c) q, = -— 2-o-.-. — - 00 (z-c)m -k(A ) we hz From (A.3) we have - JPdz g(z) = e 90

91 and a straightforward integration leads to Pi [1 - kR (z-c)] g(z) = (z-) (A.5) [1 + kR (z-c)] 2 o0 where 1 i p 1 iQ P= +- p =-+1 2 2k 2 2 2k The coefficient of w in equation (A.4) can be calculated using the explicit form of g(z). A faster and simpler process is to use equation (Ao3) in order to obtain g and g. Using equation (A.5) it is found that z zz gzz g P p2 S(z) = + p - + = - (A.6) g g 2 4 substituting for P and Q their appropriate expressions we finally have: b2 2 2 2 2(2 2 2 B(k + ) + 2k R R (k + S(z) =-. - + p 2 (z )2 1 - k2R (z-c) [1 - k R2(z-c), o 0 A.2 DETERMINATION OF THE EXPONENTS FROM THE IIIDtCIAL EQUATION S(z) can be written as: S(z) = +- + (A.7) 22 R k 1 (z-c) +o z-(c- -) -(c+- ] (~-~ k z-(CR k o 0 A o wheek z-(c- R z-(c+

92 A =A - k2 B(k2 2 A ) R2 (k2+ ) A 2 = A2 - k R - 2 A2 = R ( 4 10 2 o 2 2 o 4 The poles of S(z) are clearly at: 1 1 Z - C, z = C + - Z= -— C - kR' kR 0 0 and the coefficients of the Laurant expansion, using standard notation, are: At z = c: a = o b = -2. -1 -2 1 2 At z = c +:a = b2 1 2 At z = c - Rk a1 =0, b -2 o R k 0 Substituting in the indicial equation n(n-l) + a -n + b = o, (A.8) -1 -2 the exponents relative to the singular points are found to be: At z = c: n = 2, - 1 1 2 2 1 i At z = c + -: n -n +, n = (1 +-). - R k R2k2 2 - k o 0 A.3 QUASIGEOSTROPHIC SOLUTION OF THE BAROCLINIC PROBLEM Equation (2.11) with R = o reads o (+ )w = o w(o) = w() = o (A.9) (z-c)2

93 where 7y = B(k 2+ ) and g(z) reduces to g(z) = z-c. Two independent solutions of (A.9) are: w = e _)) W e7Z (+ 1 ) (A.10) ~1 *Y(z-c) 2, (Z-c)) as can be verified by direct substitution and a Wronskian test. The general solution, in terms of two constants to be determined from the boundary conditions, is w = A1W + A W 11 22 Application of the boundary conditions leads to 1 1 At z =o: A (l - )) + A l + ) 1 Y(-c) 2- 7(-c) (A.11) At z = 1: A (1- + A ( + )e =. The condition for nontrivial A A2 is the vanishing of the coefficients' determinant. This condition leads to the guasigeostrophic eigenvalue: 1 1 1 1/2 (A.12) c = + - + - - - cothy] (A.12) 7 Solving for A (in terms of A of course) we have A y- A 2 l-yc 1 The solution for w, up to a multiplicative constant, can now be written, recalling the factor g(z):

94 A 1 1 -7z w = (l-yc)(z-c - )e + (1 + c)(z-c + -)e7 (A.13) 7 7 where c is given by (A.12). A.4 CALCULATION OF NONGEOSTROPHIC PRESSURE AND TEMPERATURE The system (2, lb)...(2, 5b) can be written: ikR (z-c) -1 ik 0 u -R2w o 0 1 ikR (z-c) ii 0 v o (A. 14) ik i 0 0 p -R w o z L 0 1 0 ik(z-c) p Bw where the hydrostatic relation has been omittted and should be verified when (A.14) is solved. Let ikR (z-c) -1 ik 0 o 1 ikR (z-c) i O0 / = o~t = (A.15) ik iQ 0 0 1 O ik(z-c) = -R k2(k2+2 )(z-c)2 ikR (z-c) -1 ik -R w o o 1 ikR (z-c) i O0 A = (A.16) ik iQ 0 -R w o z 0 1 0 Bw = iR k(z-c)[k R (z-c)2 - l]w - [iQR +k 2R(z-c)]w} 0 0 Z 0 0

95 ikR (z-c) -1 ik -R w o o 1 ikR (z-c) iQ 0 A = (A.17) ik iQ 0 -R w o z 0 1 0 Bw = iR k[iR + B(k2+2 )(z-c)]w + [1 - iQR (z-c)]w 0 0 0 Z then A i i~l - 2 2]A 2 2 p = P = 2 {[1 - k R (z-c) ]w +liR + k R (z-c)]w} (A.18) A 2 2 o z o o k(k +i )(z-c) P i 22 -P = - = — { [[iQR + B(k2+Q )(z-c)J- +[1 - ilR (z-c)]w j p 22 2 o o z P k(k2+ 2)(z-c)2 0 (A.l9) and the relation P +p = o can be verified to hold. The temperature, in the z Boussinesq approximation, is given by: p = p (1 - a(T - T ). (A.20) o o Since we are interested in normalized fields and since the density and temperature are linearly related in this model, the density can be used to represent the temperature distribution. A.5 QUASIGEOSTROPHIC PRESSURE AND TEMPERATURE Let R = o in equations (A.18) and (A.19). Applying the remark at the end o of the last section we have i A p 2 2 w (A,21) k(k +Q )(z-c)

96 1T = -,, i,- [W + B(k2+2)(z-c)W] (A.22) k(k2++2)(z-c) z A A Substituting for IT and w their values as determined by equation (A.13) of the appendix, the system (2.27) is obtained.

BIBLIOGRAPHY Arnason, G., et. a. 1967. A case study of the validity of finite difference approximnat.ions in solving dynamic stability problems. J. Atmos Sciences 24, p. 10 17. Case K.M. 1960. Stability of inviscid plane couette flow. Phys. Fluids. 3, p. 145. Charn:ey J.G. 194-7. The dynamics of long waves in a baroclinic westerly current. J. of Met 4, No. 5. Eady, ET. 1949. Long waves and cyclone waves. Tellus 1 No. 3. Gulf Stream Summaries, U.S. Naval Oceanographic office. Hloward, L.N. 1963. Neutral curves and stability boundaries in stratified flow, J. Fluid Mech. 16, p. 3335. Ifce, E. L Ordinary Differential Equations, Dover publication. Jacobs, S.J. 1967. An asymptotic solution of the tidal equations. J Fluidc Mech. 20, p. 417. Jacobs, S.J. and A. Wiin-Nielsen 1966. J. Atmos Sci. 23 pp. 682-687. Kno, H. L. 1949. Dynamic instability of a two-dimensional nondivergent flow ir1 a baiLotropic atmosphere. J. Meteor, 6 105-1.22. Lipps, FoB. 19635 Stability of jets in a divergent barotropic fluid. J. Atrrn:. Sci. 20, p. 120. Mc:lintyre, M.E. 1965. A separable nongecstrophic baroclinic stability problerm. J. Atmos. Sci. 22, p. 730. Miles, J.W 1961.. On the stability of heterogeneous flows. J. Fluid Mech. 10 p. 496.:Pedlc sky, J. 1964. The stability of currents in the atmosphere and the ocean. Part I. J. Atmos, Sci. 21 p 201 Pecdiosky, J. 1965. The stability of currents in the atmosphere and the ocean. Part II, J. Atmos, Sci 21, p. 3L42. 97

98 Phillips, N.A. 1964. An overlooked aspect of the barotropic stability problem. Tellus 16, 2, p. 268. Phillips, N.A. 1967 Thermal convection colloquium, NCAR, Technical notes TN 24, p. 32. Rossby, C.G. 1939. Relation between variation in the intensity of the zonal circulation of the atmosphere and the displacement of semi-permanent centers of action. J. Mar. Res., 2, 38-55. Starr, V.P. 1954. Studies of the atmospheric general circulation. Final report, Part I, General Circulation project contract No. AF-19-122-53, Dep. of Meteor, M.I.T., p. 535. Stern, M. 1961. The stability of thermoclinic jets. Tellus, 13, p. 503. Stommel, H. 1965. The gulf stream. University of California press. Wiin-Nielsen, A. 1961. On short- and long-term variations in quasibarotropic flow. Mon. Weather Rev. 89, p. 461.