Shallow katabatic flow on a non-uniformly cooled slope

We examine katabatic flow driven by a non-uniformly cooled slope surface but unaffected by Coriolis acceleration. A general formulation is given, valid for non-uniform surface buoyancy distributions over a down-slope length scale L≫δ0\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$L\gg \delta _0$$\end{document}, where δ0=ν/(Nsinα)1/2\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\delta _0=\nu /(N\sin \alpha )^{1/2}$$\end{document} is the slope-normal Prandtl depth, for a kinematic viscosity ν\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\nu$$\end{document}, buoyancy frequency N and slope angle α\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\alpha$$\end{document}. We demonstrate that the similarity solution of Shapiro and Fedorovich (J Fluid Mech 571:149–175, 2007) can remain quantitatively relevant local to the end of a non-uniformly cooled region. The usefulness of the steady similarity solution is determined by a spatial eigenvalue problem on the L length scale. Broadly speaking, there are also two modes of temporal instability; stationary down-slope aligned vortices and down-slope propagating waves. By considering the limiting inviscid stability problem, we show that the origin of the vortex mode is spatial oscillation of the buoyancy profile normal to the slope. This leads to vortex growth in a region displaced from the slope surface, at a point of buoyancy inflection, just as the propagating modes owe their existence to an inflectional velocity. Non-uniform katabatic flows that detrain fluid to the ambient are shown to further destabilise the vortex mode whereas entraining flows lead to weaker vortex growth rates. Rayleigh waves dominate in general, but the vortex modes become more significant at small slope angles and we quantify their relative growth rates.


Introduction
We consider a planar sloping boundary, immersed in a fluid that is density stratified. If the slope is cold relative to the stratified ambient fluid at the same horizontal level, then a down-slope 'katabatic' flow of denser slope-adjacent fluid is to be expected. Such downslope flows are found commonly in applications with cooling of sloping terrain in a stratified atmosphere, see for example the review [1]. Similar weak flows can also be obtained for a stably salt-stratified fluid adjacent to a sloping boundary, driven by a zero-flux condition at the boundary, as discussed by Philips [2] and Wunsch [3].
If a constant potential temperature difference between the slope and the ambient fluid is maintained at all heights, then the buoyancy at the slope surface is spatially uniform and a solution exists due to Prandtl [4]. This solution will be explicitly discussed below and is given by Eq. (6); this is a parallel flow, with no motion in the slope-normal direction. In terms of the velocity field, the Prandtl solution shares some of the characteristics of a wall jet or gravity current, with a maximum downwards flow near the slope surface but the down-slope velocity profile is bi-directional. There is a maximum downwards flow near the slope surface with a return flow predicted at larger distances from the surface. The downslope velocity profile and the buoyancy profile both oscillate and decay in the slope-normal coordinate, as illustrated in Fig. 1.
Later extensions to the Prandtl state included the effects of a height dependent diffusivity profile; see for example the formulation of Rao and Snodgrass [5] and later discussion of Grisogono and Oerlemans [6]. However, observations made in otherwise calm conditions remain reasonably close to both the Prandtl prediction and the slower decaying profiles of Rao and Snodgrass [5] but with some scatter particularly in the far field as shown in Fig. 3 of [7]. Depth-integrated models require assumptions to be made for entrainment and profile shapes, but have also been obtained by (for example) Manins and Sawford [8] or Kondo and Sato [7]. are aligned with the downslope and slope-normal directions, with corresponding non-dimensional coordinates (x, z). The ambient fluid is stratified in the vertical Z * direction, with buoyancy frequency N. The buoyancy distribution (red) varies across a layer of width 0 near the surface driving a downslope velocity (black). Both the buoyancy and velocity profile oscillate about zero as they decay into the ambient. Entrainment/detrainment into the layer is dependent on a down slope variation in the surface buoyancy on the L scale Reduced models of this form have been applied with some success to jets/plumes, however simple extensions to unsteady flows are known to be problematic, in some cases leading to illposed models as discussed in [9,10]. It is worth noting that the role of Coriolis forcing in katabatic flow is not a trivial one. It is known, e.g. [11], that even for down slope flows induced by uniform surface buoyancy the effects of background rotation lead to unsteady thickening of the near-surface layer. This behaviour is akin to that observed in Ekman boundary layers on sloping surfaces in a stratified fluid; see for example the 'growing boundary layer' states identified by Duck et al. [12] and Hewitt et al. [13]. In what follows we will assume that Coriolis effects remain negligible throughout.
More recently Xiao and Senocak have described the stability of Prandtl solutions for katabatic [14], anabatic [15] and katabatic flows including ambient forcing [16], assuming a constant surface-buoyancy flux. In general the stability of these states is governed by three parameters associated with slope angle, Prandtl number and a measure of buoyancy forcing at the slope surface. By constructing numerical solutions of the linearised stability equations at finite values of a forcing parameter and at fixed slope angles and Prandtl numbers, critical values of buoyancy forcing are obtained for temporally neutral down-slope propagating waves or downslope aligned vortices. Numerical results close to onset indicated that the propagating wave modes appear first for steep slopes, whilst vortex modes occur first for shallow slopes. We will return to these modes of instability in this work too.
Of some practical relevance are cases of katabatic flow with a non-uniform distribution of surface buoyancy. Such flows must be qualitatively different, leading to a "non-parallel" development (ie. the layer evolves with the down-slope coordinate) and necessitating a corresponding detrainment/entrainment of fluid between the ambient and the near-slope flow, required by mass conservation. The work of Shapiro and Fedorovich [17] considered spatially non-uniform forcing of down-slope flows under the assumption that the surface buoyancy increases/decreases linearly with the down slope coordinate. Their implicit assumption, and an issue we will explicitly address here, was that this simplified flow could find application to slowly varying surface buoyancy profiles or transition zones between regions of uniform surface buoyancy.
In what follows we begin by formulating a general boundary-layer approximation similar to that of Shapiro and Fedorovich [18] but applied to an arbitrary surface buoyancy forcing over a down-slope lengthscale that is substantially larger than the Prandtl depth. The resulting bi-directional parabolic problem for the two-dimensional katabatic flow is solved in Sect. 2 and results are compared to the self-similar predictions made in [17]. Having obtained the two-dimensional base flows we then go on to consider their stability in Sect. 3. By tackling the stability problem in the strongly-driven limit we can pick out the dominant mechanisms associated with the buoyancy-driven vortex perturbations and the inflectional wave perturbations. This approach also allows for a quasi-local treatment, allowing us to determine how the local velocity and buoyancy profiles in the non-uniformly cooled two-dimensional base state affect the stability compared to the uniform case of Xiao and Senocak [14]. A discussion is given in Sect. 4, with particular attention paid to describing the relative importance of the two main instability modes.

Formulation
We begin with the formulation of Shapiro and Fedorovich [17], using a Cartesian (X, Z) coordinate system that is aligned with a slope of angle ( Fig. 1), where X is measured down the slope and Z is normal to the slope; in addition there is a vertical coordinate Z * over which the ambient fluid's density varies. We will denote the velocities in the directions of increasing (X, Z) by (U, W), whilst the pressure is P and buoyancy is B. In terms of the potential temperature , the buoyancy is B = g( − ∞ (Z * ))∕ r , where g is the gravitational acceleration, ∞ (Z * ) is the corresponding ambient value and r is a constant reference value. Similarly we work with a pressure that is measured relative to the ambient hydrostatic pressure, The dimensional equations governing conservation of energy, momentum and mass are where , are a constant kinematic viscosity and thermal diffusivity, ∇ 2 is the two-dimensional Laplacian and subscripts denote differentiation. Here is the buoyancy frequency, measured with respect to the vertical coordinate Z * and is assumed to be constant. The surface buoyancy B(X, Z = 0, T) takes a value B 0 at some reference location, and the spatial variation of B occurs over a known length scale L, down slope of this point. The coordinate system of (1) is then made non-dimensional via We begin with solutions that are independent of the cross slope Y coordinate, but will return to more general behaviour in the context of the stability problem of Sect. 3. The slope-normal length scale, 0 , is the depth of the Prandtl solution [4] when surface buoyancy is spatially uniform.
For the buoyancy, velocities and pressure we introduce where the given U 0 is the appropriate velocity scale of the Prandtl solution assuming an O(1) Prandtl number; this solution is discussed more below and given explicitly by (6).
Here is an aspect ratio formed by the (slope-normal) Prandtl thickness ( 0 ) relative to the down slope lengthscale (L) over which the surface buoyancy is varying. In this choice of non-dimensionalisation we use the absolute value of the reference surface buoyancy |B 0 | , and a negative sign to drive a katabatic flow will be imposed explicitly in the surface boundary condition for B.
The governing dimensionless equations obtained via (1) and (3) become Here Fr is an internal Froude number based on the Prandtl velocity scale U 0 = |B 0 |∕N and the vertical lengthscale L sin , Pr is the corresponding Prandtl number ( ∕ ) and

Boundary and initial conditions
Our focus is initially on steady solutions to (4), where the boundary conditions are These correspond to no-slip, impermeability, a known buoyancy on the slope surface, and no buoyancy perturbation far from the surface. The constraint u sin − w cos → 0 is consistent with no buoyancy perturbation, as can be seen from the form of (4a), and corresponds to a requirement for purely horizontal along-isentrope flow in the far field. The function b s (x) is any general distribution of surface buoyancy. If b s ≡ −1 for all x, then a non-entraining ( w = 0 ) steady solution exists, which is the classical x-independent Prandtl [4] This parallel solution is independent of Fr and any variation with the other parameters is already captured by (3). Whilst (4) can be solved for arbitrary our interest here is restricted to shallow flows for which ≪ 1 , and in this case the condition for u in (5) simplifies to u → 0 in the far field. In the shallow limit the equations are independent of the slope angle; any variation with is captured by the scaling factors of (3), which is consistent with the observations of Shapiro and Fedorovich [18]. Analogues of this ≪ 1 boundary-layer system have been used in other related configurations, for example the Ekman layer on a sloping boundary in a stratified fluid as discussed in [12,19], but in those cases driven by rotation rather than surface buoyancy.

Steady shallow katabatic flow
We now present solutions to the equations (4) subject to conditions (5) for shallow flows with ≪ 1 . Rather than uniform surface conditions associated with the Prandtl solution (6) we allow for more general surface buoyancy distributions. To this end, we consider the piece-wise distribution Here the unit length of the non-uniformly cooled region is set by our choice of non-dimensionalisation of x based on L, and we shall restrict attention to down-slope flows with b s (x) ≤ 0 by choosing ≤ 1.
In the shallow-flow limit ( ≪ 1 ) the governing system (4) becomes parabolic. However the problem is bi-directional with the down slope velocity u having both positive and negative regions. Therefore, the computational approach requires both up-slope and downslope conditions, and in this case we assume that as |x| → ∞ the solution reverts to the appropriate x independent Prandtl state (6) with surface buoyancy values of −1 far up the slope and −1 + far down the slope. Any entrainment/detrainment where w ≠ 0 as z → ∞ therefore occurs only in the vicinity of the transition zone defined by (7), decaying both up ( x → −∞ ) and down ( x → ∞ ) slope.
We obtain steady solutions to the nonlinear two-dimensional problem (4) by a secondorder central finite-difference scheme and Newton iteration for the nonlinear terms. Solving for all x positions simultaneously in this way avoids issues with the bi-directionality, alternatively this can be dealt with by (for example) the zig-zag scheme [20]. For a nonuniformly spaced computational mesh with N x × N z nodes, this leaves each Newton iteration to invert a 3N x N z × 3N x N z sparse linear system which is achieved using the MUMPS library [21]. The results we present are independent of N x , N z and the domain truncation in both x and z.
In Fig. 2 we show computational results for the x-dependent entrainment (negative) and detrainment (positive) velocity w(x, z → ∞) , as obtained for a range of and Fr = 1 .
Here defines the down-slope gradient of surface buoyancy as given by the piece-wise profile (7). When 0 < < 1 the bulk (u) motion remains down-slope but decreases with x leading to a detrainment from the surface layer; if < 0 the down-slope flow increases with x and there is a corresponding entrainment from the ambient fluid. Sufficiently far up/down slope the flow returns to the parallel Prandtl solution (6) which has no entrainment or detrainment. Figure 2 shows that decay back to a Prandtl solution as x increases becomes slower for increasingly large negative values of and we address this behaviour below in Sect. 2.1.
A similarity solution to (4) was obtained by Shapiro & Fedorovich [17]. Their solution assumes an unbounded linear variation of the buoyancy at the surface of the slope, which replaces (7) with b s (x) = (x − x 0 ) at all x locations, for some constant . A solution to (4) can then be sought in the form where x 0 is a constant associated with the choice of a reference position. As noted in [17], the subscript l terms give a standard stagnation-point form of solution; the (u 0 , b 0 , p 0 ) terms are not present in the leading-order boundary-layer ≪ 1 limit. Similar stagnation-point assumptions have been applied to problems of rotating stratified flow as in [13,22]. The solution for (u l , w l , p l , b l ) only depends on the parameter combination Fr because the lengthscale L is absent in the unbounded self-similar formulation. The similarity solution leads to a constant x-independent entrainment/detrainment, which is w l (z → ∞) and represented in Fig. 2 by dashed horizontal lines. We can see that the similarity solution does give an effective quantitative prediction provided that Fr > 0 is not too large. Figure 3 presents the 'peak' flow in/out of the katabatic layer (measured over all x positions) as a function of Fr (solid line) and compares this to the simpler similarity solution (data points). Again the similarity solution remains predictive but diverges rapidly from the two-dimensional computations in the region of Fr ≈ 1.
To explain the observed divergence from the self-similar prediction (as shown in Fig. 3) we consider further corrections to the assumed self-similar form (8) in the shallow limit via with similar expressions for b and p. Linearisation of the tilde quantities gives a one-dimensional eigenvalue problem for . Examination of the spectrum of eigenvalues for increasing values of Fr shows that the real part of changes sign from positive to negative at Fr ≈ 0.9 when Pr = 0.71 . The correction terms therefore grow near x 0 when r < 0 and this delineates the region of applicability of one-dimensional similarity solutions to the two-dimensional computational results.
In addition, the similarity solution has unsteady periodic states [17], which arise through a temporal instability via a supercritical Hopf bifurcation. However, this (self-similar) unsteady instability appears to have no relevance to the full problem defined by (4). Peak entrainment/detrainment Fig. 3 Peak levels of non-dimensional entrainment/detrainment, as measured by the maximum/minimum of w(x, z → ∞) taken over all x positions. The dimensional entrainment is U 0 w(x, z → ∞) as given by (3). Solid line: ≪ 1 solutions of (4) subject to (7) with Fr = 1 and variable . Data points: the similarity solution (8), which only depends on the combination Fr . There is quantitative agreement in general, but this diverges rapidly in the region of Fr near to unity Instead, katabatic flows become temporally unstable to modes of wavelength comparable to the short lengthscale 0 and that grow on a much faster timescale than that defined in (3). We return to describe these instabilities in Sect. 3 below.

Strongly decreasing surface buoyancy
When the surface buoyancy decreases rapidly in the down-slope direction ( | | ≫ 1 with negative) the flow takes larger distances to equilibrate to the Prandtl solution. We can approach this limit directly, by noting that as → −∞ the steady form of (4) (with ≪ 1 ) can be re-scaled via This results in a slightly modified system, to be solved subject to the same homogeneous boundary conditions, but with the re-scaled surface buoyancy B (x,Z = 0) = x for x ∈ [0, 1] and B (x > 1,Z = 0) = −1 . The only difference in this limit is the loss of the forcing term in the energy equation (10c), which is highlighted above by the O(Ū∕| | Fr) term.
Crucially, this leading-order solution differs from the general case in that Ū ≥ 0 and B ≤ 0 for all Z ≥ 0 and there is no oscillation with Z as one approaches the ambient conditions. As we shall describe in the next section, this has consequences for the stability properties of solutions (10) compared to typical solutions of (4), and even (6).
An asymptotic description of (10) is possible (related to the solution of Pohlhausen [23]), which is valid sufficiently far down the slope ( x ≫ 1 ), past the non-uniformly cooled region. This solution thickens in the down-slope direction, increasing in magnitude and with a corresponding (decreasing) entrainment from the ambient fluid. The reduced equations are such that where prime indicates differentiation with respect to the re-scaled coordinate Ẑ .
We cannot neglect the term Ū ∕| |Fr in (10c) at all x positions and when sufficiently far downstream from the non-uniformly cooled region it must be re-included for the solution to eventually recover down slope to a Prandtl state. The lengthscale at which this adjustment occurs can be determined from re-balancing the previously neglected term for large x, for example in (10c) we impose the balance where from (11) . Therefore, eventually the neglected term must be reintroduced, which occurs when x = O(| |Fr) , at which point a down-slope x-invariant Prandtl solution (6) is recovered.

Instability in the shallow-flow limit
We can now consider the stability of the two-dimensional katabatic flows described in the previous section, to temporally growing small-amplitude perturbations. We are free to look for disturbances of any length scale, but if this disturbance length scale is comparable to the slope-normal scale 0 rather than the down slope scale L (with 0 ≪ L ) then the resulting stability problem is determined locally at each x position rather than globally. In addition, if we can find such unstable short waves, they have an associated fast time scale which means they grow rapidly. We therefore introduce These perturbations are applied to a steady two-dimensional solution of (4) which varies on the slower downslope x scale of the non-uniform buoyancy forcing. We then seek a solution in the form where Δ ≪ 1 for small amplitude perturbations and û,v,ŵ,p,b are shape functions that depend on the slope-normal coordinate z. The base flow is two dimensional but we allow for a cross slope flow ( v ≠ 0 ) in the stability problem.
After some rearrangement, the perturbation equations for Δ ≪ 1 become Here we have neglected x derivatives of the base katabatic flow since they are on a length scale much longer than the wavelength of the perturbation. Whilst the base flow equations (4) do not explicitly depend on the slope angle, this is not true for the linearised stability equations (14). The combination ∕Fr can be interpreted as an inverse Reynolds number based on the Prandtl solution velocity scale U 0 (3d) and thickness 0 = L (3a), so in what follows we use Apart from some notation differences, the linear stability problem (14) was presented in the recent work of Xiao and Senocak [14] applied to the classical (linear) one-dimensional solution of Prandtl (6), subject to a zero buoyancy flux condition at the boundary. Their approach determined a critical buoyancy flux (here this would be a critical value for Re ) for instability. In this work we will focus on the shallow-flow limit ( ≪ 1 , Fr = O(1) ) which corresponds to Re ≫ 1.
For simplicity, given the bi-directional nature of the base flow, we will restrict attention here to temporal stability, treating ̂ in (13) as an eigenvalue to be determined. We also restrict attention to two broad classes of disturbance (i) vortices aligned with the down-slope direction ( k x = 0 ) and periodic in the cross-slope direction ( k y ≠ 0 ), and (ii) waves that propagate up/down the slope ( k x ≠ 0 ) and are independent of the cross-slope coordinate ( k y = 0). Downslope vortices In the limit Re ≫ 1 , requiring a balance of the unsteady/advection terms, pressure gradient and the buoyancy forcing term ( b ∕ tan ) of (14d) reveals that and we can reduce (14a) and (14c)-(14e) to where the down-slope x coordinate only appears parametrically via the local buoyancy gradient. This problem is to be solved subject to w = 0 at the boundary ( z = 0 ) and in the far field ( z → ∞).
For large wavenumber k y ≫ 1 the resulting small-wavelength vortex becomes spatially localised and offset from the slope surface. This localisation occurs in a thin layer around z = z 0 which is of width O(k −1∕2 y ) , and z 0 is the point at which the buoyancy gradient term − b∕ z is maximised. In this layer the second-derivative term in (17) is therefore O(k y ) and the leading-order growth rate can be found directly from the other O(k 2 y ) terms. We therefore find that where z 0 is an inflection point in the buoyancy profile with a negative buoyancy gradient.
To determine the leading-order vortex shape and corrections to this growth rate requires continuation to higher order in the large k y expansion. For the special case of a uniform surface buoyancy, the associated Prandtl solution given by (6) has a growth rate bounded above by and the vortices are located at a distance z 0 = √ 2 Pr −1∕4 from the boundary. In Fig. 4 we solve (17) at successive x positions, and present the variation of Ω i √ tan for a range of non-uniform surface buoyancy profiles, as defined by (7) for varying . In this Re ≫ 1 problem the dependence on the slope angle is entirely described by this √ tan factor. If the surface buoyancy increases with x ( > 0 ) the detraining slope flow becomes locally more unstable as seen in Fig. 4a before eventually readjusting to a slower and more stable down-slope Prandtl state. However, for surface buoyancy that decreases with x the entraining flow becomes locally more stable, as seen in Fig. 4b, before returning to a faster and more unstable Prandtl solution for large x.
Gradual stabilisation of the vortex mode for entraining flows is explained by the results of Sect. 2.1. The accelerating solution given by (10) is crucially free (at leading order) from any oscillation with the slope-normal coordinate z. In such a flow there is a monotonic decay of buoyancy flux from the cold surface and it is therefore stable to vortex modes.
In general the vortex is displaced from the boundary at a local minimum of the buoyancy gradient. For the special case of the Prandtl solution (6), this z-location is also where the down-slope flow first changes sign from positive (down the slope) to negative (up the slope). For non-uniform surface buoyancy, the vortex z-location varies with down-slope In Fig. 5 we show an example of the vortex location for increasing x, in a non-uniform surface buoyancy solution of the type shown in Fig. 2. The colour contours show the downslope flow (red for down-slope and blue for up-slope) whilst the vector field is the vortex mode's velocity. The solid/dashed horizontal lines show where the katabatic flow changes sign ( u = 0 ) and where the buoyancy flux is a max/min ( 2 b∕ z 2 = 0 ). At lower values of x the vortex is located around u = 0 (Fig. 5a), but at x = 0.8 (Fig. 5c) the effects of the nonuniform buoyancy forcing are felt and the bulk of the vortex is found in the up-slope flow u < 0 , before returning back to where u = 0 further down the slope (Fig. 5d).
For larger values of k y , viscous effects in the thin vortex layer eventually act to restabilise the flow. The reader is referred to the related discussions in the centrifugal literature, for example Hall [24] or Denier et al. [25]. This restabilisation occurs when k y = O( Re 1∕4 ). Rayleigh waves In the inviscid limit Re ≫ 1 with k y =v = 0 and k x = O(1) we recover a standard Rayleigh problem for the growth rate â gain to be solved subject to ŵ = 0 at the boundary z = 0 and in the far field.
Sufficiently far down the slope we always recover the Prandtl solution (6) which is more/less unstable for increased/decreased down-slope flows associated with < 0 and > 0 respectively. For these Prandtl states the temporal growth rate ̂i peaks at wave numbers near k x ≈ 0.25. Figure 6 shows the variation of the wave's growth rate for a variety of solutions to (4). We typically find a monotonic increase/decrease of the growth rate, however, we must qualify this statement by noting that there is a very slight non-monotonic adjustment in the growth rate near to x = 0 just visible in Fig. 4b.
The vortex modes are strongly stabilised for accelerated/entraining down slope flows (Fig. 4b) because of a reduced oscillation in the buoyancy profile. However the Rayleigh waves ultimately owe their existence to inflectional down-slope velocity profiles, which persist even in the limiting solutions (10)-(11).

Conclusion
Two-dimensional solutions have been obtained for a katabatic flow down a non-uniformly cooled slope, assuming that the down-slope length scale L over which the surface buoyancy varies is long compared to the boundary-layer thickness 0 = √ ∕N sin . For a cooled slope surface, the main flow is down the slope, but there are regions of up-slope return flow. Both the buoyancy and down-slope velocity profiles approach the ambient conditions via a decaying spatial oscillation, which makes the parabolic boundary layer for shallow flows bi-directional.
We apply a numerical scheme to deal with this bi-directionality and compare our two-dimensional solutions with the one-dimensional similarity assumption of Shapiro and Fedorovich [17]. These similarity states are associated with an unbounded surface buoyancy that varies linearly with down-slope coordinate x. Nevertheless, the similarity solutions give good quantitative estimates for the entrainment/detrainment associated with the non-uniform buoyancy region. These estimates become poor beyond limiting value of Fr , where is the down-slope gradient of the surface buoyancy and Fr is an internal Froude number. This limiting value varies with Prandtl number and can be found by an eigenvalue problem associated with the similarity states. The temporal stability of two-dimensional non-uniformly cooled states (to smallamplitude three-dimensional disturbances) is assessed under the assumption that Fr ∕ ≫ 1 . In this limit there are (broadly) two types of (linear) temporal instability, (i) down-slope aligned stationary vortices that are periodic in the cross-slope direction, and of amplitude proportional to (ii) down-slope propagating Rayleigh waves of amplitude proportional to where Ω i , ̂i are the constants presented in Figs. 4 and 6 (for example), T is the dimensional timescale and N the stratification frequency. In this limit we are able to explicitly obtain the dependence of the growth rates on the slope angle ( ), something that is not possible for Fr ∕ = O (1).
Unstable vortex modes are obtained even though the slope surface is being cooled because the buoyancy profile of katabatic flow is oscillatory (Fig. 1). This oscillation means that an unstable region still exists, but it is displaced away from the boundary (Fig. 5).
As can be seen from the above expressions, the wave modes always grow more quickly if = O(1) (steep slopes) but the growth rates are comparable on shallow slopes if The growth rate results of Figs. 4 and 6 show that Ω i can be an order of magnitude larger than ̂i making this balance plausible even when Fr ∕ is large enough for the boundarylayer approximation to be useful. Of course combinations of the two modes in the form of oblique propagating waves are also possible. The presence of both vortex and wave modes at finite values of Fr ∕ (rather than the limiting behaviour described above) was demonstrated by the DNS results of Xiao and Senocak [14]. Of particular interest for katabatic flow is that the behaviour of fully developed turbulent states remains similar to the simple predictions of the Prandtl states (6) despite these instabilities. The role of nonlinearity for these disturbances is clearly a topic for future study.
Finally, while the vortex and wave modes develop with length scales comparable to the katabatic thickness 0 and corresponding time scales 0 ∕U 0 , there are also be non-growing perturbations on the (gravity wave) time scale associated with the stratification component 1∕N sin . These slower scale oscillations decay in entraining flows (and indeed for the Prandtl solution) but can be sustained in detraining flows, leading to a propagating and diffusing oscillation moving into the ambient as described in the appendix below.
This work was motivated by the time REH spent as a post-doc with Prof. PA Davies in the late 90s. Peter, your expertise, hospitality and humour over the intervening years are greatly appreciated! A leading-order balance obtained from substitution of (A3) into (A1) gives and at next order Eliminating the subscript-1 correction terms gives an eigenvalue problem for : where b 0 → 0 as → ±∞ . This has both symmetric and antisymmetric modes, where A is a free constant determined by the history of the flow. Therefore a (linear) disturbance exists with = 0 , for which the Kummer function 1 F 1 simplifies, giving a streamfunction (A1) that is non-decaying: an internal-wave analogue of the inertial-wave behaviour discussed in [12]. This solution corresponds to a perturbation of amplitude Ā , advected with the detraining flow, oscillating with dimensional frequency N sin and diffusing in z. Such solutions require w ∞ > 0 , which is associated with a katabatic flow that decreases in magnitude down the slope owing to non-uniform cooling; that is > 0 in (7).
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http:// creat iveco mmons. org/ licen ses/ by/4. 0/.