On the propagation of waves in the atmosphere

The leading-order equations governing the unsteady dynamics of large-scale atmospheric motions are derived, via a systematic asymptotic approach based on the thin-shell approximation applied to the ellipsoidal model of the Earth’s geoid. We present some solutions of this single set of equations that capture properties of specific atmospheric flows, using field data to choose models for the heat sources that drive the motion. In particular, we describe standing-waves solutions, waves propagating towards the Equator, equatorially trapped waves and we discuss the African Easterly Jet/Waves. This work aims to show the benefits of a systematic analysis based on the governing equations of fluid dynamics.


Introduction
Waves play a fundamental role in the development and evolution of our atmosphere; it is generally accepted, for example, that Rossby waves are the most important large-scale waves (controlling the weather in the midand higher latitudes), although gravity (buoyancy) waves and Kelvin waves are also significant. Of course, there are many other ingredients that contribute to the motion of the atmosphere, both locally and globally; here, we focus on the more familiar and well-documented wave-like motions. The complexities of the unsteady flows, including wave motions, suggest-apparentlythat these cannot be systematically extracted from the full set of governing equations. Indeed, the familiar way forward is to simplify, to extreme levels, both the geometry and the assumed background state (which is  is the eccentricity. The coordinate system (ϕ, β, z ), z being the vertical distance up from the surface of the ellipsoid, is associated with the ellipsoid, which is rotating about its polar axis with (constant) angular speed Ω ≈ 7.29 × 10 −5 rad s −1 . In this system, the unit tangent vectors at the surface of the ellipsoid are (e ϕ , e β , e z ); e ϕ points from West to East along the geodetic parallel, e β from South to North along the geodetic meridian and e z points upwards (figure 1).
The system is valid throughout the space, except along the direction of the polar axis, where this description fails because e ϕ and e β are not well-defined at the two Poles. Note that the direction of apparent gravity (seen by an observer in the rotating frame) is normal to the surface of the ellipsoid, which is the isosurface of the geopotential that defines zero elevation; the geopotential Away from the polar axis, we represent a point P in the atmosphere using the hybrid spherical-geopotential rotating coordinate system (ϕ, θ , z ), obtained from the spherical system (e ϕ , e θ , e r ) and the geopotential system (e ϕ , e β , e z ). Here, ϕ and θ are the longitude and geocentric latitude of P, respectively, β is the geodetic latitude of the projection P * of P on the ellipsoidal geoid and e z points upwards along the normal P * P to the geoid (which intersects the equatorial plane in the point P e ). The unit vectors (e θ , e r ) are obtained by rotating the unit vectors (e β , e z ) by the angle (β − θ ), in the plane of fixed longitude ϕ. (Online version in colour.) being the sum of the Newtonian (gravitational) potential and the centrifugal potential due to the Earth's diurnal rotation. The details of this geometric model are fully described in [3]. In summary, we transform from spherical coordinates, (ϕ, θ, r ), to the hybrid spherical-geopotential coordinates (ϕ, θ, z ), and couple this with a transformation of the velocity vector and the gravity term. Thus the components of these vectors, which are defined normal and tangential to the surface of the ellipsoid, are decomposed into components in the spherical system. As pointed out in [3], the advantage of using this spherical-geopotential hybrid rotating coordinate system rather than the conventional spherical potential approximation (described in [12]) is that the formulation retains the details of the curved-space geometry of the Earth and the leading-order (geometrical) correction terms apply to the background state of the atmosphere but do not interact (in the leading-order perturbation) with the dynamics of the atmosphere. The relation between r and z is given by the identity with z = εz (see below). Furthermore, we must use a non-dimensionalization that enables us to put the equations in a form that is relevant for a discussion of atmospheric flows. To do this, we introduce a scale length, which we take to be the maximum height of the troposphere, H (which is about 16 km at the Equator), and an associated speed scale Ω H (approx. 1.2 m s −1 ). Thus we define in spherical coordinates is (u , v , w ) and (u, v, w) is a non-dimensional velocity with w normal, and (u, v) tangential, to the ellipsoidal geoid, where β − θ = e 2 sin θ cos θ + e 4 sin θ cos 3 θ − εe 2 z sin θ cos θ + O(e 6 , εe 4 , ε 2 e 2 , ε 3 ).
The constant k measures the size of the velocity component normal to the ellipsoid, and its choice controls the type of problem to be examined. In the case of thin-shell theory (i.e. small ε, as we have here), the choice consistent with the equations (and kinematic boundary conditions, for example) is k = ε; we use this identification hereafter.
The pressure (p ) and the temperature (T ) are correspondingly non-dimensionalized: where R ≈ 287 m 2 s −2 K −1 is the universal gas constant and (Ω d E ) 2 /R ≈ 800 • K. The dynamic eddy viscosity is assumed to vary only in the radial direction (this being the choice that was made in [3]); the transformation from r to z then shows that, for the viscous terms that appear in the Navier-Stokes equation, we have We are now in a position to quote the governing equations, suitably approximated for small ε and small δ = e 2 , this latter parameter measuring the effects of small deviations of the ellipsoid from the spherical. The equations follow directly from those derived in [3]; here, though, we retain the time dependence; the errors in each equation are recorded and we write out each component of the Navier-Stokes equation. We obtain The additional parameters that we have introduced here, based on the (constant) average acceleration of gravity at the surface of the Earth (g ≈ 9.8 m s −2 ), the specific heat of air (c p ≈ 1500 m 2 s −2 K −1 ) and the thermal diffusivity of predominantly dry air (ν ≈ 2 × 10 −5 m 2 s −1 ), are all O(1), i.e. held fixed as the limiting process is performed. Note that the appearance of large Reynolds number, R e , and small κ, suggests that additional approximations could be performed, but these are unnecessary (and would then require the introduction of viscous and thermal boundary layers, for example). We note that equation (2.6) is the first law of thermodynamics, with Q representing the (non-dimensional) heat sources (or sinks), expressing the change of total energy due to any heat exchanges. The second law of thermodynamics, setting limits for the transformations between heat energy and the sum of kinetic and potential energies (see [13,14]), will be not be of direct concern here: in the general discussion we do not, nor do we need to, specify the sources of the heat input since we do not address the issue of the climate being in a non-equilibrium thermodynamical state (see [15][16][17] for the challenges encountered in trying to exploit the second law to determine a lower bound to the entropy production). Furthermore, we choose to apply our results to regions of the Earth where the topography is relatively unchanging; the inclusion of significant orography, which will affect the bottom boundary condition and the local generation and structure of waves, is a refinement relegated to a future investigation.

Asymptotic structure of the solution
We seek a solution of (2.1)-(2.6) as an asymptotic expansion which begins where q (and correspondingly q n and q n ) represent each of the variables u, v, w, p, ρ, T. Furthermore, we assume that the boundary and initial conditions are consistent with these expansions, and that the heat-source term also follows this pattern. (The lack of any evident singularities in the velocity, temperature or pressure fields, providing that the neighbourhood of the Poles is avoided, and in the small-eccentricity corrections, indicates that this asymptotic structure is otherwise uniformly valid, i.e. higher-order terms remain small; more details are given in [3].) The leading-order problem (which is time-independent, being driven by ∂ρ 0 /∂t = 0 from (2.4)), then gives note that the third equation above appears twice in the set (2.1)-(2.6). There is a solution of this system (see [3]) in which the temperature decreases linearly with height, independently of the velocity field (which nevertheless does appear in the equations at this order): the classical adiabatic model of the lower atmosphere. This corresponds to an external heat source of zero (Q 0 ≡ 0): the heating is supplied by heat transfer from the surface of the Earth upwards into the atmosphere. (It should be noted, in our approach, that we determine the heat sources that are consistent with the motions that we describe; clearly such heat sources then drive the well-behaved flow properties. There is therefore no need to add general, auxiliary conditions on the background heating that would ensure the existence of solutions.) The resulting solution is We impose the boundary condition the constant B is fixed by knowing the pressure (or density) at ground level at some θ. We have chosen to work with a background state which is steady; other choices are possible, e.g. via the velocity field (u 0 , v 0 , w 0 ) or by allowing variations on a suitable, long time scale. As formulated here, we may choose to use the boundary conditions appropriate to the region of the Earth centred on (z = 0, θ = θ 0 ) and to a particular season. This solution describes the background state, at leading order in ε and δ, for the troposphere, and it is only this section of the atmosphere that we shall discuss here; a related description for higher levels in the atmosphere can be found in [3].
The dominant effects of the ellipsoidal-geometry correction, as these distort the otherwise spherical solution, arise at O(δ); the equations at this order are these equations have been simplified by using the equations that describe the leading order. The solution of this system simply produces a (uniformly valid) small adjustment to the background state, as described in [3]. However, the solution for these temperature, pressure and density perturbations is time-independent, and so the first law of thermodynamics (at this order) shows that, in order to maintain this structure, a heat source that moves with the fluid is required, e.g. some latent-heat input. We do not pursue the details of this contribution here because it uncouples from the leading-order time-dependent element of the motion (e.g. waves) that we want to investigate. The important time-dependence appears at O(ε), which is also where the dominant contribution to both the dynamics and thermodynamics is evident; this general structure, for steady flow, is discussed in [3]. The resulting equations are and with Q 0 ≡ 0 for the troposphere. For time-independent flows, equations (3.1)-(3.6) recover precisely those obtained in [3], but with (3.6) replaced by a more complete version of the first law of thermodynamics, which provides a mechanism for describing the heat sources associated with the motion. This difference arises because the leading terms in the time derivative vanish, so we must, perforce, go to the next order for the description of a steady dynamic-thermodynamic balance. The considerations in [3] describe the combined dynamic and thermodynamic properties of the steady troposphere; this further analysis has now demonstrated that, on time scales measured by 1/Ω , time dependence appears at this same order. We can expect, therefore, to find unsteady solutions which may not be wave-like, but also (superimposed on any such motions) wave solutions, both being the same size (in the ε sense) and both appearing in the dominant dynamic-thermodynamic structure of the atmosphere and, significantly, not as a small perturbation to a constant state. Our main aim-to put the theoretical study of the motions of our idealized model of the atmosphere on a solid foundation based on the general equations of fluid mechanicsis motivated, in part, by the reported observations and data from field studies (e.g. satellite observations, regional and global weather data, together with weather forecasts, climate modelling and reanalyses). Indeed, atmospheric waves-we are interested here in stable wavesare typically identified experimentally from satellite observations, after having subtracted the background flow field via appropriate filtering methods. Most studies use spectral analysis to isolate various wave modes, seeking agreement with structures revealed in theoretical studies; see the discussion in [18]. The system (3.1)-(3.6) shows that a careful asymptotic derivation can accommodate a dynamic-thermodynamic balance, which encompasses more behaviour (and is readily accessible) than simply perturbing a constant zonal flow (as described in [1,19,20]). The latter oversimplification of the dynamics was driven, presumably, by the need to use a very restrictive form of heat forcing, coupled to the limitations of a simplistic model of the background state.

Development of the unsteady problem at leading order
The main aim in this initial investigation is to examine the unsteady problem associated with the dominant dynamic-thermodynamic balance, which arises at O(ε) in our idealized model of the atmosphere. (This model excludes, therefore, the role of turbulence, instabilities, etc.; it is planned to extend our approach by adding more extensive attributes of the atmosphere in subsequent studies.) In particular, we seek solutions-hopefully in a manageable form-of two types: unsteady, but not wave-like, and wave solutions. However, before we proceed, it is expedient to make some small adjustments to the presentation of equations (3.1)-(3.6); this will ease the later analysis.
First, writing equation (3.6) can be integrated to give where A 1 (ϕ, θ, z) is an arbitrary function, fixed by the initial data on the perturbation temperature and pressure. Now from equation (3.3), we see that and B 1 (ϕ, θ, t) is an arbitrary function determined, for example, by the perturbation pressure on the ground. The expression for Π 1 measures the total heat input, weighted with respect to the background temperature, from the bottom of the atmosphere upwards, and over time. The complete description of the thermodynamic properties of the atmosphere, at this order, is then obtained from (3.5) and (3.6) in the form where p 1 is given by (4.2).
It is now convenient to introduce and first we show how the heat forcing (4.5) can be related directly to a familiar descriptor used in atmospheric studies. To do this, we introduce the potential temperature, T = Tp −1/c p , and the square of the (non-dimensional) Brunt-Väisälä frequency, N, where see [19]. (We have introduced Ω d E /H for the non-dimensionalization leading to N.) This is used as a measure of the static stability of the atmosphere, the adiabatic perturbation of a fluid parcel about its equilibrium position being governed by for small vertical displacements ζ . Considering only the unsteady perturbation, i.e. the O(ε) terms (because the O(δ) terms uncouple from any wave-like motion and contribute only to the background state), we have which gives and then invoking (4.3)-(4.4) we obtain, at this order Therefore, we may identify It is usual to assume that N 2 is a constant (see [19]), in which case the general explicit solution of equation (4.7) follows directly, the necessary and sufficient condition for stability being N 2 > 0.
The assumption of constant N 2 is clearly an oversimplification. While an explicit, general formula for the solution of equation (4.7), for N 2 varying in the vertical direction, is not available, the comparison method shows that stability holds if N 2 > 0 (and N = 0 implies neutral stability, while instability corresponds to N 2 < 0); see [21]. and respectively, together with from equation (3.4). Moreover, because the background state is a function of only, it is useful to transform from (ϕ, θ , z, t) coordinates to (ϕ, θ, ζ , t); equations (4.9)-(4.11) then become . This completes the development of the system that we plan to discuss here, namely, equations (4.12)-(4.14), in conjunction with (4.2) and (4.4).

Solutions of the time-dependent problem
We now turn to the main thrust of this presentation: the construction of some explicit solutions, all derived from our overarching set of general equations that arise at O(ε). In particular, we examine certain types of wave-like motion; among the plethora of unsteady motions, this is a reasonable testing ground for the applicability of our equations. Although we strongly advocate an investigation, using numerical techniques, of the complete system (4.2), (4.4) and (4.12)-(4.14), for realistic forcing functions and flow properties, e.g. suitable, variable eddy viscosity, our aim here is to produce analytical results. This will, we believe, both confirm and explain some familiar theories (but now placed in a wider context, both in terms of minimal simplifying assumptions and precise approximations and errors) and also generate new observations about the largescale structures of the moving atmosphere. There are clearly two main avenues to explore: either assume viscous flow or take the inviscid limit; the former produces closed-form solutions only in special cases (e.g. M = constant), but the latter offers many possibilities. Indeed, because of the very large Reynolds number (see §2), it might be argued that, in any event, the role of viscosity can be ignored. So we look, firstly, at some inviscid flows. A comment about viscous flows will be presented in §7.
The inviscid limit of equations (4.12)-(4.14) can be written as  (5.4) and so the explicit dependence on the background state is suppressed. This system provides a complete description of the (inviscid) velocity field at this order, for a given forcing F; this F can, in turn, be used to produce the perturbation to the thermodynamic state (given by p 1 , ρ 1 , T 1 ), together with the identification of any required heat sources. It is of some significance that equations (5.1) and (5.2) take a very familiar form (see [19]), but the forcing here drives the leading-order velocity field, not a perturbation of it about some uniform state. On the other hand, if we assume that W = 0 (which is not a consistent choice within the thin-shell approximation), then equation (5.3), combined with the other two, produces a necessary constraint on F in order for such a solution to exist: so, we must expect that we can allow only very special heat sources to drive the atmospheric flow, in this case. We are now in a position to explore, in some detail, what equations (5.1)-(5. 3) tell us about unsteady motions in the atmosphere: we outline the type of results that can be obtained by presenting a few examples, together with some general observations and a case study.

(a) Solution harmonic in time
The simplest solution to seek is one in which all the variables are proportional to e −iωt , for some real constant ω; so we write where each of the 'hatted' functions depends on (ϕ, θ, ζ ). From equations (5.1) and (5.2), we find directly that It should be noted that, if this horizontal velocity field is to be defined for − π 2 < θ < π 2 (so avoiding the Poles; see below), we must have |ω| > 2. We now use (5.5) and (5.6) in (5.3), and so obtain the equation relating W and F: Since W = 0 at ground level, from equation (5.7) we may determine the vertical velocity component, w 0 = W/ρ 0 , given the forcing function F (which, in turn, produces the horizontal velocity components and the perturbation to the thermodynamic state, and can be used to identify the associated heat sources). Note that the neighbourhood of the Poles must be avoided, because this is where W is undefined, according to equation ( We can investigate a little further: equation (5.7) admits solutions for F harmonic in ϕ and ζ (but not in θ ); we set where k and l are constants (and we are interested in real k, so that we have a component of propagation in the azimuthal direction, but l may be complex-valued); we do not impose the same structure on W since it would amount to W ≡ 0, given that W = 0 at ground level. The resulting equation relatingF and W is

(i) Reappraisal of the classical dispersion relation
While the dependence on θ in equation (5.8) makes any general further development challenging, we can be confident that we have the correct representation of the behaviour in the meridional direction. It is therefore this equation that should be examined if we require a consistent representation in the neighbourhood of particular θ; this we now do by seeking a harmonic solution close to θ = θ 0 . To accomplish this, we evaluate all the coefficients in equation (5.8) on θ = θ 0 , and then seek a zonally harmonic solution: where m is a constant, not necessarily real. The resulting 'dispersion' relation becomes F 0 2kl sin 2 θ 0 cos θ 0 + 8ωm sin θ 0 cos 2 θ 0 ω 2 − 4 sin 2 θ 0 − ωm sin θ 0 − ωl(ω 2 − 4 sin 2 θ 0 ) cos θ 0 + i ωk 2 cos θ 0 + ωm 2 cos θ 0 + ωml sin θ 0 cos 2 θ 0 − 2k cos θ 0 − 16k sin 2 θ 0 cos θ 0 ω 2 − 4 sin 2 θ 0 We should note that, if we seek such a solution, harmonic in (ϕ, ζ , t) and locally harmonic in θ , which is then evaluated for zero heat input, i.e. F 0 = 0, then there is no motion at this order. Equation (5.9) is not a dispersion relation in any conventional sense, unless some very special assumptions are made about ∂ W/∂ζ /F 0 ; rather, this equation determines W, given F 0 and all the wavenumbers and the frequency. However, our approach enables us to explain how a familiar result can be recovered. In equation (5.9), we ignore all the terms that are generated by the existence of a background state or that arise from the spherical geometry other than the term involving ∂ W/∂ζ ; in addition, we assume that there is no contribution from a heat source, which then leaves us with F 0 ωk 2 cos θ 0 + ωm 2 cos θ 0 − ig cos θ 0 (ω 2 − 4 sin 2 θ 0 ) ∂ W ∂ζ = 0.
Even though θ 0 may take any value, we must insist that cos θ 0 = 1 in this context; finally, we also assume that there is an auxiliary condition (which is not available at this order), which gives and so ω{l 2 ω 2 − [(k 2 + m 2 )N 2 ] + 4l 2 sin 2 θ 0 ]} = 0, the standard result involving the Brunt-Väisälä frequency, N; see, for example, [19]. In other words, the familiar result has no bearing on, or relevance to, the movement of the atmosphere: it cannot be recovered in any systematic analysis of the governing equations. (The errors arise as follows: (i) it is inconsistent to replace cos θ by 1 in the evaluation of a term such as (∂/∂θ){cos θ (∂/∂θ)(v cos θ )}, (ii) there is no background state and (iii) there is no equation available, at this order, to allow the introduction of N 2 .) We now return to the general time-harmonic version of our equations, (5.5)-(5.7), and seek solutions that are exponential in ϕ and ζ : where α and γ are constants (which may be chosen to be real or complex, as appropriate, in specific cases). Thus equations (5.5) and (5.6) become respectively; correspondingly, equation (5.7), after an integration in ζ and imposing the boundary condition W = 0 at the bottom of the atmosphere (i.e. on ζ = − 1 2 cos 2 θ), becomes iω cos θ df dθ + 2αf sin θ e γ ζ − e − 1 2 γ cos 2 θ e αϕ . (5.13)
This solution represents a standing wave, i.e. oscillating in time but not propagating. We note that, even though the expressions derived above exhibit no variation in the vertical direction, (5.4) shows that there is a vertical structure, by virtue of ρ 0 , in the horizontal velocity field (u 0 , v 0 ).

(iii) Modulated travelling wave
Waves propagating towards the Equator (for example, as shown in figure 2) can be obtained from our system (  where A is an arbitrary (complex) constant and λ some real constant. Then, from (5.11) and (5.12), and reinstating the harmonic time-dependence, we obtain a modulated wave motion described by U = 2Ae i(λθ−ωt)+γ ζ sin 2 θ and V = −iωAe i(λθ−ωt)+γ ζ sin θ, where γ may be chosen to be complex-valued (to produce an exponential behaviour with height, together with a tilt in the direction of propagation). Indeed, we could simplify further by imposing the condition |γ | → 0, which produces the simpler solution This solution, we observe, describes W varying linearly with height, but our development makes clear that many other choices are available, for which suitable heat sources can be identified. The above analysis demonstrates that we can reasonably expect to match what is observed. Of course, the identification of flow structures and processes from the various available satellite data is rather subjective [22], but such exercises nevertheless provide major sources of insight. In the case of the wave pattern visible in figure 2, it is observed, during May, that the surface temperature of the north Indian Ocean becomes the highest among the world's ocean surface temperatures, often in excess of 28 • C over the Arabian Sea, where evaporation largely exceeds precipitation; see [23]. The sea-surface temperature is comparable with the daily maximal land temperature in the coastal areas adjacent to the Arabian Sea, with a significant drop in temperature at night over the land. These changes in the temperature over space and time can give rise to wave patterns   These two examples demonstrate that our single formulation can capture different modes and types of wave propagation in the atmosphere, avoiding the need for ad hoc modelling.

(b) Equatorially trapped waves
Due to the nature of their energy sources, as well as the small Coriolis forcing, large-scale equatorial atmospheric flows present specific structural features. Synoptic-scale disturbances of the atmosphere outside the tropics are mainly driven by horizontal temperature gradients, the primary energy source being the potential energy associated with the latitudinal temperature gradient [19]. In the tropics, however, horizontal temperature gradients are very small, and the energy harnessed by atmospheric flows comes from diabatic heating due to latent heat release, mostly occurring in association with moist convection. The cloudiest regions on Earth (with the oceans significantly cloudier than land) are the tropics, and the higher temperatures in these regions enhance water evaporation into the air, so that the atmosphere becomes very humid. The rapid drop of temperature and pressure with altitude causes some of the water vapour in the air to condense, forming clouds (on average 19 • K colder than the ground) and, if enough water condenses, the cloud droplets can become large enough to fall as rain. There is a considerable release of heat during precipitations, which alters the thermal equilibrium and triggers an adjustment process by means of fast propagating waves (with speeds of the order of 20 m s −1 ; see [26]). This local reaction of the atmospheric flow may also induce remote responses through the excitation of equatorial waves. In particular, equatorially trapped waves are related to the heat forcing that occurs over the Indonesian 'maritime continent', an area straddling the Equator in the middle of the warmest body of water (between the Indian and western Pacific warm pools), whose land-sea heat-capacity differences (due to the mixture of ocean and islands with long coastlines) generate the world's largest rainfall on diurnal cycles, with enormous energy release by convective condensation. With all this in mind, we investigate how our system of equations might model some aspects of these phenomena.
We start by imposing a vanishing meridional velocity component (a reasonable assumption in equatorial regions), and then the system (5.1)-(5.3) reduces to 2U sin θ = − ∂F ∂θ (5.15) and We now seek a solution generated by a heat forcingF(η, θ, ζ ), where η = ϕ − ct, which is azimuthally localized (that is, it vanishes for |η − η 0 | ≥ α, for suitable η 0 and some α > 0, which is taken to be fairly small); furthermore, this function is to be symmetric about the Equator and describes a heat source that moves in the zonal direction at constant speed, c. From (5.15), we obtain (η, θ, ζ ), (5.17) and the meridional symmetry ofF ensures (by l'Hospital's rule) that U is not singular along θ = 0. Combining (5.14) with (5.17) yields the differential equation The fact that W = 0 at ground level ζ = −(1/2) cos 2 θ rules out the possibility that ∂W/∂ζ < 0 somewhere at ground level. Since (5.20) ensures the existence of regions where A > 0, we conclude from (5.22) that A (and therefore the heat forcing) must feature variations in the vertical direction, with ∂A/∂ζ < 0 at ground level. Field data (see fig. 12 in [27]) confirms that the heat forcing typically decreases in the lower and in the upper tropical troposphere, and it increases with respect to height at mid-altitudes (roughly between 2 and 8 km). The classical theoretical approach for describing equatorially trapped waves (see [19]) concentrates on their horizontal structure, using a shallow water model (for a fluid system of mean depth h e in a motionless basic state) and the equatorial β-plane approximation. While the predictions made by relying on these approximations (non-dispersive, eastward propagating trapped waves) are somewhat similar to ours, there are significant differences. Firstly, approximating sin θ and cos θ by the first two terms in their Taylor expansions near θ = 0 leads to an exponential meridional decay rate; ours is not so severe because we have not approximated the θ -dependence. Secondly, the wave speed is taken to be c = gh e , i.e. that for shallow water gravity waves, although an adequate choice of h e is somewhat elusive. In contrast to this, since F = ρ 0 F 1 , we deduce from (4.5) that the value of F at ground level is determined by that of the perturbation pressure, which therefore, due to (5.19), determines the speed c. The third aspect is more fundamental: the classical approach ignores variations of the heat forcing in the vertical direction, while the above discussion of the implications of (5.22) shows that this is an essential aspect of the dynamics. In this context, a simple example of type (5.19)-(5.20) that captures the main observed features is (5.23) where, in analogy to the considerations discussed in [28], and where H(ζ ) is the restriction of a cubic polynomial in ζ to the interval [−(1/2), g], which decreases for small and large values of ζ and increases in the middle part of the interval (thus capturing the previously discussed typical monotonicity with height of the heat forcing throughout the tropical troposphere). Such a polynomial can be easily obtained by means of Lagrange interpolation, using field data to identify the relevant values at the top and bottom of the troposphere, and at the two critical points in-between, where the monotonicity of the heat-forcing-profile changes.
whose general solution is given by the variation-of-constants formula (see [29]) where C and D are arbitrary functions, and The fundamental solution of the associated homogeneous constant-coefficient system d dt can be found directly as Corresponding solutions for W can be obtained from equation (5.3). Furthermore, we see that this describes the solution that arises when the external forcing is zero (F ≡ 0): it is unsteady but of a very specific form. This formulation, therefore, enables the solution to be found for any given forcing, the choice of which can be guided by the available data; this is clearly an area for extensive investigation.

Case study: the African Easterly Jet/Waves
During the Northern Hemisphere summer, strong heating of the Saharan region in North Africa, and the relatively cool and moist air to its south (in the Gulf of Guinea), creates a situation in which the usual north-south horizontal temperature gradient is reversed in the lower troposphere above the Sahel region, while in the upper troposphere the insolation-induced horizontal decrease of temperature with increasing latitude persists. Observations show the appearance of a strong westward jet, called the African Easterly Jet (AEJ), whose core on the western coast of Africa is near 15 • N, at a height of about 4 km (where the reversal of the temperature gradient in the middle troposphere occurs); see [30,31]. Denoting by a bracket the temporal mean (recall that 3 1 2 h corresponds to a unit interval for the non-dimensional variable t), in the time-periodic setting, we see that (5.1)-(5.3) yield the steady system and Since W = 0 at the ground level (ζ = −(1/2) cos 2 θ ), integration of equation (6.3) gives The relations (6.1)-(6.3) enable us to infer the basic dynamics of the AEJ from the observed behaviour of the heat forcing F. The heat-flux contribution from clouds is negligible, the main source being the ground-level temperature gradient, distributed in the lower troposphere by dry convection and diffusion (rather than by condensational heating); see the discussion in [30]. Going north from the Guinean coast towards the Sahara, in the region from 10 • N to 20 • N, the surface temperature increases gradually by about 10 • K (see the data in [30]). Thus ∂ F ∂θ > 0 and (6.1) predicts a westward zonal flow component for the AEJ. On the other hand, the surface temperature in the lower-troposphere Sahel region between 15 • W and 10 • E increases (weakly) eastwards (see the data in [30]). The fact that ∂ F /∂ϕ > 0 in the lower troposphere yields, from (6.2), that the meridional flow velocity of the AEJ is northward. Field data (see [30]) also show that near 15 • N, between 15 • W and 10 • E, the positive longitudinal gradient of the temperature increases with height near the ground and then starts to decrease towards the mid-troposphere altitude, where it reverses sign: the classic Hadley circulation coexists with a second but shallower overturning circulation in the lower part of the troposphere (see [32]). Due to (6.2), at a fixed latitude θ , we can track the variation of ∂ F /∂ϕ with height by the corresponding variation of the meridional velocity component, V. Field data for the AEJ (see [30,33]) show that V increases linearly with height in the lower third of the troposphere, from about 1-2 m s −1 near the ground, to a maximum of about 6-7 m s −1 at an altitude of about 4 km. Given that one non-dimensional speed unit corresponds to Ω H ∼ 1.16 m s −1 , this means that the vertical slope of ∂ F /∂ϕ near the ground is less than 15. Since sin(15 • ) ∼ 0.258, we see that, for the AEJ, 1/sin 2 θ exceeds the vertical slope of ∂ F /∂ϕ in the lower third of the troposphere. Consequently, (6.3) confirms that uplifting (W > 0) occurs near the ground. We have presented the main qualitative features of the AEJ, but detailed data about the temperature (which is readily available, but its use is outside the scope of this initial investigation) would provide the basis for a quantitative analysis, producing a more comprehensive description.
African easterly waves (AEWs) are westward moving oscillatory disturbances of the AEJ, initiated by mesoscale convective systems over Central Africa and propagating in the lower troposphere; see [34].  troposphere, thermally induced by the contrast between the warm Sahara Desert area to the north and the cool Gulf of Guinea to the south-upward to 6 km above the ground, with the core at about 3 km above the ground (where the vertically sheared AEJ is most intense; see [35]). The shorter waves propagating northward at about 3 m s −1 arise as perturbations of an ageostrophic meridional circulation in the upper troposphere, mainly in the region of 8-12 km above ground; see [36]. Note that the westward travelling African Easterly Waves, arising as lower-tropospheric threedimensional synoptic-scale disturbances of the AEJ (see [37]), with a maximal speed of about 11 m s −1 (attained around 4 km altitude) and wavelengths of about 2000 km, are not visible in the photograph. These more intricate wave patterns are identified by using theoretical predictions to interpret field data, for example, by performing a wavenumber-frequency spectrum analysis to isolate statistically significant spectral peaks that correspond to the available dispersion relations (see [38]). (Online version in colour.) the mean heat forcing F . With the mean flow determined by (6.1)-(6.3), the waves are therefore solutions of the system ∂ζ ∂t , whose general solution was described in §5c. Waves with shorther wavelengths than the AEWs are also observed (figure 4); the above approach also applies in this case.

Viscous flow
We return to the original set of equations that define the velocity field at this order, namely equations (4.12)-(4.14). An analysis of this system is essential if viscous effects are thought to play a role in the atmospheric flows of interest, which is certainly the case within the atmospheric boundary layer (see the discussion in [39]). The general procedure is exactly as before: solve the first two for the horizontal components, and then the third gives the vertical component (or provides an auxiliary condition, if we choose to set w 0 = 0). The simplest way to proceed, as was done in [3], is to introduce the complex velocity: equations (4.12) and (4.13) then give This can be rewritten in a slightly more convenient form by setting By means of a Liouville substitution, equation (7.2) can be transformed into a (complex-valued) Bessel equation (after the removal of a time-dependent factor), whose general solution can be expressed in terms of Bessel functions of the first kind (at least, for some reasonable choices of M); see [3]. It is clear that the development of these solutions is quite involved but readily accessible, and the results should provide important insights into these flows; this is another area that is left for future investigation.

Discussion
The ideas developed in this paper have shown that the theory presented in [3] can be extended to include time dependence (on a suitable time scale, of course). In particular, therefore, we invoke an asymptotic method that is driven by the thin-shell approximation, which represents a minimal assumption for the atmosphere enveloping the Earth. (It is expedient also to take the shape of the Earth as ellipsoidal, and use almost-spherical coordinates; the role of this approximation is simply to adjust the nature of the background state, without affecting the leading-order unsteady description of the atmosphere: the uncoupling at this order is an important property of the underlying system. Thus this aspect of the problem is included in the general formulation, but ignored in the main thesis that we propound here.) As we found in our earlier work [3], the real surprise is that the thin-shell approximation, without recourse to any other assumptions about the nature of the flow, leads to a complete description of the dominant dynamic and thermodynamic elements that are needed to describe the atmosphere. This earlier work showed how standard steady-state models that are used for specific problems in the atmosphere (e.g. Ekman and geostrophic flows, Hadley cell circulation) can be recovered from one set of simple governing equations, and with many accurate interpretations that underlie the flow configurations (such as the precise nature of the heating that is needed to drive the motion). In the current study, our set of governing equations, (4.12)-(4.14), with (4.2) and (4.4), provide the basis for a systematic study of unsteady motion and, most particularly, of various modes of wave propagation in the troposphere. The development presented here follows the path originated in [3]: a steady background state of the atmosphere (described by equations in which the velocity field appears at this same order) is perturbed by using the thin-shell parameter, with time included in the perturbation system. The time scale-a few hours-is that associated with the rotation of the Earth, which is the appropriate choice for reasonably large-scale wave-like motions of the atmosphere. The upshot is that we have a perturbation state that is unsteady, contains the main contributors to the dynamics and thermodynamics, as well as viscosity, compressibility, the Earth's rotation and the underlying spherical geometry. In addition, the formulation allows for the identification of the heat sources that are needed to initiate and maintain the motion, and this connection is altogether transparent. One main advantage of the careful and systematic construction of a set of equations with a robust pedigree is that we can hope to develop numerical solutions of this simpler system, driven by reliable data. Our systematic approach is based on theoretical first principles, augmented by phenomenological insights (regarding the nature of the heat forcing), and avoids ad hoc assumptions. This conceptually coherent modelling brings to light the significant dependency of the atmospheric flow on the heat forcing, avoiding adjustments of model parameterizations that often mask underlying deficiencies of models, which are not derived systematically from the governing equations, a procedure that may fail to capture the relevant physical processes even if, due to parameter calibration, they might exhibit agreement between simulation and data.
Our equations admit solutions that are harmonic in time, and so the majority of our examples are based on solutions that possess this property. A familiar model for simple, oscillatory flows involves the Brunt-Väisälä frequency (N), so we started by investigating how such flows appear in our system. On the one hand, we see that N 2 is directly and simply related to our forcing function (F 1 ) and so its form for various flows can be interpreted in terms of the associated heat sources. However, the standard analysis, which leads to the usual dispersion relation involving N 2 , cannot be obtained in any systematic way (as explained in §5a(i)): we conclude that any such 'derivations', and conclusions based on the results, remain highly suspect. Fortunately, there are many other avenues that we can explore which lead to relevant and reliable results, with applications to the motion of the atmosphere.
The main thrust of the examples that we have presented here is for inviscid flows; viscous effects have been included in the formulation, but are weak (and relevant only in thin boundary layers) and are therefore set aside in much of what we present. The examples are intended to show what is possible; clearly, many other choices can be made, suitable for other flow configurations. (It should be noted that, in all cases, the existence of a background state, which varies in ζ = gz − (1/2) cos 2 θ , ensures that there is always a variation in the vertical direction, with an associated distortion in the meridional direction.) The system that describes the dynamic-thermodynamic balance, with heat sources, admits solutions that are harmonic in t, ϕ and ζ , but not in θ: the dependence on θ necessarily involves some modulation in the meridional direction. With these points in mind, we have obtained a standing-wave solution, which oscillates in time only, as well as solutions that describe waves propagating in the meridional direction, and waves in the azimuthal direction. The former provides a model for waves that propagate towards the Equator (over the Arabian Sea); the latter is a solution of the type observed to be moving in an equatorial direction over the Indian Ocean. These solutions, and all similar ones, possess a vertical structure, provided in part by the existence of the background state, and all satisfying the condition W = 0 at the bottom of the atmosphere. We expect that many other examples, relevant to other modes and directions of propagation, with varying complexities of modulation in the meridional direction, can be identified; this is an area that needs further investigation.
Another example, which shows in greater depth the advantages of a general, all-embracing system of equations, arises when we consider waves that are trapped in the neighbourhood of the Equator. In contrast to the familiar approach, we need make no assumptions about the behaviour in the meridional direction: f -or β-plane approximations are altogether unnecessary. Simply by seeking a solution that imposes no motion in the meridional direction, together with a suitable heat forcing, leads to a solution that is restricted to a neighbourhood of the Equator, the trapping being represented by a power-law behaviour (and not exponential as in the simple theories). This example, perhaps more so than the previous ones described above, shows one of the significant advantages of our formulation: we can identify and interpret the heat sources required to maintain the motions.
One final example, discussed in some detail, shows how a time-averaged version of the equations, for periodic solutions, can be used to provide some general observations without recourse to the construction of explicit solutions. In particular, this procedure was applied to a description of the AEJ, showing how the temperature profile in the troposphere necessarily produces the observed properties of this atmospheric flow. Indeed, this calculation can be There are many different types of unsteady motion in the atmosphere (including various instabilities) see [40]; we have chosen, as simple examples, a few modes of wave propagation (of which there are many) in this initial investigation. These have been chosen to demonstrate how details of the motion can be extracted and interpreted and, where appropriate, how general properties can be identified without the need for explicit solutions. In addition to these examples, we have also shown how the general time-dependent problem can be formulated, avoiding any assumptions about the nature of the time-like behaviour. Indeed, if the heat forcing is given, then the complete dynamical structure of the atmosphere can be calculated (at leading order). Furthermore, as with all the examples discussed here, the specification of the heat forcing, F 1 , can be used to provide the thermodynamic elements associated with the motion and lead to the interpretation of the required heat sources. All this, we emphasize, is based on a single set of simplified equations that are robustly connected to the underlying, governing equations for a fluid. Finally, we have outlined how the corresponding viscous problem can be solved, although the details are less readily accessible; more information on the viscous problem can be found in [3].
What we have developed and described here should be, we submit, the basis for the study of unsteady atmospheric flows. On the one hand, we have a set of equations derived carefully (using precise asymptotic methods) incorporating minimal simplifying assumptions; this system provides a reliable starting point for investigation, analysis and interpretation. On the other hand, particularly with the ready availability of extensive and reliable data, these equations can be used to generate numerical solutions, which, in turn, can become the seed for numerical studies of the original, full set of equations. In either case, we have shown how familiar results, and new results, can be obtained directly and systematically from the underlying governing equations. Furthermore, any simpler model should be tested against this system: can such models be obtained by making additional assumptions consistent with these equations, working altogether systematically? We submit that the only reliable validation of model equations is their (asymptotic) derivation from a set of general governing equations, along the lines that we have presented here.
Data accessibility. All data are provided in full in the paper. Authors' contributions. The authors contributed equally to this study. The authors contributed equally to this study, conceiving its scientific content, collecting the data and drafting the manuscript. Both authors gave final approval for publication and agree to be held accountable for the work performed therein.