Generation of internal solitary waves by frontally forced intrusions in geophysical flows

Internal solitary waves are hump-shaped, large-amplitude waves that are physically analogous to surface waves except that they propagate within the fluid, along density steps that typically characterize the layered vertical structure of lakes, oceans and the atmosphere. As do surface waves, internal solitary waves may overturn and break, and the process is thought to provide a globally significant source of turbulent mixing and energy dissipation. Although commonly observed in geophysical fluids, the origins of internal solitary waves remain unclear. Here we report a rarely observed natural case of the birth of internal solitary waves from a frontally forced interfacial gravity current intruding into a two-layer and vertically sheared background environment. The results of the analysis carried out suggest that fronts may represent additional and unexpected sources of internal solitary waves in regions of lakes, oceans and atmospheres that are dynamically similar to the situation examined here in the Saguenay Fjord, Canada.

I nternal solitary waves 1-3 are hump-shaped nonlinear and nonhydrostatic gravity waves that propagate horizontally along density waveguides in stratified geophysical fluids. Although patchy and episodic, these waves are known to be widespread and energetic enough to participate in global-scale air and water mass modifications, as they may break into stratified turbulence when encountering destabilizing conditions 4 . Although these waves, once formed, are routinely observed in lakes 5 , oceans 3,6 and the atmosphere 7 , their origin often remains elusive. The current understanding is that these waves can be generated by stratified flows over topography [8][9][10][11][12][13][14] , by subharmonic interactions 15 , by internal tide, seiche and Kelvin wave steepening 5,16,17 , by tidal beams impinging on the pycnocline 18,19 or, and importantly for this paper, by the gravitational collapse of mixed fluids, equivalently referred to as intrusive gravity currents [20][21][22][23][24] .
While the gravitational collapse of mixed fluids mechanism has received considerable research attention from the theoretical and experimental points of view since the 1980s, there has been little clear field evidence to support whether this mechanism occurs in geophysical fluids in general, and in oceans in particular. One exception for the ocean is the reporting of internal solitary waves generated by the Columbia River plume flowing into the stratified coastal Pacific 25 ; an upside-down analogous situation to that studied in lock-exchange laboratory experiments where dense water is suddenly left to flow underneath a two-layer still ambient 20,22 . In the atmosphere, the Morning Glory phenomenon observed over North-Western Australia is also thought to arise from the generation of internal solitary waves by an intruding bottom gravity current in the stratified ambient 26,27 . Bottom density currents produced by thunderstorm outflows have also been reported as a potential source of atmospheric internal solitary waves 28 .
The generation of internal solitary waves by the gravitational collapse of mixed fluids is not limited to cases of surface (for example, river plumes) or bottom (for example, Morning Glory) density currents but may also arise, under some circumstances, from interfacial gravity currents or intrusions, in two-layer fluids, as reproduced in several lock-exchange laboratory experiments [20][21][22][23][24]29,30 . One important result of laboratory experiments is that large-amplitude internal solitary waves that are able to detach and lead the interfacial gravity current head are only generated when the initial motionless condition is said to be in non-equilibrium, that is, when the centre of mass of the lock water rises or falls when the gate is removed [21][22][23]31 . For full-depth lock-exchange experiments, that is, experiments where the lock water initially occupies the full tank depth 31 , where the two layers of the ambient fluid are of equal thicknesses, equilibrium situations exist when the density of the intrusion equals the mean density of the ambient, a condition easily achieved by a complete mixing of the two-layer ambient water on the lock side of the gate. To obtain non-equilibrium initial conditions, either salt is added to the mixture or a certain volume of the mixed lock water is replaced by an equivalent volume of freshwater 22 .
The generation of internal solitary waves by non-equilibrium interfacial gravity currents is not yet fully understood partly because the wide parameter space that defines this problem has only been partially explored 31 . For example, the case where the lock water initially occupies only a fraction of the total water depth, a situation referred to as partial-depth intrusion, remains little studied, although likely of greater environmental significance. In any case, whether or not non-equilibrium interfacial gravity currents are an important source of internal solitary waves in geophysical fluids remains to be demonstrated, especially since it has been suggested that non-equilibrium geophysical situations are rare and that equilibrium situations are 'notably of greater geophysical relevance' 31 .
Here we report on geophysical observations of the generation of internal solitary waves by what can be interpreted as being non-equilibrium stratification, partial-depth 31 and frontally forced interfacial gravity currents intruding into a quasi two-layer and vertically sheared ambient. This is a complex natural situation that, to our knowledge, has not been fundamentally studied theoretically or experimentally nor clearly observed in geophysical fluids before. For example, laboratory experiments on intrusive gravity currents have only considered still ambient, and non-equilibrium situations have been studied for full-depth lock-exchange problems, but hardly for partialdepth 31 . We conclude that fronts may be important sources of internal solitary waves in geophysical fluids that are dynamically similar to the situation examined here in the Saguenay Fjord, Canada (Fig. 1).

Results
Observations. The essence of our echo-sounder, current and density measurements is summarized by two repeated transect lines carried out in a region of the Saguenay Fjord 32 (Fig. 2) known to host internal solitary waves of unknown origin 33 . These data were acquired B1 h before the time of low water (20:30 UTC on 10 July 2015 at Tadoussac) during a neap ebb tide of range DZ ¼ 2.9 m. The first transect (Fig. 2a-c) shows a prominent solitary wave of amplitude a ¼ 7.4 m and maximum vertical velocities w s ¼ ± 0.2 m s À 1 propagating northward along a pycnocline of buoyancy frequency N pyc ¼ 0.25 s À 1 located at depth h ¼ 7.6 m in the undisturbed state ahead. The total depth H along this transect varies between 100 and 150 m (not shown). Idealized numerical simulations (presented below) indicate that such a solitary wave in this stratified and sheared background environment has horizontal lengthscale l ¼ 50 m and induces a perturbation mechanical energy E 0 ¼ 1.9 MJ m À 1 (see Methods for details). This wave is therefore highly nonlinear 3 (a/h ¼ 1), long relative to the surface layer thickness (l/h ¼ 7) but fairly short relative to the bottom layer thickness (H À h) that varies over the transect in the range 92-142 m such that 0.4ol/(H À h)o0.5. The fact that the wave is short relative to the bottom layer thickness will become important below for interpreting these observations. This wave was approximately the fourth wave seen in this area, and the previous ones seen had greater amplitude, up to a ¼ 10 m (Fig. 3).
The wave of Fig. 2a is followed by what can be interpreted from the echogram as being an intermediate turbulent layer, centred B10 m depth and being roughly 20 m thick, populated with Kelvin-Helmholtz billows near the top edge. Eleven density overturns were detected (see Methods for details) throughout this transect (Fig. 2a), of which seven were within the billowing layer behind the solitary wave. These seven overturns are characterized by mean Thorpe scales L T ¼ 0.92 m (s.d. ¼ 0.36 m). Assuming a 1:1 ratio between Thorpe and Ozmidov scales 34 gives an order of magnitude for the mean dissipation rate O E ð Þ $ 10 À 5 W kg À 1 within the overturns, consistent with intense ocean turbulent conditions 17 and supporting the echogram interpretation that this layer is turbulently mixing.
These observations may suggest that the turbulent layer was produced by the passage of the wave in a manner comparable to that seen in idealized numerical simulations of unstable solitary waves where turbulent fluid is left behind as the breaking wave keeps moving forward and restabilizes. This can happen for example when the underlying bottom topography shallows sufficiently to cause wave breaking 35,36 . However, this hypothesis is eliminated here because the wave remains short along this transect relative to the bottom layer thickness (0.4ol/(H À h)o0.5) such that it is virtually unaffected by the underlying changing bottom topography. This is supported by a numerical simulation we have carried out (not shown) in which the propagation of a comparable wave over the changing bottom topography was simulated. As expected, the result showed negligible change in wave shape and behaviour relative to a flat-bottom simulation.
Another possibility to explain the overturning and billowing fluid seen behind the wave could be that it was produced by shear instability that had previously developed within the solitary wave 37 . Qualitatively, parallels could well be made between the observations of Fig. 2a and the numerical results presented in Lamb and Farmer 37 (their Figs 13,14 and 18). However, there are reasons to also reject this hypothesis. First, the vorticity seen in the observed billowing layer is opposite to that expected from shear instability that would arise from such a northward propagating internal solitary wave (Fig. 4). For example, the billow centred around y 0 ¼ À 275 m and at 6 m depth rotates anti-clockwise, while billows arising from wave-induced shear instability would be expected to rotate clockwise 37 . Another reason to reject this hypothesis is that the reduced stratification within the wave is N 2 À S 2 /440 (equivalent to N 2 /S 2 41/4) 38 such that there is no indication that this wave may have been dynamically unstable.
An alternative and more plausible interpretation is that the Kelvin-Helmholtz billows arose from the vertical shear associated with an intrusive layer that lies just above the pycnocline and characterized by northward velocity of up to v ¼ 0.5 m s À 1 (Fig. 2b,e). This interpretation is supported by the fact that the reduced stratification averaged over the length of the intrusion ðN 2 À S 2 =4Þ is negative near the top edge of the sheared interface (Fig. 5). The Reynolds number Re ¼ 5 Â 10 6 is around three orders of magnitude higher than typical laboratory settings of intrusive gravity currents in two-layer fluids 22 . This near-surface intrusion differs from deeper intrusions found in other fjords [39][40][41][42] and in the Saguenay 32 by the fact that it intrudes through the The observations suggest that this dynamically unstable intrusion was forced by a small-scale convergent front found at the southernmost end of the transect (y 0 ¼ À 680 m; Figs 2 and 4). This front is characterized by steep isopycnals angled at 25°, strong backscatter intensity indicative of bubble entrainment typical of such fronts 43 , a large downward vertical velocity w front ¼ À 0.2 m s À 1 within the first 10 m or so of the water column and of width LE20 m, a surface density jump Dr ¼ 4 kg m À 3 and Rossby number where f ¼ 1.08 Â 10 À 4 s À 1 is the Coriolis parameter and vE0.5 m s À 1 . This front falls into the type of submesoscale fronts that are ubiquitous in the coastal ocean 44,45 .
By the second transect ( Fig. 2d-f), the leading solitary wave previously seen at 19:20:52 UTC had detached from the head of the intrusion and was captured again 635 m northward at 19:45:12 with an amplitude a ¼ 4.4 m. This wave had travelled against the surface current at an average ground speed c ¼ 0.4 m s À 1 . The reason for its reduced amplitude may be due to radial spreading loss as there is no evidence that this wave had lost amplitude due to wave breaking. Evidence of internal solitary wave radial spreading can be seen in this environment from shore-based cameras (for example, Fig. 1). A few other internal solitary waves of smaller amplitudes were also seen trailing behind, as if these waves had somehow been ejected from the head of the intrusion (Fig. 2d), an hypothesis that will be supported below with numerical simulations. During this second transect, the intrusion was still clearly defined and more Kelvin-Helmholtz billows were visible near the top of its sheared interface (Fig. 2d, inset).
Shore-based photography suggests that the phenomenon was phase-locked with the semi-diurnal tide as surface patterns indicative of the concomitant and proximate existence of fronts  ) and (c,f) upward velocity w (ms À 1 ). The black contours on b and c are isopycnals 2 kg m À 3 apart based on sorted density profiles. The 'Z' symbols on a mark the centre of detected density overturns and their size is proportional to the logarithm of the available potential energy of the fluctuation P found in the range 4 Â 10 À 6 À 3 Â 10 À 4 W kg À 1 . The ' þ ' symbols between a and b mark the positions of the 20 temperature-salinity profiles collected. The inset highlights Kelvin-Helmholtz billows within the dashed box. These transects correspond to the grey track presented in Fig. 1. and internal solitary waves were seen on the 6, 7, 9 and 10 of July in the same area and around the same tidal phase, that is, an hour or so before the time of low water at Tadoussac. July 8 was too windy to capture any interpretable surface patterns. The clearest camera and underwater evidence of this process in action was captured on the 6 ( Fig. 1). Note that at the time these data were acquired, we were not aware that this process was taking place. It is by chance that we happened to carry out an echo-sounder transect at that place and time, and it is only during the subsequent data analysis phase of this research that we realized we had sampled something unusual. This also explains why the echo-sounder transect does not extend across the front on that day. At the time of sampling, we were not preoccupied by fronts. Fronts and internal solitary waves will however be the focus of future field experiments.
Numerical simulations. Next, we carried out numerical simulations to determine whether this unusual series of concomitant phenomena (Kelvin-Helmholtz billows, internal solitary waves, intrusion and front) seen on Fig. 2 was simply coincidental or whether they were tied together by a single process, possibly related to intrusive gravity currents. Given the number of different nonlinear and nonhydrostatic processes involved, we based our analysis on the two-dimensional incompressible Euler equations. These are the simplest equations that can explicitly deal with all salient features of this problem (see Methods for details).
A first numerical simulation was carried out in a 150 m deep flat-bottom open channel initialized with a three-layer idealization of the observed front (Figs 5 and 6a,d,g; Methods for details). Importantly for the rest of this demonstration, the ambient is not initially still but is characterized by a convergent flow at the front, and the background environment is sheared, a complex natural situation not examined before in laboratory experiments. This configuration is equivalent to a partial-depth intrusion because   the intrusive fluid of density r 1 does not initially occupy the entire water depth 31 . For partial-depth and initially motionless intrusions, the degree of non-equilibrium can be characterized by the following non-dimensional density ratio 31 where g 0 ab ¼ g    The blue curves represent the control run based on an idealization of the field observations (see Fig. 6 for a more explicit representation of this run). Other curves represent sensitivity tests carried out using the parameters listed in Table 1. Blue, control run #1; black, run #2; magenta, run #3; green, run #4.
The result of this control run shows remarkable similarities with the field observations (Fig. 6). Such a convergent front forces a dynamically unstable intrusion comparable in thickness and velocity to the observations and characterized by roughly 5 m tall Kelvin-Helmholtz billows developing on the upper interface. Internal solitary waves are generated in both directions. A few southward propagating ones with amplitudes ao4 m are generated near the front during the first few minutes. The leading wave induces an energy perturbation E 0 ¼ 2.0 MJ m À 1 . These southward propagating waves may be similar to those generated by river plumes 25 with the difference that, here, wave generation does not require the flow to be initially supercritical. We cannot confirm whether such smaller amplitude southward propagating waves were also emitted in the Saguenay Fjord as sampling did not extend far enough southward of the front.
The most remarkable solitary waves are generated at the head of the intrusion and propagate northward. The waves are rank-ordered with the leading wave having, at t ¼ 1,100 s and y ¼ 901 m, amplitude a 1 ¼ 14.1 m, horizontal lengthscale l 1 ¼ 70 m, phase speed c 1 ¼ 0.91 m s À 1 and energy perturbation E 0 ¼ 7.8 MJ m À 1 . The head of the intrusion continuously emits waves northward with one wave being emitted roughly every six buoyancy periods (t ¼ 2p/N 2 ¼ 25 s). This provides a continuous northward depth-integrated perturbation energy flux F 0 that reaches an instantaneous value of max(F 0 ) ¼ 79 kW m À 1 for the leading wave and then decreases with time (Fig. 7). The total wavetrain perturbation energy induced during about 30 min ðE 0 ¼ R 50 min 18 min F 0 dtÞ, which is a lower bound estimate of the events duration observed in the field, is E 0 ¼ 37.4 MJ m À 1 .
Three other simulations were carried out to explore the sensitivity of the resulting internal solitary waves to whether the initial condition is in equilibrium or not and/or the stillness of the ambient ( Table 1). The first of these tests examines the situation where the initial density field is identical as the control run discussed above (that is, x ¼ À 0.18) but where the ambient is still (run #2 in Table 1). This corresponds to a non-equilibrium situation where, according to laboratory experiments, largeamplitude solitary waves leading the intrusion head are expected to be generated. This is indeed what the numerical results show but the northward propagating waves generated are at least an order of magnitude less energetic than the control run (Fig. 7, black). The second of these tests examines the still ambient with the initial density field in equilibrium (that is, x ¼ 0, run #3 in Table 1). Virtually, no waves are generated, in agreement with laboratory experiments (Fig. 7, magenta). Finally, we examine the case where the initial density field is in equilibrium but the ambient fluid is convergent and sheared (that is, x ¼ 0, run #4 in Table 1). Interestingly, this equilibrium situation generates very large-amplitude internal solitary waves (Fig. 7, green), much larger than the non-equilbrium and still simulation (Fig. 7, black), and almost as large and as energetic as the non-equilibrium control run (Fig. 7, blue). These results show that intruding gravity currents, when forced, can generate large-amplitude internal solitary waves even when the initial state at rest is in equilibrium.

Discussion
The model results support the hypothesis that the observed concomitant phenomena reported here (front, intrusion, Kelvin-Helmholtz billows, internal solitary waves; Fig. 2) are not simply randomly coincidental but arise from a single physical process that could be called a forced gravitational collapse of mixed fluids or a frontally forced intrusion. We conclude from this research that partial-depth, convergent fronts in stratified, sheared but initially dynamically stable (Ri41/4) and subcritical (Fro1) geophysical turbulent fluids (Re410 6 ) can force dynamically unstable intrusions that trigger Kelvin-Helmholtz billows, leading to strong turbulent dissipation and mixing, while continuously radiating very large-amplitude (a/hB1) internal solitary waves through the head of the intrusion. While non-equilibrium stratifications (xa0) favour the generation of larger amplitude waves, it is not a condition for generating largeamplitude internal solitary waves in frontally forced intrusions.
As this is a rather new problem, it still remains to be determined how the intrinsic properties of the intrusion, Kelvin-Helmholtz billows, internal solitary waves and turbulence vary across the complex and rich parameter space that defines this problem (for example, Re, Ri, Fr, Ro, x and D). Such relationships are required before turbulence and internal wave fluxes arising from dynamically similar geophysical intrusions and fronts as reported here could be parameterized and incorporated into large-scale models. Furthermore, historical observations should be re-examined in the light of this finding as in some circumstances fronts and forced intrusions could explain elusive field observations of internal solitary waves.

Methods
Observations. A field experiment took place during 5-11 July 2015 near Ansede-Roche in the Saguenay Fjord (Fig. 1). Measurements were collected from an 8 m research boat equipped with a towed 120 kHz Biosonics echo-sounder, a towed 500 kHz Teledyne RD Instruments 5-beam acoustic Doppler current profiler (ADCP) and a Sea-Bird Electronics SBE-19plus conductivity-temperature-depth profiler. The ADCP data were acquired in 5 s ensemble averages and 0.5 m vertical bin size. The data were then averaged to 1.0 m vertical resolution and 15 s temporal resolution. The ADCP had no bottom-tracking feature. The flow measurements were referenced to ground using global positioning system (GPS) data acquired every second. A 20-megapixel Canon EOS 6D camera with a 60 degree field of view was programmed to take one image every minute during daytime from the balcony of Mrs Simard at 48°13.439 0 N, 69°52.513 0 W and altitude of 89 m, pointing towards the southwest. Images were then georectified using landscapes and the boat as ground control points 33 .
The in situ data presented here were collected on 6 July between 16:11 and 16:20 UTC (Fig. 1) and on 10 July 2015 between 19:15 and 19:48 UTC (Fig. 2). The first transect presented (Fig. 2a-c) was carried out while our boat was freely drifting southward with surface currents, while 20 temperature-salinity casts were collected as quickly as possible down to a depth of B40 m to obtain the highest temporal resolution across the pycnocline. This sampling was nevertheless insufficient to fully resolve high-frequency internal solitary waves given that at a sampling fall speed of B1 ms À 1 it took roughly the same time to carry out one profile as the buoyancy period (E30 s). Once this sampling finished, we steamed back along more or less the same transect line while towing our echo-sounder and acoustic current profiler without performing en-route temperature-salinity casts (Fig. 2d-f), explaining the absence of density measurements on the second transect. En-route current measurements were subject to greater uncertainties due to movements of the towed body and bubble entrainment while steaming. This explains the sparsity of the current data on the second transect. This sampling sequence was guided by visual observations of sea surface patterns. Note that this sequence did not visually appear to be any different than the previous tens of internal solitary waves we had sampled during this field experiment. We had no clear real-time indications that we were sampling something unusual such that we called it a day at the end of the second transect.
The buoyancy frequency N 2 ¼ À (g/r 0 )qr/qz, where r is water density, g ¼ 9.81 m s À 2 is Earth's gravitational acceleration, r 0 ¼ 1,023 kg m À 3 is a reference density and z the upward axis, was calculated from gravitationally sorted density profiles and decimated at 1 m vertical resolution. The vertical shear was calculated as S 2 ¼ (qu/qz) 2 þ (qv/qz) 2 , where u and v are the eastward and northward velocity components, respectively.
Overturns were detected within density profiles 34 , where a minimum run length of 4 was found to distinguish unstable features from noise, which were then screened for temperature and salinity co-variance. Only overturns larger than order 1 m could therefore be detected. A total of 11 were found, for which their associated dissipation rates were estimated as the available potential energy of the fluctuations multiplied by N calculated on the sorted density profiles. Given that turbulence quantities are log-normally distributed, of the mean dissipation rate of turbulence kinetic energy was estimated as E ¼ exp m þ s 2 =2 ð Þ , where m and s 2 are the arithmetic mean and variance of ln(E/(W kg À 1 )), respectively 46 .
Model. The numerical model used here solves with second-order finite differences the following two-dimensional, non-rotating Euler equations in the vertical y-z where t is time, y the meridional axis (positive northward), z the vertical axis (positive upward), v and w the horizontal and vertical velocity, respectively, Z the position of the free surface, r the density, r 0 ¼ 1,023 kg m À 3 a constant reference density, and p the pressure. The Poisson equation needed for solving the nonhydrostatic pressure is solved using the Pardiso solver [48][49][50] . Implicit numerical viscosity and diffusion associated with the second-order flux-limiter advection scheme 47 insures numerical stability. Simulations were carried out for a flat-bottom rectangular domain of depth H ¼ 150 m, which roughly corresponds to the depth at the sampling transect. The initial density r and horizontal velocity v fields that define the convergent front were set as the following three-layer system with where with (r 0 , r 1 , where with (v 01 , v 02 , v 1 , v 2 ) ¼ ( À 0.4000, À 0.3088, 0.3000, À 0.8000) m s À 1 . This initial condition is such that The initial vertical velocity w is numerically computed by integrating over z the continuity equation 5 given v(y, z, This is an idealization of the observed conditions found on each side of the front (Figs 5 and 6a,d,g). The square of the buoyancy frequencies and shears at the pycnocline on each side of the front are, respectively, N 2 1 ¼ À (g/r 0 ) (r 1 À r 0 )/(2d 1 ) ¼ 0.023 s À 2 , N 2 2 ¼ À (g/r 0 ) (r 2 À r 0 )/(2d 2 ) ¼ 0.061 s À 2 , Þis the reduced gravity on each side of the front, respectively. Given that these Froude numbers are o1 indicates that long linear internal waves could escape and propagate in both directions on each side of the front.
These parameters listed above are for the control run. Three other simulations were carried out to test the sensitivity of the results to some change in the parameters. The parameters used for these additional runs are listed in Table 1.
The horizontal numerical domain extends in the interval y ¼ ±200 km. The grid size varies with highest resolution Dy ¼ 0.5 m in the [ À 1 km, 1 km] interval and the resolution steadily decreases to Dy ¼ 1 km outside the central domain of interest. This long coarse resolution domain acts as a sponge layer that numerically absorbs barotropic waves that escape during the initial adjustment and ensures that the solution is not affected by open boundary conditions. The grid size also varies in the vertical, with highest resolution Dz ¼ 0.5 m in the [0 m, 30 m] interval that steadily decreases to Dz ¼ 5 m below. The number of grid points is 4,663 Â 103 grid points. The time step Dt is adjusted at every iteration to be a third of the Courant-Friedrichs-Lewy condition. Two other simulations were carried out to test the sensitivity of the results of the control run (run #1 in Table 1) to either increasing or decreasing the spatial resolution in both directions by a factor of 2, and therefore correspondingly increasing/decreasing the time step by roughly a factor of 2. Qualitatively, the low-resolution simulation appears slightly more diffuse, but the main features such as the generation of large-amplitude internal solitary waves and Kelvin-Helmholtz instability are equally well reproduced. Quantitatively, the largest amplitude leading wave produced by the low-resolution run differs by 4% from the control run. The main conclusions drawn from these simulations are therefore considered to be insensitive to grid resolution. The depth-integrated perturbation energy density flux radiating away from the front in both directions at y ¼ À 1 km (' À ' subscript) and y ¼ 1 km (' þ ' subscript) is calculated as 52,53 and F þ ¼ where p 0 is the sum of the perturbation hydrostatic pressure and nonhydrostatic pressure and P is the available potential energy relative to the reference density profile taken either to the far left (for F À ) or far right (for F þ ) of the domain 54 . These energy fluxes contain not only high-frequency fluctuations due to internal waves but also lower-frequency changes due to the varying background conditions caused by the intrusion. These trends are removed with an envelope filter before integrating the flux to obtain the perturbation energy E 0 ¼ R F 0 dt where F 0 is the high-pass filtered flux (Fig. 7).
The horizontal lengthscale of an internal solitary wave is computed as 55 where a is the wave maximum vertical displacement, Z is the position of the interface, and y 1 and y 2 are chosen to include the entire wave.
Code availability. The code that supports the findings of this study, written in FORTRAN 90, is available from the corresponding author on reasonable request.
Data availability. The data that support the findings of this study are available from the corresponding author on reasonable request.