Helmholtz Billows Research Articles

Energy transfer from internal solitary waves to turbulence via high-frequency internal waves: seismic observations in the northern South China Sea

Abstract. The shoaling and breaking of internal waves (IWs) are critical processes in the ocean's energy cascade and mixing. Using seismic data, we observed high-frequency internal waves (HIWs), which were primarily distributed in the depth range of 79–184 m. Their amplitude scale is O (10 m), with half-height widths ranging from 154 to 240 m. The shoaling thermocline and gentle slope with a low internal Iribarren number suggest that observed high-frequency internal waves are likely a result of fission. The remote sensing data support this point. Instability estimations showed that, due to the strong vertical shear, the Richardson number (Ri) in the range of 20–30 km was less than 0.25, and Kelvin–Helmholtz (KH) billows can be found in the seismic transect, suggesting that these waves were unstable and might dissipate rapidly. We used the seismic data to estimate diapycnal mixing, and we found that the HIWs can enhance diapycnal mixing, averaging 10−4 m2 s−1. The maximum mixing value is up to 10−3 m2 s−1, and it is associated with the breaking of IWs caused by the strong shear. The results show a new energy cascade route from shoaling internal solitary waves (ISWs) to turbulence, i.e., the fission of ISWs into HIWs, which improves our knowledge of ISW energy dissipation and their roles in improved mixing in the northern South China Sea.

The Effect of Domain Length and Initialization Noise on Direct Numerical Simulation of Shear Stratified Turbulence

Direct numerical simulation (DNS) has been employed with success in a variety of oceanographic applications, particularly for investigating the internal dynamics of Kelvin–Helmholtz (KH) billows. However, it is difficult to relate these results directly with observations of ocean turbulence due to the significant scale differences involved (ocean shear layers are typically on the order of tens to hundreds of meters in thickness, compared to DNS studies, with layers on the order of one to tens of centimeters). As efforts continue to inform our understanding of geophysical-scale turbulence by extrapolating DNS results, it is important to understand the impact of model setup and initial conditions on the resulting turbulent quantities. Given that geophysical-scale measurements, whether through microstructures or other techniques, can only provide estimates of averaged TKE quantities (e.g., TKE dissipation or buoyancy flux), it may be necessary to compare mean turbulent quantities derived from DNS (i.e., across one or more complete billow evolutions) with ocean measurements. In this study, we analyze the effect of domain length and initial velocity noise on resulting turbulent quantities. Domain length is important, as dimensions that are not integer multiples of the natural KH billow wavelength may compress or stretch the billows and impact their energetics. The addition of random noise in the initial velocity field is often used to trigger turbulence and suppress secondary instabilities; however, the impact of noise on the resulting turbulent energetics is largely unknown. In this study, we conclude that domain lengths on the order of 1.5 times the natural wavelength or less can affect the resulting turbulent energetics by a factor of two or more. We also conclude that increasing the amplitude of random initial velocity noise decreases the resulting turbulent energetics, but that different realizations of the random noise field may have an even greater impact than amplitude. These results should be considered when designing a DNS experiment.

Mixing of a cylindrical gravity current in a stratified ambient

Simulation results from direct numerical simulations (DNSs) of cylindrical gravity currents propagating into a linearly stratified ambient with stratification strengths ranging from 0 to 0.8, at a moderate Reynolds number of Re=3450 are presented. The simulations are based on the incompressible Navier–Stokes equations, assuming small density differences such that the Boussinesq approximation is valid. The density of the ambient fluid increases linearly from the top (ρ0) to the bottom (ρb) of the domain. A comparative analysis is conducted between the stratified cases and the unstratified case to investigate the impact of ambient stratification strength on the mixing behaviour of cylindrical gravity currents. The energy conversion processes are analysed using the mechanical framework proposed by Winters et al. (1995). The energy budget is formulated by considering gravitational available potential energy, and the evolution of potential energy due to reversible stirring and irreversible diapycnal mixing in the system. The findings reveal a decrease in available potential energy and kinetic energy with increasing stratification strength, indicating lower energy exchange for gravity currents propagating in the stratified environment. Instantaneous and cumulative mixing efficiency calculations during the slumping phase indicate that Kelvin–Helmholtz billow play an important role in stirring the heavy fluid and causing irreversible mixing with the ambient fluid.

Long-wave instabilities of sloping stratified exchange flows

We investigate the linear instability of two-layer stratified shear flows in a sloping two-dimensional channel, subject to non-zero longitudinal gravitational forces. We reveal three previously unknown instabilities, distinct from the well-known Kelvin–Helmholtz instability and Holmboe wave instability, in that they have longer wavelengths (of the order of 10 to $10^3$ shear-layer depths) and often slower growth rates. Importantly, they can grow in background flows with gradient Richardson number $\gg 1$ , which offers a new mechanism to sustain turbulence and mixing in strongly stratified flows. These instabilities are shown to be generic and relatively insensitive to Reynolds number, Prandtl number, base flow profile and boundary conditions. The nonlinear evolution of these instabilities is investigated through a forced direct numerical simulation, in which the background momentum and density are sustained. The growth of long unstable waves in background flows initially stable to short wave causes a decrease in the local gradient Richardson number. This leads to local nonlinear processes that result in small-scale overturns resembling Kelvin–Helmholtz billows. Our results establish a new energy exchange pathway, where the mean kinetic energy of a strongly stratified flow is extracted by primary unstable long waves and secondary short waves, and subsequently dissipated into internal energy.

Growth and Breakdown of Kelvin–Helmholtz Billows in the Stable Atmospheric Boundary Layer

The development and breakdown of Kelvin–Helmholtz (KH) waves (billows) in the stable atmospheric boundary layer (SABL) and their impact on vertical transport of momentum and scalars have been examined utilizing large eddy simulations. These simulations are initialized with a vertically uniform geostrophic wind and a constant potential temperature lapse rate. An Ekman type of boundary layer develops, and an inflection point forms in the SABL, which triggers the KH instability (KHI). KHI develops with the kinetic energy (KE) in the KH billows growing exponentially with time. The subsequent onset of secondary shear instability along S-shaped braids leads to the turbulent breakdown of the KH billow cores and braids. The frictional ground surface tends to slow down the growth of KE near the surface, reduce the KH billow core depth, and likely suppress other types of secondary instability. KH billows induce substantial down-gradient transport of momentum and sensible heat, which can be further enhanced by the onset of secondary shear instability. Although the KHI-induced strong transport only lasts for around 10–20 min, it reduces vertical shear and stratification in the SABL, enhances surface winds, and results in a 2–3-fold increase in the SABL depth.

Statistical assessment of a Doppler radar model of TKE dissipation rate for low Richardson numbers

Abstract. The potential ability of VHF or UHF Doppler radars to measure turbulence kinetic energy (TKE) dissipation rate ε in the atmosphere is a major asset of these instruments because of the possibility to continuously monitor turbulence in the atmospheric column above the radars. Several models have been proposed over the past decades to relate ε to half the Doppler spectral width σ, corrected for non-turbulent contributions, but their relevance remains unclear. Recently, Luce et al. (2023) tested the performance of a new model expected to be valid for weakly stratified or strongly sheared conditions, i.e., for low Richardson (Ri) numbers. Its simplest expression is εS=CSσ2S, where CS∼0.64 and S=|dV/dz| is the vertical shear of the horizontal wind V. We assessed the relevance of this model with a UHF (1.357 GHz) wind profiler called WPR LQ-7, which is routinely operated at the Shigaraki Middle and Upper Atmosphere (MU) observatory (34.85∘ N, 136.10∘ E) in Japan. For this purpose, we selected turbulence events associated with Kelvin–Helmholtz (KH) billows, because their formation necessarily requires Ri<0.25 somewhere in the flow, a condition a priori favorable to the application of the model. Eleven years of WPR LQ-7 data were used for this objective. The assessment of εS was first based on its consistency with an empirical model εLout=σ3/Lout, where Lout has the dimension of an outer scale of turbulence. It was found to compare well in a KH layer with direct estimates of ε from in situ measurements for Lout≈70 m. Some degree of equivalence between εS and εLout was confirmed by statistical analysis of 192 KH layers found in the height range 0.3–5.0 km indicating that Lout≈LH/0.64, where LH=σ/S is the Hunt scale defined for neutral turbulence. The degree of equivalence is even significantly improved if Lout is not treated as a constant but depends on the depth D of the layer. We found Lout≈0.0875D or equivalently LH∼0.056D, which also means that σ is proportional to the apparent variation in the horizontal velocity (S×D) over the depth of the KH layer. Consequently, εS=0.64σ2S and εLout=σ3/0.0875D would express the same model for KH layers when Ri remains low. For such a condition, we provide a physical interpretation of εLout, which would be qualitatively identical to that for neutral boundary layers.

Kelvin–Helmholtz Billows in the Rising Turbulent Layer During Morning Evolution of the ABL at Dome C, Antarctica

Kelvin–Helmholtz billows (KHBs) within a rising turbulent layer during the transition period from stable to unstable stratification occurring in the morning hours in summertime at the interior of Antarctica (Dome C, Concordia station) are examined in this study. The wave pattern captured by high-resolution sodar echograms from November 2014–February 2015 exhibits regular braid-like structures, associated with Kelvin–Helmholtz shear instabilities. This phenomenon is observed in more than 70% of days in the selected period. Two main regimes of the morning evolution with KHBs are identified roughly, distinguished by the presence or absence of turbulence in the preceding night-time. The weather and turbulent conditions favouring the occurrence of these regimes are analyzed. Also, two distinct patterns of KHBs are identified: (i) quasi-periodical (with periods ≈ 8–15 min) trains containing 5–10 braids, (ii) about continuous series lasting 20–90 min containing 20–80 braids. A composite shape of KHBs is determined. The periodicity of these waves is estimated to be between 20 and 70 s, and their wavelength is estimated roughly to be 100–400 m. The vertical thickness of individual braids at the wave crests ranges between 5 and 25 m. The total depth of a rising turbulent layer containing these waves varies between 15 and 120 m, and the ratio of the wavelength to the depth of the wave layer varies from 3 to 12 with a mean value ≈ 8.2. The morphology of the turbulence structure in the ABL is studied as a function of both temperature and wind field characteristics retrieved from an instrumented 45-m tower and an ultrasonic anemometer-thermometer at 3.5 m. The observational results highlight the necessity of considering the interaction between convective and wave processes when occurring simultaneously.

Shear Instabilities and Stratified Turbulence in an Estuarine Fluid Mud

Abstract This study presents field observations of fluid mud and the flow instabilities that result from the interaction between mud-induced density stratification and current shear. Data collected by shipborne and bottom-mounted instruments in a hyperturbid estuarine tidal channel reveal the details of turbulent sheared layers in the fluid mud that persist throughout the tidal cycle. Shear instabilities form during periods of intense shear and strong mud-induced stratification, particularly with gradient Richardson number smaller than or fluctuating around the critical value of 0.25. Turbulent mixing plays a significant role in the vertical entrainment of fine sediment over the tidal cycle. The vertical extent of the billows identified seen in the acoustic images is the basis for two useful parameterizations. First, the aspect ratio (billow height/wavelength) is indicative of the initial Richardson number that characterizes the shear flow from which the billows grew. Second, we describe a scaling for the turbulent dissipation rate ε that holds for both observed and simulated Kelvin–Helmholtz billows. Estimates for the present observations imply, however, that billows growing on a lutocline obey an altered scaling whose origin remains to be explained.

Diapycnal diffusivities in Kelvin–Helmholtz engendered turbulent mixing: the diffusive–convection regime in the Arctic Ocean

Recent progress in the direct measurement of turbulent dissipation in the Arctic Ocean has highlighted the need for an improved parametrization of the turbulent diapycnal diffusivities of heat and salt that is suitable for application in the turbulent environment characteristic of this polar region. In support of this goal we describe herein a series of direct numerical simulations of the turbulence generated in the process of growth and breaking of Kelvin–Helmholtz billows. These simulations provide the data sets needed to serve as basis for a study of the stratified turbulent mixing processes that are expected to exist in the Arctic Ocean environment. The mixing properties of the turbulence are studied using a previously formulated procedure in which the temperature and salinity fields are sorted separately in order to enable the separation of irreversible Arctic mixing from reversible stirring processes and thus the definition of turbulent diffusivities for both heat and salt that depend solely upon irreversible mixing. These analyses allow us to demonstrate that the irreversible diapycnal diffusivities for heat and salt are both solely dependent on the buoyancy Reynolds number in the Arctic Ocean environment. These are found to be in close agreement with the functional forms inferred for these turbulent diffusivities in the previous work of Bouffard & Boegman (Dyn. Atmos. Oceans, vol. 61, 2013, pp. 14–34). Based on a detailed comparison of our simulation data with this previous empirical work, we propose an algorithm that can be used for inferring the diapycnal diffusivities from turbulent dissipation measurements in the Arctic Ocean.

A comparative study of turbulent stratified shear layers: effect of density gradient distribution

Direct numerical simulations are performed to compare the evolution of turbulent stratified shear layers with different density gradient profiles at a high Reynolds number. The density profiles include uniform stratification, two-layer hyperbolic tangent profile and a composite of these two profiles. All profiles have the same initial bulk Richardson number (Ri_{b,0}); however, the minimum gradient Richardson number and the distribution of density gradient across the shear layer are varied among the cases. The objective of the study is to provide a comparative analysis of the evolution of the shear layers in term of shear layer growth, turbulent kinetic energy as well as the mixing efficiency and its parameterization. The evolution of the shear layers in all cases shows the development of Kelvin–Helmholtz billows, the transition to turbulence by secondary instabilities followed by the decay of turbulence. Comparison among the cases reveals that the amount of turbulent mixing varies with the density gradient distribution inside the shear layer. The minimum gradient Richardson number and the initial bulk Richardson number do not correlate well with the integrated TKE production, dissipation and buoyancy flux. The bulk mixing efficiency for fixed Ri_{b,0} is found to be highest in the case with two-layer density profile and lowest in the case with uniform stratification. However, the parameterizations of the flux coefficient based on buoyancy Reynolds number and the ratio of Ozmidov and Ellison scales show similar scaling in all cases.

Dynamics of a buoyant gravity current propagating in a linearly stratified medium

In this study, we investigate partial- and full-depth buoyant gravity currents propagating along the top surface in a linearly stratified medium. Two- and three-dimensional numerical simulations are performed to study the effect of stratification and initial current depth, on the front speed, internal wave field, and turbulence characteristics. The stratification is varied through a non-dimensional parameter R=ρ0−ρCρb−ρ0, ranging between 0.04 and 85, where ρC is the constant bulk density of the current fluid and ρ0,ρb represent the densities of the ambient fluid at the top and bottom surfaces, respectively. For large values of R (ρ0−ρC≫ρb−ρ0), we observe that the resulting Froude number (Fr=UNH) is greater than 1/π, and the flow is characterized as supercritical, where the front speed exceeds the long wave speed. In the supercritical regime, Kelvin–Helmholtz billows are prominently seen along with an internal solitary wave, which propagates with the density front. As the R value decreases, the relative strength of the ambient stratification increases when compared to the horizontal density difference at the top surface, leading to a subcritical flow regime where the front speed is smaller when compared to the long-wave speed. The Kelvin–Helmholtz billows and the solitary wave gradually disappear, and vertically propagating high-mode internal waves are prominently seen for R< 1. Quantification of the Froude number for various values of R and h/H shows that it follows a power law, Fr∝(h/H×R)1/2, with the proportionality constant 0.72. This scaling works well for all the partial-depth cases considered in this study, i.e., h/H= 1/8, 1/6, 1/4, and 1/3, while a slight deviation is observed for the full-depth gravity currents that correspond to h/H= 1.

Stationary internal hydraulic jumps

This is a theoretical and laboratory study of stationary internal hydraulic jumps. These jumps are rapid transitions between internally supercritical flow, generated by placing a sill on the bed of a horizontal rectangular channel, and internally subcritical flow, generated by installing a downstream contraction. This contraction generates an approximately uniform flow downstream of the jump; thus mimicking barotropically driven two-layer flows, as found in tidally driven flows over underwater sills, and flows over mountain ranges driven by large-scale pressure gradients. Upstream of the jump a train of Kelvin–Helmholtz billows forms on the interface between the layers. Upper layer fluid is entrained into these billows, which are subsequently advected into the lower portion of the jump. These billows are broken down by the turbulence of the jump, and the entrained upper layer fluid is mixed with lower layer fluid. Downstream of the jump the upper layer remains homogeneous, the density step at the interface is weakened, the upper portion of the lower layer is approximately linearly stratified, and the lower portion of the lower layer is undisturbed. This altered density profile is the downstream conjugate state of the jump. When the contraction is narrowed the jump moves upstream and ‘drowns’ part of the train of billows, reducing the amount of entrainment. Thus, while the jump is responsible for mixing fluid from the upper layer into the lower layer, it is the position of the jump relative to the upstream train of billows that determines the amount of entrainment.

Impact of Elevated Kelvin–Helmholtz Billows on the Atmospheric Boundary Layer

AbstractThe impact of Kelvin–Helmholtz billows (KHBs) in an elevated shear layer (ESL) on the underlying atmospheric boundary layer (BL) is examined utilizing a group of large-eddy simulations. In these simulations, KHBs develop in the ESL and experience exponential growth, saturation, and exponential decay stages. In response, strong wavy motion occurs in the BL, inducing rotor circulations near the surface when the BL is stable. During the saturation stage, secondary instability develops in the ESL and the wavy BL almost simultaneously, followed by the breakdown of the quasi-two-dimensional KH billows and BL waves into three-dimensional turbulence. Consequently, during and after a KH event, the underlying BL becomes more turbulent with its depth increased and stratification weakened substantially, suggestive of significant lasting impact of elevated KH billows on the atmospheric BL. The eventual impact of KHBs on the BL is found to be sensitive to both the ESL and BL characteristics.

Convection experiments with the exponential time integration scheme

We present the construction and preliminary tests of the exponential time integration schemes for the Euler equations. The equations are written in terms of potential temperature and Exner function. The time stepping is accomplished with the EPI2 and EPI3 schemes which have been used in the past for very accurate, long time step integration of the shallow water equations on the sphere. The stability of exponential integrators of this class is ensured by approximation of the oscillatory term using an absolutely stable exponential formula. Consequently, the schemes are not subject to the CFL condition imposed by the linear oscillatory part of the system and they also eliminate the phase errors associated with all the known implicit time stepping algorithms. All calculations are performed using the phipm algorithm based on the Krylov space methods. It is shown that the spectral radius of the Jacobian and consequently the Krylov space dimension are very sensitive to the presence of spatial noise. Efficiency and accuracy increase after the Shapiro filter is applied every two time steps. The method is investigated using several standardized and well accepted benchmark tests, including convective bubbles and a cold density current proposed by Straka and collaborators. Numerical experiments show that a stable integration is obtained for the wave Courant number of the order of 60 to 250. We present a detailed discussion of the development of a cold density current capturing the spiraling Prandtl vortex and the subsequent appearance of the Kelvin–Helmholtz billows. Our results from the exponential method compare very well with those obtained from the accurate schemes reported in the recent literature. The algorithm accurately handles both advection as well as acoustic and gravity waves. These results warrant the further development of the scheme for modeling atmospheric processes on a wide range of scales.

Vertical Nitrate Flux Induced by Kelvin–Helmholtz Billows Over a Seamount in the Kuroshio

Kelvin–Helmholtz (KH) billows can facilitate microscale turbulent mixing around seamounts in the Kuroshio. This study sought to describe the influence of billow intensity (i.e., “intermittent and small” and “steady and large” billows) on vertical nitrate fluxes. KH billows led to turbulent kinetic energy dissipation rates [ε = O (10–7 to 10–6) W kg–1] and eddy diffusivities [Kρ = O (10–4 to 10–3) m2 s–1] that were significantly stronger than those outside the billow depths. The mean nitrate flux estimated using Kρ in the billow depths had a maximal value of 10.0 mmol m–2 day–1, which was much higher than estimates for the open ocean. The nitrate flux associated with the shallow KH billows at two vertical levels contributed to enrich the subsurface phytoplankton maximum, while the deeper billows closer to the summit were found to induce a large amount of nitrate flux from the deeper to the subsurface water. This study showed that KH billows make important contributions to seamount ecosystems, particularly in the cycling and vertical mixing of nutrients to make them available for potential downstream transport.

Large-Aspect-Ratio Structures in Simulated Ocean Surface Boundary Layer Turbulence under a Hurricane

AbstractA large-eddy simulation (LES) initialized and forced using observations is used to conduct a process study of ocean surface boundary layer (OSBL) turbulence in a 2-km box of ocean nominally under Hurricane Irene (2011) in 35 m of water on the New Jersey shelf. The LES captures the observed deepening, cooling, and persistent stratification of the OSBL as the storm approaches and passes. As the storm approaches, surface-intensified Ekman-layer rolls, with horizontal wavelengths of about 200 m and horizontal-to-vertical aspect and velocity magnitude ratios of about 20, dominate the kinetic energy and increase the turbulent Prandtl number from about 1 to 1.5 due partially to their restratifying vertical buoyancy flux. However, as the storm passes, these rolls are washed away in a few hours due to the rapid rotation of the wind. In the bulk OSBL, the gradient Richardson number of the mean profiles remains just above (just below) 1/4 as the storm approaches (passes). At the base of the OSBL, large-aspect-ratio Kelvin–Helmholtz billows, with Prandtl number below 1, intermittently dominate the kinetic energy. Overall, large-aspect-ratio covariance modifies the net vertical fluxes of buoyancy and momentum by about 10%, but these fluxes and the analogous diffusivity and viscosity still approximately collapse to time-independent dimensionless profiles, despite rapid changes in the forcing and the large structures. That is, the evolutions of the mean temperature and momentum profiles, which are driven by the net vertical flux convergences, mainly reflect the evolution of the wind and the initial ocean temperature profile.

Near- and far-field structure of shallow mixing layers between parallel streams

Abstract

Turbulence, Low-Level Jets, and Waves in the Tyrrhenian Coastal Zone as Shown by Sodar

The characteristics of the vertical and temporal structure of the coastal atmospheric boundary layer are variable for different sites and are often not well known. Continuous monitoring of the atmospheric boundary layer was carried out close to the Tyrrhenian Sea, near Tarquinia (Italy), in 2015–2017. A ground-based remote sensing instrument (triaxial Doppler sodar) and in situ sensors (meteorological station, ultrasonic anemometer/thermometer, and net radiometer) were used to measure vertical wind velocity profiles, the thermal structure of the atmosphere, the height of the turbulent layer, turbulent heat and momentum fluxes in the surface layer, atmospheric radiation, and precipitation. Diurnal alternation of the atmospheric stability types governed by the solar cycle coupled with local sea/land breeze circulation processes is found to be variable and is classified into several main regimes. Low-level jets (LLJ) at heights of 100–300 m above the surface with maximum wind speed in the range of 5–18 m s−1 occur in land breezes, both during the night and early in the morning. Empirical relationships between the LLJ core wind speed characteristics and those near the surface are obtained. Two separated turbulent sub-layers, both below and above the LLJ core, are often observed, with the upper layer extending up to 400–600 m. Kelvin–Helmholtz billows associated with internal gravity–shear waves occurring in these layers present opposite slopes, in correspondence with the sign of vertical wind speed gradients. Our observational results provide a basis for the further development of theoretical and modelling approaches, taking into account the wave processes occurring in the atmospheric boundary layer at the land–sea interface.

Sodar Observation of the ABL Structure and Waves over the Black Sea Offshore Site

Sodar investigations of the breeze circulation and vertical structure of the atmospheric boundary layer (ABL) were carried out in the coastal zone of the Black Sea for ten days in June 2015. The measurements were preformed at a stationary oceanographic platform located 450 m from the southern coast of the Crimean Peninsula. Complex measurements of the ABL vertical structure were performed using the three-axis Doppler minisodar Latan-3m. Auxiliary measurements were provided by a temperature profiler and two automatic weather stations. During the campaign, the weather was mostly fair with a pronounced daily cycle. Characteristic features of breeze circulation in the studied area, primarily determined by the adjacent mountains, were revealed. Wave structures with amplitudes of up to 100 m were regularly observed by sodar over the sea surface. Various forms of Kelvin–Helmholtz billows, observed at the interface between the sea breeze and the return flow aloft, are described.

Kelvin–Helmholtz billows above Richardson number

We study the dynamical system of a two-dimensional, forced, stratified mixing layer at finite Reynolds number $Re$, and Prandtl number $Pr=1$. We consider a hyperbolic tangent background velocity profile in the two cases of hyperbolic tangent and uniform background buoyancy stratifications, in a domain of fixed, finite width and height. The system is forced in such a way that these background profiles are a steady solution of the governing equations. As is well known, if the minimum gradient Richardson number of the flow, $Ri_{m}$, is less than a certain critical value $Ri_{c}$, the flow is linearly unstable to Kelvin–Helmholtz instability in both cases. Using Newton–Krylov iteration, we find steady, two-dimensional, finite-amplitude elliptical vortex structures – i.e. ‘Kelvin–Helmholtz billows’ – existing above $Ri_{c}$. Bifurcation diagrams are produced using branch continuation, and we explore how these diagrams change with varying $Re$. In particular, when $Re$ is sufficiently high we find that finite-amplitude Kelvin–Helmholtz billows exist when $Ri_{m}>1/4$ for the background flow, which is linearly stable by the Miles–Howard theorem. For the uniform background stratification, we give a simple explanation of the dynamical system, showing the dynamics can be understood on a two-dimensional manifold embedded in state space, and demonstrate the cases in which the system is bistable. In the case of a hyperbolic tangent stratification, we also describe a new, slow-growing, linear instability of the background profiles at finite $Re$, which complicates the dynamics.