Transformation and disintegration of strongly nonlinear internal waves by topography in stratified lakes

Abstract. For many lakes the nonlinear transfer of energy from basin-scale internal waves to short-period motions, such as solitary internal waves (SIW) and wave trains, their successive interaction with lake boundaries, as well as over-turning and breaking are important mechanisms for an enhanced mixing of the turbulent benthic boundary layer. In the present paper, the evolution of plane SIWs in a variable depth channel, typical of a lake of variable depth, is considered, with the basis being the Reynolds equations. The vertical fluid stratification, wave amplitudes and bottom parameters are taken close to those observed in Lake Constance, a typical mountain lake. The problem is solved numerically. Three different scenarios of a wave evolution over variable bottom topography are examined. It is found that the basic parameter controlling the mechanism of wave evolution is the ratio of the wave amplitude to the distance from the metalimnion to the bottom d . At sites with a gentle sloping bottom, where d is small, propagating (weak or strong) internal waves adjust to the local ambient conditions and preserve their form. No secondary waves or wave trains arise during wave propagation from the deep part to the shallow water. If the amplitude of the propagating waves is comparable with the distance between the metalimnion and the top of the underwater obstacle ( d ~ 1), nonlinear dispersion plays a key role. A wave approaching the bottom feature splits into a sequence of secondary waves (solitary internal waves and an attached oscillating wave tail). The energy of the SIWs above the underwater obstacle is transmitted in parts from the first baroclinic mode, to the higher modes. Most crucially, when the internal wave propagates from the deep part of a basin to the shallow boundary, a breaking event can arise. The cumulative effects of the nonlinearity lead to a steepening and overturning of the rear wave face over the inclined bottom and to the formation of a turbulent jet propagating upslope. Some time later, after the breaking event, a new stable stratification is formed at the site of wave destruction. The breaking criterion of ISWs is discussed. Key words. Oceanography: general (limnology; numerical modeling) – Oceanography: physical (internal and inertial waves)


Introduction
Internal waves play an important role in water dynamics of stratified lakes. They intensify the momentum and energy transport, enhance the turbulence and water mixing through the water column, and in turn, facilitate the ventilation of a thermocline, as well as energize the vertical diffusion of oxygen, nitrate, hydrogen sulfide and other chemical elements between near-surface and near-bottom waters (Eckert et al., 2000). The last process has a very important influence on the environment and biological productivity of lakes and prevents many ecological disasters (Vinni et al., 2000). This is why internal waves in lakes have been studied extensively over the past century and their description has received considerable attention in the literature (for detailed reviews see, for instance, Hutter, 1984aHutter, , b, 1986Grimshaw, 1986;. It is a common belief, that the energy influx to the lakes' depth is supplied by wind acting on the free surface. It drives the surface water and generates internal waves in the form of basin-scale standing waves (Mortimer, 1952;Mortimer and Horn, 1982;Hutter, 1991) or propagating nonlinear waves (Farmer, 1978). After storms, internal waves in lakes may take the form of an internal surge or packets of internal solitons, generated by the nonlinear steepening of a basin-scale finite-amplitude wave (Hunkins and Fliegel, 1973;Thorpe et al., 1972;Thorpe, 1977;Farmer, 1978;Wiegand and Carmack, 1986). Since these solitons are much shorter in length than the wind-induced initial large-scale thermocline displacements, their generation results in a transfer of energy within the internal wave field from large to small scales.
Considerable progress has been made recently in the mathematical modeling of this mechanism of generation (e.g. Hollan and Simons, 1978;Hutter, 1983Hutter, , 1984aStoker et al., 1987;Hodges et al., 2000;. At the same time surprisingly little works were devoted to the clarification of the problem of how internal wave energy is transferred along a spectrum from large to small scales. From theoretical studies it is known that wind stresses, applied to the free surface, excite the low-frequency basinscale internal (Kelvin and Poincaré) waves which, according to the dispersion relation, have a discrete spectrum of wavelengths. In situ measurements, nevertheless, show that the internal wave field in lakes has a continuous spectrum ranging from low-frequency basin-scale to the maximum buoyancy-frequency waves (Thorpe et al., 1996;Hodges et al., 2000;Imberger, 1998, 2001). Moreover, according to experimental data, basin-scale internal waves excited by wind, attenuate remarkably quicker than predicted by models. As was found by Stevens at al. (1996) on the basis of almost 3 years of measurements at Kootenay Lake (British Columbia, Canada), over 80% of the potential energy of the basin-scale internal waves is lost in long, narrow lakes in one internal wave period. On the other hand, threedimensional numerical models of the nonlinear shallow water equations are fraught with implicit numerical diffusion that basin wide internal oscillations are damped out before they may even have developed. Considerable computational expertise is needed to overcome these difficulties Umlauf et al., 1999).
There are several reasons for the decay of internal waves in stratified lakes. One of the more probable sources is dissipation due to the viscosity in layers with large velocity shears (basically in the bottom boundary layer and metalimnion, where baroclinic modes are changed dramatically). Such sites are also potential places for the development of shear (Kelvin-Helmholtz) instabilities, which can also extract energy from the mean flow to the small-scale motions. The importance of the dissipative mechanism of internal wave attenuation in lakes was discussed, for instance, by Imberger and Ivey (1991). Estimations of the energy decay by the viscous effects performed by Monismith (1985) had shown that it is possible to lose the potential energy to heat in a time comparable with the wave period.
Another possible mechanism of the energy transfer from a background wave field to the small-scale motions can be nonlinear steepening and disintegration of long internal waves to packets of short-period waves and solitons. Manifestation of such strong nonlinearity of internal waves in lakes was reported by Hunkins and Fliegel (1973), Thorpe (1977), Farmer (1978), Wiegand and Carmack (1986).
An important role in the wave-energy transfer along the spectrum is played by the dynamical processes taking place in nearshore regions. For instance, observations (Hamblin et al., 1994) show significant transverse motions introduced by the variations of the basin width. We speculate that constrictions and the sudden widening of the lake width are likely places of degeneration of basin scale waves and transfer of its energy to small-scale motions.
Finally, all waves that are generated in lakes must encounter inevitably a variable bottom topography and interact with it simply because lakes are closed basins. So, one can expect that the basic part of the internal wave energy is transferred to the small-scale motions, turbulence and mixing in nearshore regions where the water depth is small. Under conditions of variable bottom topography, the shoaling of such waves and their transformation, as well as overturning and breaking over sloping boundaries lead to the loss of most of their energy (e.g. Helfrich andMelville, 1986, Helfrich, 1992;Michalet and Ivey, 1999;Grue et al., 2000;Hutter, 2001, 2002). In the past, such processes of lake dynamics were rather poorly theoretically investigated. This is the case because the commonly used models exploit the hydrostatic pressure approximation. They are obviously not appropriate for strong, nonhydrostatic effects, such as nonlinear steepening of internal waves, their fission, and subsequent disintegration into a sequence of solitary waves, the strong interaction of waves with bottom topography, as well as overturning, breaking and successive water mixing. These processes provide a flux path from the wind to the turbulent benthic boundary layer. To take all of these effects into account in any reasonable form, some efforts were undertaken (even in the framework of the hydrostatic approach) to parameterise the nonhydrostatic effects of the evolution, propagation and dissipation of internal solitary waves (Horn et al., 1999). This is one of the possible ways; we suggest an alternative approach, which incorporates all the effects of strong nonlinearity, as well as nonhydrostaticity.
In the present study, we apply the mathematical numerical model developed on the basis of the Reynolds equations in the Boussinesq approximation, to examine the influence of the bottom topography and nonlinearity of the wave processes on the mechanism of the energy transfer from the large-to the small-scale motions and even farther -to turbulence and mixing, as extreme situations. The model takes into account nonhydrostatic pressure, continuous fluid stratification, full nonlinearity (without implementing the commonly used assumption of weak nonlinearity) and variable bottom topography. We shall concentrate on identifying the mechanism of energy transfer from large-scale internal waves to short scales. As an example for the model application, topographic features from Lake Constance are used. Several different scenarios of the nonlinear wave evolution over the bottom topography are investigated. They indicate possible important sinks of wave energy from the large scale motions to the small scales and to turbulence. Even though they are obtained for the conditions of Lake Constance, these results stand as examples but apply equally to many other lakes.

Formulation of the model
As was mentioned above, the process of disintegration and dissipation of basin-scale internal waves includes several successive stages. Here we examine the final stage of the degeneration process at which the large-scale internal wave due to the nonlinear steepening has already been transformed into a sequence of solitary internal waves (SIW) (Thorpe et al., 1972;Hunkins and Fliegel, 1973;Farmer, 1978;Saggio and Imberger, 1998;Horn et al., 2001). We focus attention on the ultimate fate of such waves when they encounter the bottom topography. We consider shoaling and breaking of solitary waves which is an important energy sink for the internal wave field; both play a major role in the driving mixing processes in stratified lakes.
For short-period internal waves like solitons, the influence of the Earth's rotation is a secondary effect and will here be neglected (Diebels et al., 1992). Therefore, one can use a 2-D hydrodynamic model. With such prerequisites the twodimensional system of Reynolds equations in the Boussinesq approximation, written in Cartesian x-z coordinates in which the x-axis lies at the undisturbed free surface and the z-axis is vertical opposite to the direction of gravity, take the form for the stream function ψ(x, z, t) (u = ψ z , w = −ψ x , where u, v are the horizontal and vertical velocity components) and vorticity ω(x, z, t) (ω = ψ xx + ψ zz ). Here, ρ(x, z, t) is the water density, k v and k h are the coefficients of vertical and horizontal turbulent diffusion, a v and a h are the coefficients of vertical and horizontal eddy viscosity, J (a, b) = a x b z − a z b x is the Jacobian operator and the subscripted indices x, z, t denote partial derivatives with respect to the subscripted variable. There is very little known about the appropriate parameterization of diffusivities of momentum and energy in lakes (see, however, e.g. Hutter, 1986). We assumed one of the popular parametrizations of the coefficients of vertical turbulent exchange a v and k v Richardson-number dependences (Pakanowski and Philander, 1981), which was successfully used in different applications for the oceanic problems, viz., where the Richardson number is defined as Here, a b and k b are dissipation parameters describing background turbulence, a o , α and p are adjustable parameters and N(x, z, t) is the local Brunt-Väisälä frequency (N 2 (x, z, t) = −gρ z /ρ). These parameterisations for the vertical turbulent kinematic viscosity and diffusivity increase the coefficients a v and k v in regions with small values of Ri.
In regions with vertical density inversions, the Richardson number was taken to be 0. For the coefficients of horizontal turbulent exchange a h , k h we applied the following gradient-dependent parameterisation: in which β 0 and µ 0 are homogeneous background values, β and µ are disposal parameters and x is the horizontal step size of the grid. For more details about this, we refer the reader to the original paper by Stacey and Zedel (1986). We are only interested in baroclinic motions, and thus use the "rigid-lid" condition at z=0. In this case, at a free surface, boundary conditions can be selected in the following form: (2) The first condition says that the free surface is a streamline, the last two express the vanishing of the shear stresses and zero mass flux. The bottom line, z = −H (x), is a streamline at which the "no-slip" condition and zero mass flux across the boundary are imposed. This yields where n is the unit normal vector to the irregular bottom relief. System (1) with boundary conditions (2), (3) was solved numerically. It was, however, useful to first perform a σtransformation of coordinates (σ = −z/H (x)), which transforms the physical domain with variable depth into a rectangular computational area. This allowed us to use a rectangular grid in the development of a numerical scheme.
The splitting-up method (Marchuk, 1974) was applied for the finite-difference approximation of the equations. The vorticity transport equation is integrated in time by splitting the temporal step into two semi-steps. The spatial derivatives are approximated by second order central finite differences. At each temporal semi-step an implicit system of equations with a tridiagonal matrix is obtained that is solved using standard techniques. The stream function is then computed from the vorticity by solving the Poisson equation ω = ψ xx + ψ zz , previously transformed with the use of the σ -coordinate. For its solution a standard relaxation method is applied. Finally, the density ρ is computed by using the same methodology as for the vorticity.
In situ measurements and laboratory experiments suggest a range of values for the adjustable parameters in the a v , k v , a h and k h parameterizations. From a numerical point of view the optimal parameter selection is one that simulates the observations best and provides the stability of the numerical scheme during wave breaking.
Several numerical experiments were performed to determine the values of the eddy viscosity and turbulent diffusivity parameters a o , a b , k b , α, p, β 0 , β, µ 0 , µ. Of course, the intensity of the background mixing essentially depends on the scales of the considered oceanic processes. They can not be the same for the short-scale processes (like those considered here in the SIW) as for the large-scale oceanic circulation considered by Jones (1973) or Pakanowski and Philander (1981). In the present study, the background mixing coefficients a b , k b are as small as in the above cited papers, but for a o we took a 100 times smaller value.
It was found that by taking t = 0.5s, x = 2 m, z = 0.5 m (in the deep part of the basin), a o = 2 · 10 −4 m 2 s −1 , a b = 10 −5 m 2 s −1 , k b = 10 −6 m 2 s −1 , p = 2, α = 4, β 0 = 10 −4 m 2 s −1 , β = 1, µ 0 = 10 −5 m 2 s −1 , µ = 1 the numerical scheme was stable even during the overturning of a wave front. At the same time the selection of these parameters provided a long distance propagation of the SIW without essential attenuation (500 wavelengths and more). This is in agreement with in situ observations, which show that SIW can propagate several hundreds of kilometres from the source of generation without essential attenuation (Osborne and Burch, 1980;Apel et al., 1985). The essential (10 times) increase in the coefficients of the eddy viscosity and turbulent diffusivity does not introduce in the model results with any remarkable distinction on the initial stage of breaking, but after the breaking, the spot of mixed water at large coefficients a v , k v , a h and k h was more regular.

Initialization of the model
The model is applied to the study of the interaction of solitary internal waves with a variable bottom topography. It is assumed that a plane wave of depression with amplitude a m , whose vertical structure is defined by the density profiles, propagates from the deep part of the lake to an underwater ridge or a slope region near the shore. The objective of our modeling efforts was to investigate how the process of evolution of a solitary internal wave reacts to the variable bottom topography. Moreover, an interest is in a possible classification of bathymetric profiles into types of distinct behaviour. To be closer to real lake conditions, we selected a region of the Mainau Island in Lake Constance where a sill acts as a constriction, (the bottom topography of the region under study is shown in Fig. 1). Since Lake Constance is elongated in one direction, one can assume that the lake's axis is the preferable direction of internal wave propagation. So, the idea of this study is to investigate the evolution of an internal wave propagating from the Obersee to theÜberlinger See through the constriction near the Mainau Island. Such waves can encounter three qualitatively different typical profiles of the bottom topography. For the cross sections 1-3 depicted in Fig. 1 they are presented in Fig. 2.
Profile 1 (central part of the lake) is characterised by a gentle slope and smooth transition from the deep part of the Obersee to the Sill of Mainau (∼ 80 m deep). Internal waves propagating more closely to the northern coast of the lake (profile 2) encounter (in Sect. A-A) a higher sill (with its top at about 50 m below the surface), followed by a relatively strong increase in water depth.
A qualitatively different profile is encountered by waves propagating in direction 3 (we should also take into account this possibility because almost 40% of the cross section B-B can be described by such a bottom profile). This slope is shorter and steeper in comparison to profiles 1 and 2. Moreover, contrary to the two first counterparts, it ends on the shore. Presumably, the wave evolution here will essentially differ from the dynamical processes in the other two places. 1 The fluid stratification used in this study is shown in Fig. 3 (Ollinger, 1999). Two profiles of the buoyancy frequency presented here reflect typical vertical distributions for the summer and autumn seasons. In summer, the density jump is more pronounced and lies above the similar autumn maximum. As a consequence, the maximum of the wave displacements in this season is located closer to the free surface than  (2) (Ollinger, 1999).
in autumn. The last circumstance leads, in turn, to a conspicuous difference between the wave behaviour over the bottom features in different seasons, as we shall soon see. Amplitudes of internal waves measured in Lake Constance usually lie in the range of 0÷15 m. To recognise the influence of the nonlinearity of a wave process on the mechanism of internal wave transformation over variable bottom topography (dependence of the wave dynamics on the wave amplitudes), four different incoming solitary internal waves of depression were taken as initial conditions for every season: • a m = 2 m -weakly nonlinear case; • a m = 5 m -intermediate nonlinearity; • a m = 10; 15 m -strong nonlinearity.
The approach used here to obtain the initial fields (at t = 0) for the incident wave in this paper was the same as that used by Vlasenko et al. (2000). At the initial stage we considered a basin with constant depth of 200 m (deep part of Lake Constance). The numerical model was initialized by using the first-mode analytical solitary wave solution of the stratified Korteweg-de Vries (K-dV) equation. Such initial fields represent a stationary solitary wave solution in a weakly nonlinear, nonhydrostatic medium, but they do not satisfy system (1) for large amplitude solitary waves, for which the nonlinear parameter lies outside the region of applicability of weakly nonlinear models. Once inserted in the numerical scheme the strong, nonlinear wave will evolve in the basin of constant depth. During this evolutional process the initial large amplitude K-dV soliton will be modified, and a new stationary solitary wave will be formed at the frontal side of the wave field. The leading wave with larger phase speed than the wave tail separates from the latter at a definite stage of evolution (at the distance of 20÷30 wavelengths from its incept), and it propagates independently farther as a solitary wave. Thus, the model is run until a leading wave separates from the wave tail, and until a new stationary solution is reached at the frontal side of the wave field. This wave is then used as an initial condition for the problem of the interaction of the intense ISW with the bottom topography.
Numerical runs were performed for every SIW listed above and every profile shown in Figs. 1 and 2, for both the summer and autumn seasons. Below is the description of the three different scenarios of a wave evolution in the considered area, which are found for different bottom profiles.

SCENARIO 1: Wave adjustment
This scenario was found for profile 1 crossing the Sill of Mainau at its centre. As was already mentioned above, this profile is characterised by a gentle slope (by a weak change of the bottom depth from 200 m in the deep part of the Obersee to ∼80 m depth at the Sill of Mainau). Taking into account that for the summer and for the autumn stratifications the maximum of the SIW amplitude lies very far from the bottom (at the horizons 20 m and 35 m, respectively), internal waves under conditions of a weakly variable bottom are not destroyed and do not transform into a secondary wave packet during their evolution. These SIWs (weak, or strong) simply adjust permanently and passively to the local ambient conditions (stratification and depth) by practically preserving their initial form. One example of such an adjustment (for autumn stratification) is presented in Fig. 4. Its panels show the density and u and w velocity fields for incident waves with an amplitude of 15 m in the deep part of the basin. The right panels represent the result of adjustment of the solitary wave to the conditions of the Sill of Mainau (depth of 80 m). No secondary SIWs or wave trains arise during wave propagation from the deep to the shallower part of the water. The incoming wave preserves its initial form as a solitary wave of deflection. The single evident qualitative change only concerns the nonessential transformation of the vertical structure of the horizontal and vertical velocities, according to the vertical profile of the first baroclinic mode in shallow water. Quantitatively, the maxima of the u and w velocities become more pronounced as a consequence of the conservation of energy.
Note that a possible reflection which can take place during the interaction of the SIWs with the bottom topography does not take place due to the small inclination of the bottom topography (the bottom angle is less than 0.6 • ). Thus, backscattering of the wave is negligible, and (except for the  losses due to friction) the wave conserves its energy during propagation. So, it can be concluded that all characteristics of the medium and wave parameters change very slowly along profile 1, and therefore, the approximation of geometric optics (or the WKB approximation) can be used in this case, as well as an approach based on the Korteweg de-Vries equation.

SCENARIO 2: Wave transformation
Propagation of SIWs along Sect. 2, which crosses the local sill A-A (see Figs. 1 and 2) leads to another scenario of the wave evolution. This more representative example of a wave transformation over the sill A-A is shown in Fig. 5. Autumn stratification was used (see Fig. 3) for which the density jump is located at about 20 m depth. So, the distance between the pycnocline axis and the top of the sill was about 30 m. At the same time, the maximum of the isopycnal deflection in the centre of the incoming wave for such stratification (defined from the standard boundary value problem for internal waves) does not lie exactly in the pycnocline centre but below it. In our case it was the horizon z = −35 m. Thus, the distance between the maximum of the eigenfunction and the top of the sill equals 15 m. Taking into account that the amplitude of the incoming wave (a m = 15 m) is comparable with this distance, one should expect very strong effects of interaction. Figure  • the first stage (prior to t = T o ) is characterised by a weak transformation of the propagating SIW and its adjustment to the ambient conditions over a slowly variable depth. No dispersive effects are visible. This stage is analogous to that described above; • during the time interval T o < t < T o + 25T , the SIW continues to propagate over a slightly inclined bottom, but its vertical scope becomes comparable with the total water depth. This leads to an increase in the nonlinear effects. First, it concerns the nonlinear phase speed, which, in fact, strongly depends on the wave amplitude. The wave trough begins to propagate more slowly in comparison to the wave periphery. This effect leads to a steepening of the rear wave face; however, the frontal face becomes more and more gently sloping (see Fig. 5b). Most likely, the SIW would be destroyed further, if the water depth would continue to decrease, but beyond the point "Top" (see Fig. 5c) the SIW leaves the shallow water zone and propagates further into relatively deep water. At this moment the "dispersive" stage of the wave evolution commences.
• for t > T o + 30T the rear steep face of the SIW, due to dispersion, starts to transform and evolve into a secondary oscillating wave tail. The result of such transformation is presented in Figs. 5f and 6. One should note, that the energy of the initial wave is transferred not only to the first baroclinic mode, but also to the higher modes. Fragment D in Fig. 5f with counterphase displacements of isopycnals represents second baroclinic mode waves. Such an effect of energy transfer to a high baroclinic mode was found experimentally in Hüttemann and Hutter (2001) and described theoretically by Vlasenko and Hutter (2001). Due to the smaller phase speed the second mode wave packets (fragment D) separates effectively from the leading first mode wave train (fragment C). The reason of the generation of the second-and higher baroclinic modes is the sudden change of the bottom topography on the top of the sill (sharp edge). The initial stage of this second mode generation is marked by the dashed line E-E in Fig. 5d.
A last remark on the "dispersive" stage of evolution: on the frontal side of the wave field, instead of the single in-coming SIW, two newly-born solitary waves of depression are formed (positions A and B in Figs. 5 and 6). Their characteristics are very close to those of the Korteweg de-Vries solitons. Soliton B is weaker than A, which is why it is wider and propagates more slowly than soliton A.
In the above described mechanism of internal wave evolution, wave reflection induced by bottom topography does not play any essential role. This is understandable if one compares the scales of appreciable bottom changes (which is more than 10 km) with typical wavelengths, which are much smaller (λ ∼150 m). The wave "behaves" locally. Thus, the key role of disintegration and fission of an incoming SIWs is played by the nonlinearity of the wave process. The next example corroborate this statement. Figure 7, like Fig. 6, also represents the result of the penetration of a SIW onto a sill (behind the Sect. A-A, as shown in Fig. 1). All conditions in this run are kept the same as for the case discussed above, with one exclusion: the amplitude of the incoming wave is 3 times smaller (5 m instead of 15 m). Therefore the fission of an incoming wave into sec- ondary waves over the sill's top and behind is much weaker (compare Figs. 6 and 7). This is not surprising since in the last case the nonlinear parameter is 3 times smaller. Under such conditions the weak propagating SIW does not change so dramatically over the top of the sill, as it does for a strong wave (see Fig. 5), and it is not destroyed. As a result, the backscattered wave tail here is much weaker. This fact of strong dependence of the efficiency of wave scattering upon the nonlinearity was checked for every wave listed above both for summer as well as for autumn stratification. It was found that for very weak waves (a m = 2 m), the energy transfer from the incident wave to the wave tail does not exceed 1%.

SCENARIO 3: Wave breaking
The last series of numerical runs was performed for Sect. 3 of Fig. 1; it is characterised by steep bottom topography and very shallow water (5 m or less) near the Mainau Island. The bottom profile used in the model is depicted by the dashed line in Fig. 2.
Such a situation prevails when the pycnocline in the deep part of the basin is located below the top of the bottom topography (isopycnals run against the bottom). A very strong effect of wave-bottom interaction must be expected during the shoaling of the internal wave. An incoming wave can be completely destroyed after its overturning and breaking in a shallow water area. An example of such a destruction is presented in Fig. 8. The referenced time t = T b is defined hereafter as the moment of wave breaking when the horizontal velocity at any point exceeds the local phase speed.
Due to the nonlinear effects, the frontal face of the propagating wave becomes more gently sloping and, correspondingly, the rear face becomes steeper (Fig. 8a). Such a change in the wave profile occurs due to a smaller propagation velocity of a trough in comparison to the wave periphery and was discussed in the previous section, when the evolution of SIWs propagating along the cross section 2 was discussed. The basic difference compared to those cases is that here, the incident wave does not penetrate into a deep-water region after passing the top of the sill, but continues to prop- agate under the conditions of a narrowing waveguide. This circumstance leads to a progressive steepening of the rear wave face, and at a definite stage of evolution, the rear wave front approaches a vertical wall (baroclinic bore) and is subsequently destroyed. Figure 8 illustrates the process of overturning of a rear steep front; it is shown at that moment of time when the top of the propagating baroclinic bore outstrips the wave trough. In this situation, the heavier and denser water penetrates into the relatively light water layers and falls down to the wave trough (Fig. 8c). Thus, we may conclude that the reason for wave breaking is the kinematic instability of the wave profile. This conclusion can be seen even more clearly in Fig. 9, where the horizontal and vertical velocities during of an overturning event are presented (the left panels show the fields of density and u and w velocities before overturning and the right panels -after overturning). Before breaking, the maximum of the horizontal velocity (30.1 cm s −1 ) was located at the free surface. During the overturning event, this maximum shifts from the free surface to the deep layers and to the lateral boundary of the steep slope, where the value of the maximum velocity increases right up to wave breaking. At the moment of breaking, the onshore fluid velocities (42.1 cm s −1 ) exceed the wave phase speed (36.4 cm s −1 ).
The overturning event is also clearly seen from the analysis of the evolution of the w-field (note the appearance of The three lowermost panels in Fig. 8 represent the ultimate fate of a solitary wave after its breaking. As seen from the evolution of the density field, the water does not mix instantaneously at the location of the overturned rear face. The dense water, penetrating into the light layers, behaves as an upstream propagating jet. In Fig. 10, this jet is presented in greater detail. It is characterized by the upslope directed core, which is located as an intermediate layer between two boundary layers with opposite velocities and which is exposed by vertical oscillations. It is obvious that zones of internal wave breaking are the potential places of enhanced energy dissipation and intensive water mixing. Local coefficients of turbulent viscosity and diffusivity in such places increase dramatically (by several orders of magnitude). This circumstance leads to a fast attenuation, dissipation and overmixing of water layers and to the formation of a local zone with a new stable vertical fluid stratification. For instance, the horizontal pulsating stream described above and presented in Fig. 8e, after the time interval t = 25T , transforms into a more regular structure, almost without vertical density inversions (see fragment A-A in panel (f) of Fig. 8). The structure of the density field, together with the patterns of horizontal and vertical speeds at t = T b + 35T , is shown in more detail in Fig. 11.
Comparison of Figs. 10 and 11 obviously shows the reconstruction of the hydrophysical fields during 2075 s (25 T ): the density field becomes more regular, the maximum horizontal velocity decreases substantially (1.5-2 times) and the vertical pulsations almost completely disappear (vertical circulations remain only on the frontal side of the moving intrusion and on its rear side). The average velocity of propagation of the density intrusion between t = T b + 25T and t = T b + 35T is 18 cm s −1 .
The above results concern the evolution and destruction of a relatively large wave. Similar effects of overturning and breaking also took place for all other considered waves. As for the basic features, the mechanisms of their evolution are very similar and exhibit the same characteristic stages: nonlinear transformation, formation of a baroclinic bore, its overturning and breaking with the formation of density intrusion and its dissipation. Deviation from such a scenario is only observed for small-amplitude waves. Let us consider these differences in more detail. Figure 12 shows the evolution of a 2 m amplitude SIW in the slope-shelf area for summer stratification. Qualitatively, the processes of transformation of the small-amplitude wave at initial stages look very similar to those which were described for strong waves. At the beginning of the wave transformation, a SIW penetrates into a shallow-water zone without any visible changes in its form up to the depths where the pycnocline crosses the bottom surface (Figs. 12a and b). The  (Figs. 12c). Position A indicates the local isopycnal elevation which was formed due to the nonlinearity of the wave process.
Position B shows the place of internal wave breaking. Note that the breaking process of the relatively weak waves is not so pronounced as it is for large-amplitude waves. A relatively weak nonlinearity leads to the result that a smallamplitude wave does not overturn over the whole depth, but only near the bottom. The local zone of overturning is marked by an ellipse in Fig. 12d. The near-surface isolines do not reveal a tendency for overturning to occur. This last circumstance leads to a slightly different scenario of wave breaking, which is valid in the case of a weak nonlinearity: the internal wave overturns and is destroyed only near the bottom, where a bottom boundary layer (BBL) exists; the turbulent pulsating spot of overmixed water, which arises after overturning, disappears rapidly due to the high level of turbulent viscosity and diffusivity in the BBL. Beyond this the wave motion does not have stochastic character, but rather it is deterministic.
At the same time the near-surface layers of the SIW, which were not exposed to overturning, form a wave of elevation, which propagates upstream. Fragments of such waves are shown in the last two panels in Fig. 12 by the letters C and D. So, in the case of weakly-nonlinear waves rather than a turbulent pulsating jet, waves of elevation (or "boluses") are produced in the breaking area. This conclusion is in agreement with results of numerous laboratory experiments on the breaking of internal waves (e.g. Helfrich and Melville, 1986;Helfrich, 1992). The last remark concerning the difference between the process of wave breaking for strong and weak waves has the consequence for the breaking criterion. Figure 13 represents the breaking criterion, which was found by Vlasenko and Hutter (2002) for oceanic conditions (stratification, bottom inclinations and wave amplitudes). In this figure, the normalized amplitude a m /(H b − H m ) is plotted against the slope angle γ . The parameter H b is the water depth at the point of wave breaking (see inset in Fig. 13). The location of wave breaking, or the break point, is defined as the position where the orbital velocities begin to surpass the phase speed for the first time. Another basic parameter controlling the process of wave breaking is the ratio of the wave amplitude a m to H − H m , which is the distance from the undisturbed isopycnal of maximum amplitude (located at depth H m ) to the bottom. Thus, we define the value H m as an undisturbed position of the isopycnal line which has maximum depression a m at the centre of the wave (see inset in Fig. 13).
The following criterion of wave breaking was found in (Vlasenko and Hutter, 2002): where γ is given in degrees. This criterion was obtained theoretically in the framework of the approach, which is very similar to that described in the present paper. A study of the transformation of large-amplitude SIW of permanent form over a slope-shelf topography was considered for the conditions observed in the Andaman and Sulu Seas. Over the range of examined parameters (0.52 • < γ < 21.8 • and 10m < a m < 80m), a breaking event arises whenever the condition (4) is satisfied. This condition gives the dependence of the depth of the breaking point on the wave amplitude and bottom steepness. The breaking events obtained in the present series of numerical runs are presented in Fig. 13. One should note the good agreement of results obtained for the lake with similar data found for oceanic conditions for relatively large waves (a m > 5 m). Points of breaking events for smallamplitude waves (a m ≤ 5 m) lie above the breaking curve. This result can be connected with the different manifestations of the breaking mechanism for strong and weak waves. Some remarks about this difference were mentioned above. Of course, a definite discrepancy can also follow from the fact that the breaking criterion (4) was found for a flat inclined linear bottom topography, whereas the bottom profile 3 (Fig. 2) is not uniform and has a whole range of characteristic angles of inclination. For instance, in autumn, the SIW with amplitude of 15 m is destroyed farther away from the coastline than the other waves (with bottom inclination angle of 1.5 • instead of 3.4 • , see Fig. 13). Nevertheless, one should state that the breaking criterion that is found for oceanic conditions can also be applied for the lakes if internal waves are relatively strong (a m > 5 m).

Discussion and conclusions
Intensive winds applied at the free surface of a stratified lake usually produce a baroclinic response (topographic waves and internal seiches). This mechanism of energizing the water column is likely the basic source of energy for baroclinic motions of many lakes. Different mechanisms are responsible for the wave energy that is transformed into a spectrum from large-to short-scale motions. This transfer enhances the water mixing, and the latter has an important influence on the lakes' environment (water quality, biological productivity, distribution of pollutant and nutrient elements).
It is known that due to dissipation in the bottom boundary layer, excited long, internal waves can gradually lose their energy. At the same time, surprisingly little is known about the cascade mechanism of the nonlinear energy transfer from large-to small-scale motions. A number of experimental data obtained from laboratory and in situ measurements reveal the importance of the nonlinear effects for the degeneration of basin-scale internal waves. An important role in this mechanism is played by the nonlinear steepening of the basin-scale waves, their disintegration and fission on a sequence of solitary waves and successive dissipation during their shoaling and breaking over inclined lake boundaries.
Unfortunately, the most popular of the existing mathematical models describing internal wave dynamics cannot be used for the investigation of such processes. These models usually assume the hydrostatic approximation for the pressure and cannot be applied for the study of nonhydrostatic and strongly nonlinear effects, such as nonlinear steepening of internal waves, their fission and disintegration into solitary waves, as well as the interaction with the bottom topography, and the overturning and breaking. This paper has been devoted to the investigation of such strongly nonlinear effects which occur during the propagation of internal waves in stratified lakes. For this purpose, we developed a numerical model, based on the Reynolds equations. The model is nonlinear, nonhydrostatic and accounts for the real vertical fluid stratification, variable coefficients of vertical and horizontal turbulent exchange and variable bottom topography. Due to computational limitations, the model is restricted to two spatial dimensions; it thus excludes the effects of the rotation of the Earth and other topographic effects of three dimensionality. The basic goal of the numerical runs was to consider the interaction of solitary internal waves of different amplitudes with the bottom topography and to examine different possible scenarios of their nonlinear evolution. As an example for the model application, topographic profiles from Lake Constance were used with summer and autumn stratifications that are typical for this lake. We concentrated on identifying the mechanism of the energy transfer from solitary internal waves into short-scale waves and mixing, depending on the wave amplitude and bottom relief. Three different shapes of the bottom topography were selected for investigation: • the gently sloping transition zone between deep and shallow basins, which can model the transition between two separate parts of many lakes; • the underwater ridge, whose top is located 50 m beneath a free surface; • the steep bottom topography that represents a sloping boundary of a lake.
Three different scenarios of the nonlinear wave evolution over the described bottom features were examined. They can, respectively, be called "Wave adjustment", "Wave transformation" and "Wave breaking".
Wave adjustment is peculiar to sites with a gentle sloping bottom where propagating waves (weak or strong) adjust permanently to the local ambient conditions (stratification and bottom), thereby preserving their form. No secondary waves or wave trains arise during the wave propagation from the deep part to the shallow water. Wave backscattering is also negligible due to the small inclination of the bottom topography, so the approximation of geometric optics or weakly nonlinear theories, such as that based on the Korteweg de-Vries equation, are appropriate in such cases.
Wave transformation takes place in sites with underwater ridges and banks in situations when the amplitude of the incoming waves is comparable to the distance between the metalimnion and the top of the ridge, so that the parameter of nonlinearity is about unity. As a consequence, nonlinear dispersion starts to play a key role in the mechanism of the wave evolution. Due to the strong nonlinearity, the frontal face of the propagating wave becomes more gently sloping, but the rear face becomes steeper. The basic reason for such a transformation is that the propagation velocity of a wave's trough is smaller in comparison to the speed of the wave periphery. The increase in the rear face steepness can lead to overturning and breaking of the internal wave, according to "the wave breaking" mechanism. The conditions of their origin are described below. If the internal wave is not destroyed above the ridge top and penetrates into the upcoming deep water zone, it splits behind the obstacle into a sequence of secondary waves. The leading waves of a scattered wave field are the solitary internal waves arranged by the amplitude scales, followed by a weak oscillating wave tail. It was also found that the energy of the incoming wave at the position of the underwater obstacle is transmitted not only to the first baroclinic mode, but also to higher modes.
Wave breaking is the more crucial case (and probably most typical for many lakes), when the internal wave propagates from the deep part of the basin into a shallow water zone, where a breaking event of the internal wave can arise. The cumulative effects of the nonlinearity lead to a steepening and overturning of the rear wave face over the inclined bottom. Just before breaking, the horizontal orbital velocity at the site of instability exceeds the phase speed of the ISW. The top of the propagating baroclinic bore outstrips the wave trough, and at this site the heavier and denser water penetrates into the relatively light water layers and falls down to the wave trough.
In the following stages of evolution of the overturned hydraulic jump, the water does not mix instantaneously at the place of the overturned rear wave face. The dense water behaves as an upstream propagating pulsating jet, penetrating into the light layers. It is characterized by the upslope directed core which is located in the intermediate layers and which is exposed by vertical oscillations.
Zones of internal wave breaking are the potential places of enhanced energy dissipation and intensive water mixing. Incident internal waves are completely destroyed in such zones during the breaking event and give rise to a new stable vertical stratification at the site of destruction.
The breaking criterion of ISWs over the sloping bottom (4), which was found by Vlasenko and Hutter (2002) was also tested for lakes conditions. Good agreement was obtained for relatively large waves (a m > 5 m). The mechanisms of overturning and breaking of small-amplitude waves (a m < 5 m) into basic features goes through the same characteristic stages (nonlinear transformation, formation of baroclinic bore, its overturning and breaking, formation of density intrusion and its dissipation). These results must be contrasted with scenarios for small-amplitude waves, where the breaking processes of relatively weak waves are not so pronounced as those for large-amplitude waves. A relatively weak nonlinearity leads to the fact that a small-amplitude wave does not overturn through the whole depth, but rather only near the bottom. The near-surface isolines do not reveal a tendency for overturning to occur.
The turbulent pulsating spot of overmixed water, which arises after overturning the weakly nonlinear waves, disappears rapidly due to a high level of turbulent viscosity and diffusivity in the bottom boundary layer. At the same time, the near-surface fragments of the SIW, which were not exposed to overturning, form a wave of elevation that propagates upstream. Thus, in the case of weakly-nonlinear waves, instead of a turbulent pulsating jet, waves of elevation (or "boluses") are produced in the breaking area. This conclusion is in agreement with results of laboratory experiments by Helfrich andMellville (1986) andHelfrich (1992).
So, the results presented in this paper indicate possible important sinks of wave energy from the large-scale motions to small-scales and to turbulence. Obtained for conditions adequate for Lake Constance, they can also be applied to other lakes. We highlighted only a few possible scenarios of wave evolution in lakes and tried to indicate the parameters controlling the process of wave evolution. Of course, to some extent, such a classification is a little bit conventional, and the real manifestation of internal wave dynamics can be more complicated and can also incorporate 3-D effects like refraction, influence of variable lake width or wave-wave interaction. Their influence on wave dynamics should also be studied, which is why we plan to expand our model in this direction for long, narrow lakes. Future studies will also help us to answer some other interesting and important questions related to the water mixing in lakes due to the breaking of internal waves (for instance, multiple wave breaking, distribution of energy dissipation, generation of secondary internal waves from the spots of mixed water and so on). Nevertheless, we really believe that the results obtained above can help to understand the basic features of nonlinear energy transfer from large-to small-scales in lakes and also can be useful in many practical applications.