On the role of domain aspect ratio in the westward intensification of wind-driven surface ocean circulation

The two seminal studies on westward intensification, carried out by Stommel and Munk over 70 years ago, are revisited to elucidate the role of the domain aspect ratio (i.e., meridional to zonal extents of the basin) in determining the transport of the western boundary current (WBC). We examine the general mathematical properties of the two models by transforming them to differential problems that contain only two parameters – the domain aspect ratio and the non-dimensional damping (viscous) coefficient. Explicit analytical expressions are obtained from solutions of the non-dimensional vorticity equations and verified by long-term numerical simulations of the corresponding time-dependent equations. The analytical expressions as well as the simulations imply that in Stommel's model both the domain aspect ratio and the damping parameter contribute to the non-dimensional transport of the WBC. However, the transport increases as a cubic power of the aspect ratio and decreases linearly with the damping coefficient. On the other hand, in Munk's model the WBC's transport increases linearly with the domain aspect ratio, while the damping coefficient plays a minor role only. This finding is employed to explain the weak WBC in the South Pacific. The decrease in transport of the WBC for small-domain aspect ratio results from the decrease in Sverdrup transport in the basin's interior because the meridional shear of the zonal velocity cannot be neglected as an additional vorticity term.

1 Introduction

As was noted by Henry Stommel, in the opening sentence of his seminal 1948 study, “Perhaps the most striking feature of the general oceanic wind-driven circulation is the intense crowding of streamlines near the western borders of the oceans”. These strong and narrow poleward-directed currents, often referred to as “western boundary currents” (WBCs), counterbalance the weak and wide equatorward (Sverdrup) flow in the interior of the basin. In the North Atlantic this current is the Gulf Stream, and it was known to oceanographers and explorers for a few centuries – see Stommel (1958) for a historical review. Similar WBCs exist in other basins as well, and these include the Kuroshio in the North Pacific and the Brazil Current in the South Atlantic. These currents transport large amounts of heat from low to high latitudes, thus playing an important role in the climate system. The winds overlying, though, are easterlies along the Equator (the trade winds) and westerlies around 40^∘ N. There are no strong northward winds along the western boundaries of the ocean basins and, as is well understood now, the WBCs are not obviously correlated with the overlying wind patterns. Interestingly, two such WBCs lie in the Pacific viz. the Kuroshio and the East Australian Current (EAC). Both the Kuroshio and the EAC are centered close to 26^∘ latitude in their respective hemispheres, are driven by similar wind stresses and are adjacent to a ∼ 2000 km long coastline. Despite these structural similarities, the maximal volumetric transport of the Kuroshio Current is 55 Sv (1 Sv =10⁶ m³ s⁻¹) (Qiu, 2019), whereas that of the EAC is around 30 Sv (Archer et al., 2017). The maximum velocity that EAC attains is also substantially smaller than that of the Kuroshio (Campisi-Pinto et al., 2020).

Stommel, apparently in his first oceanography paper (Stommel, 1948, hereafter referred to as S48) was the first to formulate a simple, yet comprehensive, mathematical model of the WBCs (see for example Kunzig, 1999). S48 is now regarded as a seminal paper in theoretical physical oceanography (e.g., http://empslocal.ex.ac.uk/people/staff/gv219/classics.d/oceanic.html, last access: 10 February 2021). S48's model probably provides the simplest explanation for the existence of WBCs: in this linear and frictional model on the β plane the ocean is taken to be a flat bottom rectangle forced by a cos(latitude)-dependent zonal wind pattern. Walter Munk further extended this work to a different frictional (viscous) parameterization and a more general form of the wind stress (Munk, 1950, hereafter referred to as M50).

In the last 70 years, both models have been modified and extended to further explore the phenomenon of westward intensification in different settings or to evaluate the importance of different specific processes and terms in the governing equations (Munk and Carrier, 1950; Veronis, 1966a, b; Pedlosky, 1987; Vallis, 2017, and references therein).

As in S48 and M50, a large number of these subsequent studies employed the dimensional form of the governing equations which are the time-independent rotating linearized shallow water equations compounded by friction and forcing. These dimensional models include numerous parameters: the zonal and meridional extents of the basin; either the coefficient of linear drag (i.e., the coefficient in the Rayleigh frictional term) or the kinematic eddy viscosity (i.e., the coefficient in parameterization of the viscous term); the amplitude (and possibly meridional structure) of the wind stress; the gradient of Coriolis frequency (β effect). On the other hand, a few studies (Welander, 1976; Bye and Veronis, 1979) employed the alternate, concise approach of non-dimensionalizing the governing equation (or the vorticity equation) to investigate the depth-averaged wind-driven ocean circulation. The non-dimensional approach not only simplifies the problem by reducing the number of dimensional parameters in the model to fewer non-dimensional ones but also brings out some salient features associated with the problem which are difficult to unveil in the dimensional formulation.

By employing a non-dimensional approach, Welander (1976) successfully identified a zonally uniform regime in both S48's and M50's models of wind-driven ocean circulation, and using the same approach, Bye and Veronis (1979) derived a correction to the Sverdrup transport in S48's model. The aforementioned studies highlighted the importance of the ratio between meridional and zonal extents of the basin as one of the two fundamental parameters in both S48's and M50's models. The aim of this study is to further elaborate on the role of the domain aspect ratio (defined here as the ratio between the basin's meridional and zonal extents) in S48's and M50's models of westward intensification. In particular, we examine the role of domain aspect ratio in the transport of the WBC as was first hypothesized by Bye and Veronis (1979) in the context of S48's model: “… the tendency of north–south diffusive processes to be more significant in basins with a large [small in the present scaling] aspect ratio makes sense physically and may play a quantitative role in the transport of the western boundary current”. We also examine the relevance of our results to the observed difference in strengths of the five WBCs in the world ocean.

The paper is organized as follows. Section 2 outlines our proposed scaling (which is slightly different from the one employed in Welander, 1976; Bye and Veronis, 1979) that reduces the number of parameters in the vorticity equations corresponding to S48's and M50's models from five dimensional ones to two non-dimensional ones – one of which is the domain aspect ratio (the other is damping). The solution for the stream function in the two cases is outlined, and using this we obtain the expression for the non-dimensional transport of the WBC in both S48's and M50's models. The applicability of the analytical expression of transport for relevant values of the model parameters is validated in Sect. 3 by simulating the time-dependent equations. We discuss the results and conclude in Sect. 4. We also note that there were some typos in the expressions of zonal velocity and sea surface height (but not the stream function itself) in S48, and for completeness we list them in Appendix A. These typos do not change the scientific conclusions drawn in S48.

2 The two-parameter differential problems, their solutions and the transports of the WBC

2.1 S48's non-dimensional counterpart

S48's dimensional vorticity equation for the spatial structure of the stream function, ψ, is given by

$$\begin{array}{}\text{(1)}& r{\mathrm{\nabla }}^{\mathrm{2}}\mathit{\psi }+\mathit{\beta }{\displaystyle \frac{\partial \mathit{\psi }}{\partial x}}={\mathit{\tau }}_{\mathrm{0}}{\displaystyle \frac{\mathit{\pi }}{{\mathit{\rho }}_{\mathrm{0}}{H}_{\mathrm{0}}{L}_y}}\sin \left({\displaystyle \frac{\mathit{\pi }y}{L_y}}\right),\end{array}$$

where r is the Rayleigh friction coefficient, β is the meridional gradient of the Coriolis frequency, and τ₀ is the amplitude of wind stress. The operator ∇² is the two-dimensional Laplacian, H₀ is the mean depth of the barotropic ocean with density ρ₀, and L_y is the meridional dimension (and L_x is the zonal dimension) of the basin. The velocity components in the zonal and meridional directions, u and v, are related to the stream function via $u=\frac{\partial \mathit{\psi }}{\partial y}$ and $v=-\frac{\partial \mathit{\psi }}{\partial x}$.

We begin by scaling Eq. (1) as follows: x (the zonal coordinate) on L_x; y (the meridional coordinate) on L_y and ψ on $\mathit{\gamma }\mathit{\beta }{L}_y^{\mathrm{3}}$, where $\mathit{\gamma }={\mathit{\tau }}_{\mathrm{0}}\frac{\mathit{\pi }}{L_y}\left(\frac{L_x}{{\mathit{\rho }}_{\mathrm{0}}{H}_{\mathrm{0}}{\mathit{\beta }}^{\mathrm{2}}{L}_y^{\mathrm{3}}}\right)$ is the non-dimensional amplitude of the wind stress curl. With this scaling the non-dimensional form of S48's vorticity equation is

$$\begin{array}{}\text{(2)}& {\displaystyle \frac{\mathit{ϵ}}{{\mathit{\delta }}^{\mathrm{2}}}}{\mathrm{\nabla }}^{\mathrm{2}}\mathit{\psi }+{\displaystyle \frac{\partial \mathit{\psi }}{\partial x}}=\sin \left(\mathit{\pi }y\right),\end{array}$$

where

$$\begin{array}{}\text{(3)}& \mathit{ϵ}={\displaystyle \frac{r}{\mathit{\beta }{L}_x}},{\mathrm{\nabla }}^{\mathrm{2}}={\mathit{\delta }}^{\mathrm{2}}{\displaystyle \frac{{\partial }^{\mathrm{2}}}{\partial x^{\mathrm{2}}}}+{\displaystyle \frac{{\partial }^{\mathrm{2}}}{\partial y^{\mathrm{2}}}}.\end{array}$$

From this point onwards, both the variables and the operators in the differential equation(s) are non-dimensional, while dimensional quantities will be accompanied by an asterisk (^*). Here ∇² is the non-dimensional Laplacian, $\mathit{\delta }=\frac{L_y}{L_x}$ is the ratio of meridional and zonal extents of the basin (referred to as the domain aspect ratio), and ϵ is the non-dimensional width of the WBC (and also a proxy of the damping). The definition of the non-dimensional stream function implies that the zonal velocity and the meridional velocity are given by $u=\frac{\partial \mathit{\psi }}{\partial y}$ and $v=-\mathit{\delta }\frac{\partial \mathit{\psi }}{\partial x}$, respectively. It is evident from Eqs. (2) and (3) that the two parameters, ϵ and δ, govern the structure of the flow in the basin. The no-normal-flow conditions at the basin's boundaries mandate that the stream function ψ satisfies the boundary conditions: $\mathit{\psi }(x,\mathrm{0})=\mathit{\psi }(\mathrm{1},y)=\mathit{\psi }(x,\mathrm{1})=\mathit{\psi }(\mathrm{0},y)=\mathrm{0}$. We note that the term domain aspect ratio in the present study, $\mathit{\delta }=\frac{L_y}{L_x}$, is the inverse of the domain aspect ratio used in Bye and Veronis (1979), who derived a non-dimensional equation similar to Eq. (2).

As has been stated earlier, the non-dimensional formulation lumps the five dimensional parameters in S48's model – zonal and meridional extent of the basin, gradient of Coriolis frequency, wind stress amplitude and Rayleigh friction coefficient – into just two non-dimensional ones: ϵ and δ (both of which appear only in the first term of Eq. 2). Following S48, an explicit expression of the solution for ψ in Eq. (2) is given by

$$\begin{array}{}\text{(4)}& \mathit{\psi }(x,y)={\displaystyle \frac{{\mathit{\delta }}^{\mathrm{2}}}{\mathit{ϵ}{\mathit{\pi }}^{\mathrm{2}}}}\sin \left(\mathit{\pi }y\right)\left(p{e}^{Ax}+q{e}^{Bx}-\mathrm{1}\right),\end{array}$$

where

$$\begin{array}{ll}{\displaystyle }p& {\displaystyle }={\displaystyle \frac{\mathrm{1}-e^B}{e^A-e^B}},\\ {\displaystyle }q& {\displaystyle }=\mathrm{1}-p,\end{array}$$

and

$$\begin{array}{ll}{\displaystyle }A& {\displaystyle }=-{\displaystyle \frac{\mathrm{1}}{\mathrm{2}\mathit{ϵ}}}+{\displaystyle \frac{\mathit{\pi }}{\mathit{\delta }}}\sqrt{\mathrm{1}+{\displaystyle \frac{{\mathit{\delta }}^{\mathrm{2}}}{\mathrm{4}{\mathit{\pi }}^{\mathrm{2}}{\mathit{ϵ}}^{\mathrm{2}}}}},\\ {\displaystyle }B& {\displaystyle }=-{\displaystyle \frac{\mathrm{1}}{\mathrm{2}\mathit{ϵ}}}-{\displaystyle \frac{\mathit{\pi }}{\mathit{\delta }}}\sqrt{\mathrm{1}+{\displaystyle \frac{{\mathit{\delta }}^{\mathrm{2}}}{\mathrm{4}{\mathit{\pi }}^{\mathrm{2}}{\mathit{ϵ}}^{\mathrm{2}}}}}.\end{array}$$

As is evident from Eq. (4), the spatial structure of the stream function is controlled by both ϵ and δ. Figure 1a and c depict the stream functions for two ϵ regimes of S48's model: (i) weak damping (ϵ≤δ²) and (ii) strong damping (ϵ>δ²). For ϵ≤δ², the solution ψ given by Eq. (4) becomes linear in x and thus can satisfy only one boundary condition out of two. This solution is commonly assumed to approximate the exact solution for ψ in the frictionless interior of the basin while a different approximation applies in the narrow, frictional boundary layer adjacent to x=0. Figure 1a depicts this narrow boundary layer for ϵ=0.1δ² where the stream function first decreases fast with x for small x and then increases slowly with x for large x. For ϵ>δ², the solution, ψ, is symmetric about $x=\frac{\mathrm{1}}{\mathrm{2}}$ and can satisfy the two boundary conditions, $\mathit{\psi }(\mathrm{0},y)=\mathrm{0}=\mathit{\psi }(\mathrm{1},y)$. This is demonstrated in the symmetric stream function depicted in Fig. 1c for ϵ=10δ². The explicit expressions of ψ in the two ranges of ϵ are given in Appendix B.

Figure 1 The stream functions in different ϵ regimes of S48's and M50's models for $\mathit{\delta }=\mathrm{2}\mathit{\pi }/\mathrm{10}$: (a) and (b) weak damping (ϵ≤δ² in S48's model and $\mathit{ϵ}\le \mathrm{0.1}{\mathit{\delta }}^{\mathrm{4}/\mathrm{3}}$ in M50's model) – there exists a narrow fast flowing current along the western edge of the basin; (c) and (d) strong damping (ϵ>δ² in S48's model and $\mathit{ϵ}>\mathrm{0.1}{\mathit{\delta }}^{\mathrm{4}/\mathrm{3}}$ in M50's model) – the stream function is (nearly) symmetric about x=0.5, which indicates that there is no westward intensification.

In S48's model, we define the transport of the WBC as the product of its width, ϵ, and the average of the meridional velocity, $v=-\mathit{\delta }\frac{\partial \mathit{\psi }}{\partial x}$, between the western edge of the basin, x=0, and x=ϵ evaluated along $y=\frac{\mathrm{1}}{\mathrm{2}}$, i.e., $\mathrm{Tr}=\mathit{ϵ}\left(\frac{\mathrm{1}}{\mathit{ϵ}}{\int }_{\mathrm{0}}^{\mathit{ϵ}}-\mathit{\delta }{\left.\frac{\partial \mathit{\psi }}{\partial x}\right|}_{y=\frac{\mathrm{1}}{\mathrm{2}}}\mathrm{d}x\right)$. The integral in the definition of the transport, Tr, simplifies to the product of the domain aspect ratio and the difference in the values of the stream function evaluated at x=0 and x=ϵ along $y=\frac{\mathrm{1}}{\mathrm{2}}$, i.e., Tr $=\mathit{\delta }\left[\mathit{\psi }\left(\mathrm{0},\frac{\mathrm{1}}{\mathrm{2}}\right)-\mathit{\psi }\left(\mathit{ϵ},\frac{\mathrm{1}}{\mathrm{2}}\right)\right]$. Substituting the boundary condition $\mathit{\psi }\left(\mathrm{0},\frac{\mathrm{1}}{\mathrm{2}}\right)=\mathrm{0}$ and using the explicit solution Eq. (4) yields

$$\begin{array}{}\text{(5)}& \mathrm{Tr}={\displaystyle \frac{{\mathit{\delta }}^{\mathrm{3}}}{\mathit{ϵ}{\mathit{\pi }}^{\mathrm{2}}}}(\mathrm{1}-p{e}^{A\mathit{ϵ}}-q{e}^{B\mathit{ϵ}}).\end{array}$$

This expression will be compared below to its counterpart in M50's model and will be compared in Sect. 3 with transports calculated by numerical simulations.

Here, we note that the definition of the WBC's width is somewhat arbitrary, and for definiteness we choose it to be ϵ (as in Welander, 1976; Bye and Veronis, 1979; Vallis, 2017). However, the conclusions drawn in this study are independent of the precise definition; for instance, the width of the WBC can also be defined as the value of x at which the stream function reaches an extremum. According to this definition the WBC's width ϵ^′ equals ∼5ϵ, and the corresponding transport is given by ${\mathrm{Tr}}^{\prime }=\frac{{\mathit{\delta }}^{\mathrm{3}}}{\mathrm{5}\mathit{ϵ}{\mathit{\pi }}^{\mathrm{2}}}(\mathrm{1}-p{e}^{\mathrm{5}A\mathit{ϵ}}-q{e}^{\mathrm{5}B\mathit{ϵ}})$. Both expressions of Tr and Tr^′ yield that the transport of the WBC in S48's model varies as $\sim \frac{{\mathit{\delta }}^{\mathrm{3}}}{\mathit{ϵ}}$.

2.2 M50's non-dimensional counterpart

The non-dimensional counterpart of M50's vorticity equation, obtained by employing the scaling proposed in this study in a similar manner to that of S48 (refer to Munk, 1950 for the dimensional equation), is given by

$$\begin{array}{}\text{(6)}& -{\displaystyle \frac{{\mathit{ϵ}}^{\mathrm{3}}}{{\mathit{\delta }}^{\mathrm{4}}}}{\mathrm{\nabla }}^{\mathrm{4}}\mathit{\psi }+{\displaystyle \frac{\partial \mathit{\psi }}{\partial x}}=\sin \left(\mathit{\pi }y\right),\end{array}$$

where

$$\begin{array}{}\text{(7)}& \mathit{ϵ}={\displaystyle \frac{\mathrm{1}}{L_x}}{\left({\displaystyle \frac{\mathit{\mu }}{\mathit{\beta }}}\right)}^{\mathrm{1}/\mathrm{3}},{\mathrm{\nabla }}^{\mathrm{4}}={\mathit{\delta }}^{\mathrm{4}}{\displaystyle \frac{{\partial }^{\mathrm{4}}}{\partial x^{\mathrm{4}}}}+\mathrm{2}{\mathit{\delta }}^{\mathrm{2}}{\displaystyle \frac{{\partial }^{\mathrm{4}}}{\partial x^{\mathrm{2}}\partial y^{\mathrm{2}}}}+{\displaystyle \frac{{\partial }^{\mathrm{4}}}{\partial y^{\mathrm{4}}}},\end{array}$$

where μ is the (dimensional) horizontal eddy viscosity coefficient. We note that contrary to Eq. (2), the sign in front of the first term in Eq. (6) is negative. This dissimilarity arises because, unlike the parametrization in S48, in M50's model the damping is parametrized by the two-dimensional bi-Laplacian operator. Also, in addition to stream function vanishing at the edges of the basin another set of boundary conditions has to be specified to solve the fourth-order Eq. (6). The additional boundary conditions employed by M50 originate from the inclusion of lateral viscosity which implies that there should be no tangential flow at the basin's edges, i.e., ${\left.\frac{\partial \mathit{\psi }}{\partial x}\right|}_{x=\mathrm{0},\mathrm{1}}={\left.\frac{\partial \mathit{\psi }}{\partial y}\right|}_{y=\mathrm{0},\mathrm{1}}=\mathrm{0}$. Following the mathematical steps in M50 yields the following approximate solution of Eq. (6):

$$\begin{array}{}\text{(8)}& \mathit{\psi }=-\sin \left(\mathit{\pi }y\right)\left[\mathrm{1}-x+\mathit{ϵ}{e}^{(x-\mathrm{1})/\mathit{ϵ}}-e^{-(x/\mathrm{2}\mathit{ϵ})}\mathit{\xi }\left(\mathit{ϵ}\right)\right]\end{array}$$

and

$$\mathit{\xi }\left(\mathit{ϵ}\right)=\left[\cos \left({\displaystyle \frac{\sqrt{\mathrm{3}}x}{\mathrm{2}\mathit{ϵ}}}\right)+{\displaystyle \frac{\mathrm{1}-\mathrm{2}\mathit{ϵ}}{\sqrt{\mathrm{3}}}}\sin \left({\displaystyle \frac{\sqrt{\mathrm{3}}x}{\mathrm{2}\mathit{ϵ}}}\right)\right].$$

Figure 2 Numerically obtained, non-dimensional stream functions for ϵ=0.01 and L_x=10000 km. (a) S48's model with $\mathit{\delta }=\mathrm{2}\mathit{\pi }/\mathrm{10}$, (b) M50's model with $\mathit{\delta }=\mathrm{2}\mathit{\pi }/\mathrm{10}$. Panels (c) and (d) are the same as (a) and (b) but for $\mathit{\delta }=\mathrm{0.25}\mathit{\pi }/\mathrm{10}$; i.e., the meridional extent of the basin is one-eighth of that in (a) and (b). Note the different color bars in panels (a) and (c).

Figure 1b and d depict the stream function for small and large damping in M50's model. For large damping the stream function shown in Fig. 1d is not entirely symmetric about $x=\frac{\mathrm{1}}{\mathrm{2}}$. Also, unlike the behavior of the stream function in S48's model, the stream function in M50's model skews more towards the eastern boundary with the increase in damping. This less-than-optimal behavior of the stream function in M50's model occurs because the stream function does not vanish identically along the eastern boundary and is, instead, a function of ϵ itself (although, for small ϵ, the zonal velocity there is small compared to the rest of the basin).

We turn now to the estimation of the WBC's transport in M50's model. As was done in S48's model, this transport is also defined as the product of the boundary layer width (ϵ) and the mean meridional velocity of the current between x=0 and x=ϵ along $y=\frac{\mathrm{1}}{\mathrm{2}}$. Following the arguments laid out in the previous section (see the paragraph above Eq. 5) an expression for transport can be obtained by multiplying the domain aspect ratio by the difference of the stream function values between x=0 and x=ϵ along $y=\frac{\mathrm{1}}{\mathrm{2}}$. Furthermore, substituting the boundary condition $\mathit{\psi }\left(\mathrm{0},\frac{\mathrm{1}}{\mathrm{2}}\right)=\mathrm{0}$ yields $\mathrm{Tr}=-\mathit{\delta }\mathit{\psi }\left(\mathit{ϵ},\frac{\mathrm{1}}{\mathrm{2}}\right)$. Evaluating ψ in Eq. (8) at $\left(\mathit{ϵ},\frac{\mathrm{1}}{\mathrm{2}}\right)$ for ϵ≪1 yields the following simplified expression for the WBC's transport in M50's model:

$$\begin{array}{}\text{(9)}& \mathrm{Tr}=\mathit{\delta }\left(\mathrm{1}-e^{(-\mathrm{1}/\mathrm{2})}\left[\cos \left({\displaystyle \frac{\sqrt{\mathrm{3}}}{\mathrm{2}}}\right)+{\displaystyle \frac{\mathrm{1}-\mathrm{2}\mathit{ϵ}}{\sqrt{\mathrm{3}}}}\sin \left({\displaystyle \frac{\sqrt{\mathrm{3}}}{\mathrm{2}}}\right)\right]\right).\end{array}$$

As anticipated by Bye and Veronis (1979), the transport of the WBC in S48's model (given by Eq. 5) is governed by both damping (ϵ) and the domain aspect ratio (δ). However, in M50's model the dependence of the WBC's transport on the two parameters is strikingly different: the transport is governed primarily by δ and is weakly dependent on ϵ. In the next section we validate these claims using (dimensional) numerical simulations and then apply our results to the present-day world ocean.

3 Numerical simulations and application to the world ocean

The numerical simulations described below were carried out using the time-dependent, forced-dissipative, rotating shallow water equation (SWE) dimensional solver that was successfully used in previous studies. The solver employs the finite-difference method to solve SWEs on the β plane, and the simulations are carried out on an Arakawa C grid with leapfrog time difference scheme. Though the solver can include nonlinear terms, these terms were neglected in the present application. The reader should refer to Gildor et al. (2016) and Shamir et al. (2019) for a more detailed description of the solver.

The simulations presented here were carried out in a barotropic ocean with the same characteristics as in S48, i.e., on an equatorial β plane (f₀=0), forced by a wind stress that varies as $-{\mathit{\tau }}_{\mathrm{0}}\cos \left(\frac{\mathit{\pi }{y}^{*}}{L_y}\right)$. Three of the dimensional parameters remained fixed in all the simulations presented below – the gradient of the Coriolis frequency (given by $\mathit{\beta }=\mathrm{2}\times {\mathrm{10}}^{-\mathrm{11}}$ ${\mathrm{m}}^{-\mathrm{1}}{\mathrm{s}}^{-\mathrm{1}}$), the zonal extent of the basin (L_x = 10 000 km) and the amplitude of the prescribed forcing (τ₀=0.2 N m⁻²). The other two-dimensional parameters in the two WBC models, i.e., the damping coefficients (Rayleigh friction coefficient (r) in S48's model and horizontal eddy viscosity (μ) in M50's model) and the meridional extent of the basin (L_y), are varied to examine the effect of ϵ and δ on the transport. We note that keeping τ₀ fixed and varying L_y will yield different values of γ in the simulation; however, since we scale our ψ^* on $\mathit{\gamma }\mathit{\beta }{L}_y^{\mathrm{3}}$ and only look at the non-dimensional transport we do not have to account for the effects of changes in γ. The results are consistent with what one would obtain by keeping only β and L_x fixed and varying τ₀ along with the damping coefficients and L_y to keep γ constant. The boundary conditions are that the (dimensional) zonal and meridional velocities vanish along the basin's meridional and zonal boundaries respectively, i.e., $u^{*}{|}_{y^{*}=\mathrm{0},L_y}=v^{*}{|}_{x^{*}=\mathrm{0},L_x}=\mathrm{0}$. The numerical solver is integrated until a steady state is reached. The steady state of the time-dependent simulations is defined as the state at which the dependent variables in the SWEs (dimensional zonal velocity (u^*), meridional velocity (v^*) and sea surface height (η^*)) cease to evolve for sufficiently long time.

Figure 2a and c depict the numerically obtained, non-dimensional stream function $\left(\mathit{\psi }=\frac{{\mathit{\psi }}^{*}}{\mathit{\gamma }\mathit{\beta }{L}_y^{\mathrm{3}}}\right)$ in the steady state for the dimensional parameters as in S48's model (and the corresponding value of δ is noted above these panels) while Fig. 2b and d depict the numerically obtained, non-dimensional ψ in the steady state for the parameters relevant to M50's model (and here too the corresponding value of δ is noted above these panels). The WBC's width is kept constant in all four panels and is noted in the caption. The reader should note that the meridional extent (L_y) of the basin shown in Fig. 2a and b is 2π×1000 km, whereas, the meridional extent of the basin shown in Fig. 2c and d is $L_y=\mathit{\pi }/\mathrm{4}\times \mathrm{1000}$ km. In all four cases the shape of the stream function is very similar to the steady non-dimensional stream functions shown in Fig. 1a and b. We note that for the given values of (ϵ, δ), the values of ψ obtained by dividing the numerically calculated values of ψ^* by the corresponding values of $\mathit{\gamma }\mathit{\beta }{L}_y^{\mathrm{3}}$ are in agreement with the ψ calculated analytically for the same values of (ϵ, δ) using Eq. (4) for S48's model and Eq. (8) for M50's model.

Figure 2 shows that in both S48's and M50's models, for a fixed value of ϵ=0.01 (damping and the width of the WBC), the gradient of the stream function increases with δ. The higher zonal gradient of the stream function near the western boundary yields a larger meridional velocity, thus increasing the transport of the WBC (given by the product of width and average meridional velocity of the WBC). Clearly, δ exercises control over the transport of the WBC and hence cannot be ignored.

Figure 3 Comparison between analytically (solid lines) and numerically (dots) calculated values of transport (Tr) as a function of ϵ for different values of δ in (a) S48's and (b) M50's models. The dashed lines depict the cubic-spline-interpolated curves between the numerically calculated transports (dots).

Figure 3 compares the analytic and numerically computed values of the non-dimensional transport (Tr) of the WBC in S48's and M50's models as a function of ϵ for several values of δ. The solid lines denote the analytic value of Tr obtained from the expressions given by Eqs. (5) and (9). The “numerical transport” of the WBC is obtained by taking the product of δ and $-\frac{{\mathit{\psi }}^{*}(\mathit{ϵ},\frac{\mathrm{1}}{\mathrm{2}})}{\mathit{\gamma }\mathit{\beta }{L}_y^{\mathrm{3}}}$. Here, ${\mathit{\psi }}^{*}\left(\mathit{ϵ},\frac{\mathrm{1}}{\mathrm{2}}\right)$ is the value of the steady-state dimensional stream function at $\left(\mathit{ϵ},\frac{\mathrm{1}}{\mathrm{2}}\right)$ obtained from the results of the numerical simulation for a given set of parameters which correspond to a certain (ϵ,δ).

As is evident from Fig. 3a, the analytic and numerically calculated non-dimensional transports of the WBC in S48's model are in good agreement. Figure 3b shows that the analytic transports of the WBC in M50's model are nearly independent of ϵ and are governed primarily by δ. The numerically calculated transports of the WBC depicted by the dashed lines in Fig. 3b show a similar dependence on δ, but in contrast to the approximate analytic expression these transports vary slightly with ϵ. We also note that there is a discernible difference between the analytically estimated and numerically calculated values of transport in M50's model for nearly all values of (ϵ, δ). This is because the expression for the stream function in M50's model, (i.e., Eq. 8) only crudely approximates the actual stream function.

Figure 4 The non-dimensional transport of the western boundary current (WBC) as a function of ϵ and δ in (a) S48's model and (b) M50's model. Contours depicted are obtained from numerical simulations for fixed L_x=10 000 km, $\mathit{\beta }=\mathrm{2}\times {\mathrm{10}}^{-\mathrm{11}}$ ${\mathrm{m}}^{-\mathrm{1}}{\mathrm{s}}^{-\mathrm{1}}$, and by varying r between $\mathrm{1}/\mathrm{100}$ and $\mathrm{1}/\mathrm{10}$ d⁻¹ in S48's model (a) and μ=10 000 m² s⁻¹ in M50's model (b). The different WBCs in the world ocean are depicted with different symbols, and the error bars denote the possible variability of parameters that can occur because of an error in estimating the zonal and meridional extents of the basins that contain the WBC. The width of each of the WBCs, ϵ, is estimated using $r=\mathrm{1}/\mathrm{30}$ d⁻¹ (a) and μ=10 000 m² s⁻¹ (b) and $\mathit{\beta }=\mathrm{2}\times {\mathrm{10}}^{-\mathrm{11}}$ ${\mathrm{m}}^{-\mathrm{1}}{\mathrm{s}}^{-\mathrm{1}}$. The zonal extent, L_x, of the basin containing the particular WBC is given in Table C1. The error in ϵ is not accounted for in (b) because the WBC's transport in M50's model is (nearly) independent of ϵ. The East Australian Current's non-dimensional transport, as calculated from both S48's and M50's models, is less than those of the other four WBCs. The uncertainty in δ for the Brazil Current extends up to 0.375 in both models; the contours have been restricted to better resolve the boundary currents. The range between which the non-dimensional transport varies is similar in both the models.

Figure 4 depicts the non-dimensional transport of the WBC in S48's (a) and M50's models (b) as contours on the (ϵ,δ) plane. The contours were obtained by interpolating (using the cubic spline method) between the numerically calculated values of the WBC's transport as shown in Fig. 3 (plus a few additional values at low ϵ in S48's model). As is evident from Fig. 4a, the non-dimensional transport is a function of both ϵ and δ in S48's model. On the other hand, Fig. 4b shows that the transport of the WBC is only weakly dependent on ϵ and is governed primarily by δ (the contours are nearly parallel to the abscissa). The position of the different WBCs in the (ϵ,δ) parameter space is marked with different symbols, and the error bars account for the inaccuracies in the assigned values of the zonal and meridional extents of the basins. The details of how the irregular basins in the world ocean are approximated with rectangles are discussed in Appendix C. The error in ϵ has been omitted from Fig. 4b because the WBC's transport in M50's model is nearly independent of ϵ. Despite the large uncertainty in the damping parameters (relevant in S48's model only) and domain aspect ratios (relevant in both S48's and M50's models) of the five WBCs, the non-dimensional transport of the East Australian Current is distinctly smaller than that of the other WBCs.

4 Summary and discussions

Since the introduction of S48's and M50's models about 70 years ago, numerous theoretical and numerical investigations have been carried out to further explore the characteristics of westward intensification (Munk and Carrier, 1950; Stommel, 1958; Hogg and Johns, 1995; Pedlosky, 1987; Vallis, 2017, and references therein). Both S48's and M50's dimensional models clearly bring out the contribution of each source of vorticity, damping, planetary gradient and wind forcing, in producing the characteristic east–west asymmetry of the flow in a basin. However, it is difficult to quantify the contribution of each of the five dimensional parameters (L_x, L_y, β, τ₀ and $r/\mathit{\mu }$) to the transport of WBC, using the dimensional models. A better alternative is to combine several dimensional parameters to yield a system with fewer non-dimensional parameters as was employed by, for example, Welander (1976) to identify a zonally uniform regime in ocean circulation and by Bye and Veronis (1979) to identify the correction to the Sverdrup transport in the context of S48's original model.

In this article, we address the issue raised by Bye and Veronis (1979) regarding the effect of domain aspect ratio on the WBC's transport by providing explicit expressions of the non-dimensional transport in both S48's and M50's models. These expressions are then benchmarked against numerical simulations of the time-dependent, forced-dissipative, rotating shallow water equations. Both the analytic expressions and steady-state simulations show that the WBCs' transports depend on both ϵ and δ in S48's model; however, a change in δ has a stronger effect on Tr when compared to a change in ϵ $\left(\mathrm{Tr}\sim \frac{{\mathit{\delta }}^{\mathrm{3}}}{\mathit{ϵ}}\right)$. In contrast, the transport of the WBC in M50's model is nearly independent of ϵ and is governed primarily by δ (Tr∼δ).

In the traditional description of the S48 model the flow is decomposed into two parts: a slow, anti-cyclonic flow in the inner basin where the velocities are tiny so frictional effects can be neglected, and a return boundary flow where the frictional vorticity associated with the zonal shear of the poleward directed velocity balances the planetary vorticity advected by this velocity. According to this paradigm the WBC simply returns the frictionless equatorward Sverdrup transport of the inner basin so its transport is independent of the friction coefficient, and since the (dissipation) Laplacian term does not affect the Sverdrup interior flow, the transport of the WBC should also be independent of the domain aspect ratio. The present study demonstrates that the assumption of small damping, ϵ≪1, implies that only the term $\mathit{ϵ}\frac{{\partial }^{\mathrm{2}}\mathit{\psi }}{\partial x^{\mathrm{2}}}$ of the Laplacian in Eq. (2) can be neglected in this limit while the second term, $\frac{\mathit{ϵ}}{{\mathit{\delta }}^{\mathrm{2}}}\frac{{\partial }^{\mathrm{2}}\mathit{\psi }}{\partial y^{\mathrm{2}}}$, cannot be neglected in the interior solution when δ²∼O(ϵ). The implication of our analysis is that the Sverdrup interior flow depends on δ for sufficiently small δ and therefore so does the (return) transport of the WBC.

To appreciate this subtle issue one should compare a square basin, where δ=1, with a narrow and long “channel-like” basin where δ≪1. In a square basin, the classical approach of equating $\frac{\partial \mathit{\psi }}{\partial x}$ to sin (πy) in the inner basin works well since the north–south gradient of the zonal velocity (represented by $\frac{{\partial }^{\mathrm{2}}\mathit{\psi }}{\partial y^{\mathrm{2}}}$) is small and can be neglected from the interior solution. However, in a “channel-like” ocean this quantity is large and cannot be neglected from the balance of terms in the interior solution. An examination of the three vorticity terms in the interior [$\frac{\partial \mathit{\psi }}{\partial x}$, $\frac{{\partial }^{\mathrm{2}}\mathit{\psi }}{\partial y^{\mathrm{2}}}$ and −sin (πy)] clarifies that the transport of the WBC (as well as the equatorward transport in the interior) in a “channel-like” ocean should be smaller compared to a square ocean since the meridional shear of the zonal velocity lowers the vorticity induced by the curl of the wind stress. In an extreme “channel-like” ocean with $\frac{\mathit{ϵ}}{{\mathit{\delta }}^{\mathrm{2}}}\gg \mathrm{1}$ the only term that can balance $\frac{{\partial }^{\mathrm{2}}\mathit{\psi }}{\partial y^{\mathrm{2}}}$ is $\frac{\partial \mathit{\psi }}{\partial x}$ that implies a strong, equatorward directed velocity. Indeed, as was shown by Welander (1976) for a small domain aspect ratio, a boundary layer develops along the basin's eastern boundary in which the strong current flows equatorward.

In M50's model the vorticity balance of the interior is more involved since the bi-Laplacian dissipation operator (∇⁴) has three terms, each of which has a coefficient of different power of δ. Thus, the distinction between terms associated with the inner basin and those with the boundary solution is not as clear as in S48's model. However, under the assumption of small damping, ϵ≪1, the third term in the ∇⁴ operator (given by Eq. 7) cannot be neglected for δ⁴∼O(ϵ³). Similar to S48's model, the vorticity balance in the interior, which is determined in this limit by the interplay of three terms, $\frac{\partial \mathit{\psi }}{\partial x}$, $\frac{{\partial }^{\mathrm{4}}\mathit{\psi }}{\partial y^{\mathrm{4}}}$ and −sin (πy), yields a δ-dependent equatorward transport. This δ-dependent transport in the interior is balanced by an equal, poleward-directed transport along the western boundary. Although the vorticity associated with $\frac{{\partial }^{\mathrm{4}}\mathit{\psi }}{\partial y^{\mathrm{4}}}$ is not as intuitive as that associated with $\frac{{\partial }^{\mathrm{2}}\mathit{\psi }}{\partial y^{\mathrm{2}}}$ in S48's model, the change it entails in Sverdrup's interior solution is similar.

The results derived here highlight an important effect that was overlooked in the classical (traditional) WBC theory, namely, the effect of the domain aspect ratio on the Sverdrup solution of the inner basin which results from the meridional shear of the zonal velocity in a narrow zonal channel.

The non-dimensional formulation presented here does not alter the physical basis of the S48 and M50 models. We emphasize that the dimensional transport (calculated from the product of the non-dimensional transport and $\mathit{\gamma }\mathit{\beta }{L}_y^{\mathrm{3}}{H}_{\mathrm{0}}$) in S48's model varies linearly with the Rayleigh friction coefficient (r), while in M50's model it is nearly independent of the eddy viscosity (μ). In both models the transport is linear with the magnitude of the wind-stress curl.

The application of our results to the present-day ocean attributes the small transport of the EAC compared to the other WBCs to the geometry of the South Pacific ocean. In reality, factors other than the domain aspect ratio may also be important in determining the transport. For instance, the Brazil Current's volumetric transport is low (especially in the northern part) because the current is largely confined to the continental shelf (Stramma et al., 1990). Temperature-driven buoyancy fluxes can also affect the transport of a WBC (Hogg and Gayen, 2020).

It is highly plausible that with a different arrangement of the continents in previous geologic times, the small domain aspect ratio that persisted in the ocean in those times was unable to support a strong WBC. Thus, the higher pole-to-Equator temperature gradient might have resulted from the smaller poleward heat flux by the WBC due to the decrease in the domain aspect ratio. This hypothesis should be addressed in a future work.

Appendix A: Typos in Stommel (1948)

There are some typos in the expression for u (Eq. 21) and η (Eq. 23) in Stommel (1948). The correct expressions are given as follows:

$$\begin{array}{}\text{(A1)}& {\displaystyle }u=\mathit{\gamma }\left(b/\mathit{\pi }\right)\cos \left(\mathit{\pi }y/b\right)\left(p{e}^{Ax}+q{e}^{Bx}-\mathrm{1}\right),\text{(A2)}& {\displaystyle }\begin{array}{rl}\mathit{\eta }& =-(F/gD)\cos (\mathit{\pi }y/b)(e^{Ax}p/A+e^{Bx}q/B)\\ & -(f\mathit{\gamma }/g)(b/\mathit{\pi }{)}^{\mathrm{2}}\sin (\mathit{\pi }y/b\left)\right(p{e}^{Ax}+qB{e}^{Bx}-\mathrm{1})\\ & +(\partial f/\partial y)(\mathit{\gamma }/g)(b/\mathit{\pi }{)}^{\mathrm{3}}\cos (\mathit{\pi }y/b).\end{array}\end{array}$$

For the reader's perusal, the variables in the aforementioned equations are the same as the ones defined in Stommel (1948). Figure A1 provides excerpts from Stommel (1948), over which the corrections have been highlighted.

Figure A1 Corrections to u and h indicated over excerpts from Stommel (1948).

Appendix B: Limiting cases of stream function ψ in S48's model

In the limit ϵ≤δ², the solution ψ tends to

$$\begin{array}{}\text{(B1)}& {lim}_{\mathit{ϵ}\le {\mathit{\delta }}^{\mathrm{2}}}\mathit{\psi }(x,y)=c_{\mathrm{1}}{\displaystyle \frac{{\mathit{\delta }}^{\mathrm{2}}}{\mathit{ϵ}{\mathit{\pi }}^{\mathrm{2}}}}\sin \left(\mathit{\pi }y\right)\left(x\right),\end{array}$$

where $c_{\mathrm{1}}={lim}_{\mathit{ϵ}\le {\mathit{\delta }}^{\mathrm{2}}}A$ is a number ≪1. On the other hand, in the limit of ϵ>δ², the solution ψ becomes

$$\begin{array}{}\text{(B2)}& {lim}_{\mathit{ϵ}>{\mathit{\delta }}^{\mathrm{2}}}\mathit{\psi }(x,y)={\displaystyle \frac{{\mathit{\delta }}^{\mathrm{2}}}{\mathit{ϵ}{\mathit{\pi }}^{\mathrm{2}}}}\sin \left(\mathit{\pi }y\right)\left[p\left(e^{c_{\mathrm{2}}x}+e^{c_{\mathrm{2}}(\mathrm{1}-x)}\right)-\mathrm{1}\right],\end{array}$$

where $c_{\mathrm{2}}={lim}_{\mathit{ϵ}>{\mathit{\delta }}^{\mathrm{2}}}A=\frac{\mathit{\pi }}{\mathit{\delta }}$ and $p=\frac{e^{c_{\mathrm{2}}}-\mathrm{1}}{e^{\mathrm{2}{c}_{\mathrm{2}}}-\mathrm{1}}$. The function ${lim}_{\mathit{ϵ}>{\mathit{\delta }}^{\mathrm{2}}}\mathit{\psi }(x,y)$ is symmetric about $x=\frac{\mathrm{1}}{\mathrm{2}}$.

Appendix C: Zonal and meridional extents of the five western boundary currents in present-day world ocean

To determine the zonal and meridional extents of a basin containing a WBC, we identified the mean initiation and termination latitudes of each WBC based on the available literature. The Gulf Stream begins at the tip of Florida (∼ 25^∘ N) and runs up to ∼ 38^∘ N, where it breaks of into hot and cold rings (Hogg and Johns, 1995). The Kuroshio originates from the bifurcation of the North Equatorial Current at 12–13^∘ N, although this bifurcation point can vary between 10–15^∘ N (Qiu and Lukas, 1996); it separates from the Japan coast at 35^∘ N as a meandering current colloquially known as the Kuroshio extension which stretches as far as ∼ 38^∘ N (Kida et al., 2016). The East Madagascar–Agulhas current, in the South Indian ocean, runs from 20 to 40^∘ S (Lutjeharms et al., 1981; Gordon, 1985; Lutjeharms, 2006) – however, the current retroflects between 38 to 40^∘ S (Quartly and Srokosz, 1993). Moreover, the African continental landmass ends close to 35^∘ S. The Brazil Current begins between 10 and 12^∘ S (Peterson and Stramma, 1991; Stramma et al., 1990), but the intense current attains the intense speed characteristic of a WBC only when it crosses the Vitoria–Trindade Ridge at 20.5^∘ S (Evans et al., 1983). This current separates from the coastline at a mean value of 36^∘ S ± 1.1^∘ (Olson et al., 1988). The last of the five WBCs in the world ocean is the East Australian Current (EAC) that extends from 18^∘ S to around 35^∘ S (Boland and Church, 1981; Ridgway and Godfrey, 1994), but a characteristic southward flow is evident only when EAC crosses 22^∘ S (Ridgway and Dunn, 2003); the current usually separates from the coast at 33^∘ S (Archer et al., 2017).

Table C1 Dimensions of the gyres that contain the five western boundary currents in the present-day world ocean.

We define the meridional extent (L_y) as the distance between the initiation and termination latitudes of the WBC. On the other hand, to determine the zonal extent (L_x) we calculate the distances between the land masses at both the initiation latitude and termination latitude. The average of the two distances is defined as the typical L_x for a given WBC. For instance, the approximate initiation and termination coordinates for the Kuroshio are 13^∘ N, 125^∘ E and 35^∘ N, 140^∘ E respectively, which yields L_y≈2500 km. The distances to the opposite landmass, the North American continent (which forms the eastern boundary of the basin), as calculated from the initiation coordinate and termination coordinate are ∼ 9000 km and ∼ 15 000 km respectively. Thus, the typical zonal extent of the basin is assumed to be 12 000 km.

The mean dimensions L_x and L_y for all the five WBCs in the world ocean are given by Table 1. The “error” in L_y accounts for the variation between different references of the initiation and termination latitudes, and the error in L_x is the deviation of the measured zonal distances along initiation and termination latitude from the mean value. Based on these values of L_x and L_y a range of parameter damping (ϵ) and domain aspect ratios (δ) corresponding to every WBC was estimated, and these values of ϵ and δ were employed to distinguish between the five WBCs. Typical values of L_x and L_y for the ocean basins that contain the WBCs were also estimated using the mean streamlines in the ocean as calculated by Maximenko et al. (2009) – these values were well within the range cited in Table C1.