Variations in bottom topography may result in second-class (low-frequency) waves, similar to planetary Rossby waves (e.g. Rhines, 1969a,b, 1989; Csanady, 1982; Pedlosky, 2003). The similarity between the two wave types becomes evident upon considering the general case in which potential vorticity varies with both the Coriolis parameter and water depth. Due to this similarity, rotating tank experiments with varying bottom topography are often used in the study of many aspects of the two types of (Rossby) waves. The potential vorticity considerations make it clear that in the Northern Hemisphere, shallower depths correspond to higher latitudes and thus topographic Rossby waves propagate with the shallower depths on their right, similar to the westward propagation of planetary Rossby waves under the influence of varying Coriolis parameter.
Topographic Rossby waves are of great importance in coastal ocean and lake dynamics. There are numerous observational indications for their existence (e.g. Csanady, 1973, 1976; Saylor et al., 1980; Raudsepp et al., 2003), and their unique dynamics have provided explanation of observed features throughout the world ocean (e.g. Ball, 1963, 1965; Shilo et al., 2007, 2008). Most theories of topographic Rossby waves have concentrated on a particular depth profile (mostly linear and rarely quadratic or exponential) that are amenable to analytical solutions, while general depth profiles received little attention. Here, we propose a simple approximation for topographic wave dynamics in an arbitrary bottom profile in one direction, assuming only that the bottom topography is slowly varying. With this approximation it is possible to estimate the evolution of an initial perturbation, even when the eigenstates of the system are not known. The approximation may be used to study more general coastal, and open ocean, topographic wave dynamics.
The paper is organised as follows: in Section 2.1, we develop the new approximation in Cartesian coordinates and then compare the analytical approximation to numerical simulation results obtained by the MITgcm (Section 2.2). In Section 3, we propose a similar approximation for a polar (cylindrical) coordinate system and a similar comparison to the numerical simulation. We conclude the study with a short discussion in Section 4.
We start from the linearised shallow water equations (LSWEs) in a channel of varying depth (see Fig. 1) on the f-plane (see Gill, 1982):
By differentiating eqs. (1) and (2) with respect to t, it is possible to express the zonal and meridional velocities, u, v, as functions of η:
We non-dimensionalise the equation by scaling the variables as follows: , , and , where the ‘hat’ indicates non-dimensional variables and is the Rossby radius of deformation. To complete the scaling, we also use , although this scaling is not necessary for the mathematical treatment presented below. α is a parameter that characterises the bottom slope and can be chosen to be the maximum of ∣b′(y)∣. Note that the proposed scaling is invalid in the absence of bottom slope, in which case one expects to find only Kelvin and Poincaré waves. α is a small parameter and thus the above scaling of y and t implies ‘shrinking’ of these variables; this indicates that the underlying assumption of this scaling is long meridional wavelength (either compared with the zonal wavelength or compared with the Rossby radius of deformation) and slow wave frequency (compared with f_{0}). The idea behind the different scaling of the zonal and meridional coordinates is the asymmetry of the system, in which the bottom topography is constant in x but varies with y. As it is well known and consistent with the numerical results presented below, topographic waves propagate mainly along isobaths and hardly across bathymetry lines.
For convenience, we now drop the ‘hat’ sign and, unless indicated, the variables are non-dimensional. Under the above scaling, the non-dimensional form of eq. (7) is:
Next, we assume that the solution to eq. (8) is of the form:
To allow for a perturbation series expansion, is assumed to be of the form:
To evaluate the velocities we first apply the above scaling in eqs. (4) and (5) and scale the velocity components on , the speed of gravity waves, that is and . The non-dimensional counterparts of eqs. (4) and (5) are:
It follows from the fact that ω depends on y that the phase speed, , varies with y, depending on the local slope, b′(y), and on the water depth, ; this is true regardless of the specific bottom topography profile, b(y). Thus, different parts of an initial perturbation will propagate with different phase speeds. This behaviour will be demonstrated and explained below.
It is possible to find the specific bottom topography for which the (zeroth-order) phase speed is constant, that is . Consistent with previous studies (Buchwald and Adams, 1968; Gill, 1982), this specific profile is an exponential profile.
For future reference we note that the dimensional frequency associated with the dispersion relation (12) is given by:
The LSWEs (1)–(3) can be rescaled using the above scaling. Then, the small parameter multiplies only the time derivative of the zonal momentum equation. Thus, the proposed scaling is relevant when the time derivative in eq. (1) is small compared with the other term in this equation. Such a situation occurs when the wave frequency is small. According to eq. (20), the wave frequency is small when the bottom topography is slowly varying; we thus associate the small parameter α with the bottom slope. In addition, when both the zonal and meridional wavelengths are long compared with the Rossby radius of deformation, or when the water depth approaches zero (i.e. ), the dispersion relation can be approximated as . In this limit, the wave is non-dispersive in the zonal direction and long zonal wavelength (small k) results in a smaller frequency and, hence, more accurate perturbation expansion. This limit is the analogue of the extratropical planetary Rossby waves frequency under the long-wave approximation, for which the wave frequency is . Lastly, when the Rossby radius of deformation is much greater than the zonal and meridional wavelengths and the water depth is sufficiently large, the dispersion relation [eq. (20)] can be approximated as . In this limit, short zonal waves (large k) result in smaller wave frequency and hence better perturbation expansion. This limit is the analogue of equatorial planetary Rossby wave frequency, where .
Using eqs. (11) and (13) and writing a differential equation for each order of approximation, we obtain the following general form for :
A few comments regarding the approximations presented above are worth noting.
First, the different Ω_{n}'s contain both real and imaginary terms. The imaginary terms ( where j is odd) will lead to exponential growth or decay with time, depending on the sign of . Even if the infinite sum in eq. (13) converges, each order of approximation may contain secular terms and, thus, the approximation is valid for times such that the highest order term is smaller than its preceding term, that is .
Second, the coefficients depend on , b and ω that are functions of y. In addition, in the derivatives of are divided by , such that singular points may occur when decays to zero slower than exponential (while for an exponential decay, the ratio between the derivative of and is constant). This limitation does not affect the solution outside the vicinity of the singularity. Such singular points can be seen in the numerical example discussed below, yet, as expected, other regions are not affected by these singularities.
Third, the higher-order approximations result in a non-typical wave propagation, that is a different evolution than that of a monochromatic wave propagation where ω is constant. Thus, when computing the phase speed using, for example, longitude-time (Hovemöller) diagram, the phase speed changes with time (i.e. the crest or trough lines of the longitude-time diagram are not linear). We observe such changes of phase speed with time in the numerical simulation–this support the proposed analytical result. The reason for this cross-bathymetry phase speed dependence might be due to the set of eigenstates that make up the initial perturbation, each covering its own spatial range, that propagate in its corresponding phase speed. When starting from a true eigenstate of the system (which is hard to find for non-trivial bottom topography), there will be no change in wave frequency neither with time nor in space. The advantage of the proposed approximation is that there is no need to know the eigenstates of the system (corresponding to the bottom topography) in order to evaluate the short time evolution of an initial perturbation.
We use the Massachusetts Institute of Technology general circulation model (MITgcm, see Adcroft et al., 2002, 2004) to validate the analytical approximation. We consider a channel on the f-plane (f_{0}=10^{−4} s^{−1}) covering 120 km in the zonal and meridional directions, with periodic boundary conditions in the zonal direction and walls in the meridional boundaries. Free-slip boundary conditions and vanishing meridional velocities at the channel walls are assumed. The lateral resolution is 2 km. There is only one layer and, since the MITgcm can handle partial cells, one can consider depth variations even in a single layer ocean. The maximal depth in the domain is H_{0}=170 m, at the southern wall, and it decreases monotonically to a depth of 100 m at the northern edge, such that ɛ ≈ 3.4×10^{−7}; linear and quadratic depth profiles are considered (Fig. 1). We use these two profiles in the numerical verification of the analytical results to demonstrate the dependence of the zonal phase speed on both the depth gradient and the actual depth [eq. (20)], so the wave characteristics in the two depth profiles are different even though the overall bottom difference across the 120 km wide channel is identical. The gravitation constant used is g=9.81 m s^{−2}. The simulations were run for 1 week, a typical time of topographic Rossby waves evolution.
The initial conditions for the free surface of the two simulations are sinusoidal functions in both directions (two cycles in the zonal direction and half a cycle in the meridional direction), as shown in Fig. 2a; the initial velocities were calculated based on the geostrophic approximation (Fig. 2b, c), which is not an eigensolution of the governing equations. As expected from the analytical considerations presented above, the wave propagates to the West, such that the shallower depths are on the right relative to the direction of propagation. However, there is a noticeable difference between the two cases of linear and quadratic bottom topography: in the linear topography, the northward parts of the wave propagated faster than the southern parts, whereas the opposite occurred in the quadratic topography (Fig. 3). The agreement between the analytical approximation and the numerical results is good, both for the zeroth-order (Fig. 3c, d) and first-order (Fig. 3e, f) approximations. Similar comparisons between the numerical results and the analytical approximations for the zonal and meridional velocities are shown in Figs. 4 and 5, respectively. Here, the first-order approximation is closer to the numerical results. Singular behaviour may occur close to the channel's walls (Fig. 3e, 4e, f, 5e), as expected from eq. (23).
In Fig. 6, we present the numerically calculated latitudinal variation of the zonal phase speed in the two cases and the corresponding analytical estimates. The phase speed is negative, indicating that indeed the wave travels to the West (such that the shallower depths are on the right). The observation of faster (i.e. more negative) zonal phase speeds at higher y for the linear topography, and the vice versa for the quadratic topography is clearly evident. The agreement between the numerically computed zonal phase speed and the analytical speed is good, where, as expected, the first-order approximation yields a closer agreement to the numerical results. We note that the agreement between the numerical and analytical approximations becomes less accurate with time due to the presence of the secular term in the analytical solution; see Section 2.1.2. As predicted by the analytical approximation, the numerically calculated phase speed depends on time, and the one shown in Fig. 6 is based on the mean phase speed of the first 14 d of the simulations.
The initial states used here are eigenmodes (harmonic functions) of the quasi-geostrophic approximation. These do not preserve their structure with time, as the perturbation propagates with different speeds for different latitudes; thus, the quasi-geostrophic solution provides an inaccurate approximation for the set-up studied here. This inaccuracy is most probably due to the more significant bottom slope, which cannot be considered to be small as required in the quasi-geostrophic approximation.
The reason for the difference between the zonal phase speeds of the two runs is clear from eq. (20). For the linear topography case, only the local water depth affects the zonal phase speed and thus, as the depth decreases, the denominator becomes smaller and hence the phase speed increases. In contrast, for the quadratic topography, the local bottom slope appears in the numerator and as it decreases towards higher y, so does the phase speed. In this case, the effect of the local bottom slope in reducing the zonal phase speed towards larger y outweighs the effect of local depth that tends to increasing the zonal phase speed at larger y.
It is possible to approximate the eigenstates of LSWEs in the case of a linear depth profile using Airy Functions and Parabolic Cylinder Functions Cohen et al. [(2010); see also the planetary counterpart in Paldor and Sigalov (2008)]. These functions oscillate at large depths and decay (faster than exponential) at shallower depths. Higher eigenmodes extend farther into shallower depths, where each of them propagates with its own phase speed. These straightforward considerations can explain the different phase speeds for different y's by viewing an arbitrary initial perturbation as a collection of eigenstates, each propagating with its own speed in a different subregime of the domain, leading to different phase speeds at different latitudes.
Although we consider here the case of topographic wave within a channel, our approximation is also relevant for coastal shelf and open ocean topographic wave dynamics. The reason is that the problem that we solve is essentially an initial value problem in the meridional direction. Thus, the channel width can be as large as desired—it is possible to approximate the dynamics of initially meridionally concentrated perturbation within the channel. The perturbation can be located far from the channel walls, to mimic open ocean dynamics, or close to the shallow boundary to mimic coastal shelf dynamics. We have verified these claims by performing additional numerical simulations (results not shown) with a channel that is four times wider than the one considered in the numerical simulations discussed above. We found similar agreement between the analytical approximation and the numerical results, with better performance for the perturbation that is located far from the channel walls. As expected, the first-order approximation was more accurate than the zeroth-order one.
We performed an additional numerical simulation with a sinusoidal (two waves in the y-direction) bottom topography. Here, only the zeroth-order approximation exhibited a reasonable agreement with the numerical simulations. We speculate that the secular terms of the approximation diverge faster for rough topography profile.
It is unlikely to find exact analytical solutions of the LSWEs for arbitrary bottom topography. The analytical approximation presented here provides the wave structure and zonal phase speed for bottom topography that is uniform in the zonal direction (i.e. the long-shelf direction), but is arbitrary in the meridional (cross-shelf) direction, and the approximation is valid for small ɛ values. Note that the time of validity of the approximate solution becomes shorter as ɛ increases. In addition, the details of the wave's meridional structure are not required for calculating the (zeroth-order) dispersion relation.
The derivations of the previous section can be straightforwardly extended to basins other than the channel studied above. Since circular basins seem to be closer to the geometry of a typical lake geometry, we find it useful to discuss this geometry in details.
The LSWEs in polar coordinates on the f-plane and with varying depth H_{0} − b(r) are given by:
where v_{r} and v_{θ} are the radial and azimuthal velocities, η is the free surface displacement, f_{0} is the constant Coriolis parameter, g is the gravitation constant, H_{0} is the unperturbed water depth and b(r) is the (radial) bottom topography. The independent variables are the radial, azimuthal and time coordinates (r, θ and t, respectively). Note that we assume here, similar to the previous section, that the depth changes only in the radial direction, r.
The velocities may be expressed in terms of the free surface:
We now perform the following scaling to switch to non-dimensional (or scaled) variables: , , , and , where the ‘hat’ indicates non-dimensional (scaled) variable. α is a parameter that characterises the mean bottom slope. We now omit the ‘hat’ sign and, unless otherwise indicated, the different variables designate non-dimensional variables.
Using the above scaling, eq. (27) transforms to:
We next assume that is of the form:
It is possible to find the zero-order velocities , using the scaling and and using eqs. (25) and (26).
Since ω is a function of r, the phase speed varies as a function of r, depending on the local slope b′(r) and local water depth . There exists a bottom topography profile for which the (zero-order) phase speed is constant, that is . This particular profile is given by:
The dimensional zeroth-order frequency is given by:
In polar coordinates, the higher-order approximation is essentially the same as that in Cartesian coordinates given in eqs. (21) and (22), with the trivial replacement of y by r, and using eqs. (31) and (33).
The coefficients of the first-order approximations for Ω are:
As in Section 2.2, we used the MITgcm for the numerical verification of our analytical approximations. The spatial dimensions of the square used in the simulation are 120×120 m with a resolution of 2 km. The Coriolis parameter is f_{0}=10^{−4} s^{−1}, the gravitation constant is g=9.81 m s^{−2} and H_{0}=90 m. The (dimensional) b(r) is shown in Fig. 7, and it follows eq. (38) so that the zeroth-order frequency is constant and corresponds to a topographic wave rotation of one cycle per 10 d.
In Fig. 8, we present a comparison between the numerical simulation and the zeroth-order approximation of the free surface. The simulation is initiated with a Gaussian profile in the radial direction, multiplied by a sinusoidal function in the azimuthal direction, with a wave number, k=2 (see Fig. 8a). The initial velocities were chosen to be in geostrophic balance with the free surface and the initial perturbation rotated clockwise. The numerically calculated free surface for t=3 and 8 d is shown in Fig. 8b, c. It is evident that the initial perturbation rotates with constant angular velocity for different values of r, confirming the analytical prediction of constant ω for the profile given in eq. (38). In Fig. 8d, the zeroth-order approximation free surface for t=8 d is depicted. The orientation of the numerical and analytical free surface for t=8 d (Fig. 8c, d) is similar; yet, while the zeroth-order approximation exhibits only rotation of the initial perturbation, the numerically calculated free surface exhibits also radial propagation. This difference may be attributed to viscosity and non-linear terms that are included in the numerical model. We note that the first-order approximation exhibits dependence of wave frequency on the radial coordinate; this difference between the analytical first-order approximation and the numerical results may be linked to higher wave frequency close to the periphery of the basin, making the perturbation expansion less accurate.
We have also compared numerical simulations (results not shown) to the zeroth-order approximation for other bottom topography profiles. As expected, the numerical simulations of these cases indicated different angular rotations at different values of r, that is ω does vary with r. The remarks listed in Section 2.1.3 for a periodic zonal channel are also relevant to the current case of a circular basin.
We have studied topographic Rossby waves dynamics on the f-plane using the LSWEs. A slowly varying arbitrary bottom topography profile in one of the directions was assumed. We propose a new approximation, according to which the time and cross-bathymetry directions are considered to be slow and long, that is changes in the cross-bathymetry direction are small. The small non-dimensional parameter is the square of a number representing the slope of the bottom topography (which we choose to be the maximum of the absolute value of the slope). It then follows that the cross-bathymetry coordinate appears as a parameter, that is the underlying differential equation, for each order of approximation, does not contain derivatives with respect to the cross-bathymetry coordinate. This equation leads to a wave frequency (and hence phase speed) that explicitly depends on the cross-bathymetry coordinate where the wave propagation depends on the local slope of the bottom topography and on the local water depth. [We note that it is possible to apply an alternative expansion in which the wave number depends on y (i.e. k(y)) instead of the ω(y) dependence. We however find the latter easier to interpret and implement as an initial value problem.] The analytical approximation is valid for slowly varying arbitrary initial conditions in the cross-bathymetry direction that satisfy the boundary conditions; the analytical approximation exhibits close similarity with the numerical simulations carried out using the MITgcm. The proposed approximation may be relevant to coastal shelf and open ocean wave dynamics since the channel can be very wide and the evolution of a meridionally localised perturbation is very well approximated by the analytical theory. Both linear and quadratic bottom topographies were used in the simulations, both spanning the same overall height difference. Yet, the linear bottom topography yielded faster wave propagation towards shallower depths, while the quadratic bottom topography exhibited the opposite behaviour, due to the smaller local slope at shallower depths.
We have developed a similar approximation for a circular basin using the LSWEs in polar coordinates on the f-plane and a general bottom topography profile in the radial direction. The approximation may be found useful in the study of lakes with (nearly) circular symmetry. This case shares similar characteristics with a channel of periodic boundary conditions. We have found a specific form of the bottom topography that yields the same angular rotation as a function of the radial dimension; constant (zeroth-order) topographic wave rotation may be just clockwise rotation.
It had been shown previously that when changes in bottom topography are very moderate with respect to the changes in the wavelength, it is possible to replace the mean water depth by the local water depth (e.g. Pedlosky, 2003). In essence, our proposed approximation is based on similar arguments, formulating the approximation rigorously and developing also the first-order approximation. The main advantage of the proposed approximation over existing ones is the ability to find the dispersion relation for arbitrary bottom topography. Extension of the present work to basins with elliptical shape is possible, which may be more relevant to lake dynamics in some cases.
Adcroft, A, Campin, J.-M, Heimbach, P, Hill, C and Marshall, J. 2002. MITgcm Release Manual. MIT/EAPS: CambridgeMA. Online at: http://mitgcm.org/sealion/online_documents/manual.html.
AdcroftA. CampinJ.-M. HillC. MarshallJ. Implementation of an atmosphere-ocean general circulation model on the expanded spherical cube. Mon. Weather Rev. 2004; 132(12): 2845–2863. https://doi.org/10.3402/tellusa.v64i0.18160.
AdemJ. A series solution for the barotropic vorticity equation and its application in the study of atmospheric vortices. Tellus. 1956; 8(3): 364–372. https://doi.org/10.3402/tellusa.v64i0.18160.
Ashkenazy, Y, Paldor, N and Zarmi, Y. 2011. On the meridional structure of extra-tropical Rossby waves. Tellus A. 63(4), 817–827. https://doi.org/10.3402/tellusa.v64i0.18160.
BallF. K. The effect of rotation on the simpler modes of motion of a liquid in an elliptic paraboloid. J. Fluid Mech. 1963; 22: 529–545. https://doi.org/10.3402/tellusa.v64i0.18160.
BallF. K. Second-class motions of a shallow liquid. J. Fluid Mech. 1965; 23: 545–561. https://doi.org/10.3402/tellusa.v64i0.18160.
BuchwaldV. T. AdamsJ. K. Propagation of continental shelf waves. Proc. R. Soc. Lond. A. 1968; 305(1481): 235–250. https://doi.org/10.3402/tellusa.v64i0.18160.
CohenY. PaldorN. SommeriaJ. Laboratory experiments and a non-harmonic theory for topographic Rossby waves over a linearly sloping bottom on the f-plane. J. Fluid Mech. 2010; 645: 479–496. https://doi.org/10.3402/tellusa.v64i0.18160.
CsanadyG. T. Wind-induced barotropic motions in long lakes. J. Phys. Oceanogr. 1973; 3: 429–438.
CsanadyG. T. Topographic waves in Lake Ontario. J. Phys. Oceanogr. 1976; 6: 93–103.
CsanadyG. T. Circulation in the Coastal Ocean. Reidel: Dordrecht, 1982
GillA. E. Atmosphere–Ocean Dynamics. Academic Press: London, 1982
PaldorN. SigalovA. Trapped waves on the mid-latitude beta-plane. Tellus. 2008; 60A: 742–748.
PedloskyJ. Waves in the Ocean and Atmosphere. Springer-Verlag: Berlin Heidelberg, New York, 2003
PrimeauF. Long Rossby wave basin-crossing time and the resonance of low-frequency basin modes. J. Phys. Oceanogr. 2002; 32: 2652–2665. https://doi.org/10.3402/tellusa.v64i0.18160.
RaudseppU. BeletskyD. SchwabD. Basin-scale topographic waves in the Gulf of Riga. J. Phys. Oceanogr. 2003; 33: 1129–1140.
RhinesP. B. Slow oscillations in an ocean of varying depth. 1. Abrupt topography. J. Fluid Mech. 1969a; 37: 161–189. https://doi.org/10.3402/tellusa.v64i0.18160.
RhinesP. B. Slow oscillations in an ocean of varying depth. 2. Islands and seamounts. J. Fluid Mech. 1969b; 37: 191–205. https://doi.org/10.3402/tellusa.v64i0.18160.
RhinesP. B. Deep planetary circulation and topography: simple models of midocean flows. J. Phys. Oceanogr. 1989; 19: 1449–1470.
SaylorJ. H. HuangC. K. ReidR. Vortex modes in southern Lake Michigan. J. Phys. Oceanogr. 1980; 10: 1814–1823.
Shilo, E, Ashkenazy, Y, Rimmer, A, Assouline, S, Katsafados, P and co-authors. 2007. The effect of wind variability on topographic waves: Lake Kinneret case. J. Geophys. Res. 112, C12024.
ShiloE. AshkenazyY. RimmerA. AssoulineS. MahrerY. Wind spatial variability and topographic wave frequency. J. Phys. Oceanogr. 2008; 38: 2085–2096. https://doi.org/10.3402/tellusa.v64i0.18160.
VallisG. K. Atmospheric and Oceanic Fluid Dynamics. Cambridge University Press: Cambridge, 2006