In a fluid with uniform density stratification the internal waves obey a specific dispersion relation [1] which connects the wave frequency and the angle of inclination of the wave beam with respect to the vector of gravity but does not contain the length scale. The consequence of the dispersion relation is the specific law of reflection from the boundaries under which the wave-ray angle with respect to the vertical, but not with respect to the normal to the surface, is conserved during reflection. In its turn, such a reflection law leads to the focusing of the wave beams due to reflection at inclined wall [2, 3]. In the presence of inclined walls in the closed volume of stratified fluid the wave focusing can lead to energy concentration along closed trajectories, the so-called wave attractors [4, 5]. Numerous studies reviewed in [3, 6] are devoted to the wave attractors. In geophysical applications, tidal currents are considered as the energy sources that generate the internal wave motions [7, 8]. The interaction of the surface tides with the bottom topography leads to generation of internal tides; in this case, under certain conditions, the internal tidal motion can exist in form of attractors [9–12]. It is well known that the surface tides have a fairly complex spectrum and well-studied deterministic nature. The full-scale measurements show that the internal tides are distinguished by fairly lesser predictability. This results from both the high sensitivity of internal waves to the fine structure of stratification [13] and the development of instability [14, 15].

In order to investigate specifics of dynamics of the regimes with wave attractors, special methods of processing and visualizing the data, both experimental and numerical, were used. These methods include the time-and- frequency diagrams, the bi-spectra, the bi-correlations, and the Hilbert transforms with filtering the frequency range and the wave vectors [16, 17]. The first attempts of two-dimensional numerical simulation demonstrated qualitative correspondence of the low-amplitude regimes of wave attractors [18, 19]. The investigation of three-dimensional effect of wave attractors carried out in [20, 21] showed the adequate qualitative and quantitative correspondence of the results of direct numerical simulation and experiments and the difference between the results of previous two- and three-dimensional numerical calculations was explained. It should be noted that almost all the literature in the field of simulation of internal tides and wave attractors is devoted to the investigation of monochromatic external forcing [3, 6]. In the present study, the linear and nonlinear dynamics of stratified liquid are investigated numerically in a model system which admits the existence of internal wave attractors under the impact of a biharmonic perturbation.

1 MATHEMATICAL FORMULATION OF THE PROBLEM

We will consider a tank of trapezoidal form with the greater base at top filled with a uniformly stratified liquid in the gravity field (Fig. 1). This geometric structure is fairly well studied in literature, namely, there is a detailed classification of the geometric configurations of attractors within the framework of the ray path theory [5], the stream functions are constructed within the framework of the inviscid liquid model [22], the linear [18] and nonlinear [23] scaling for the width of attractor’s wave beams is investigated, and the energy cascade in the nonlinear regime is described [24]. The stratification can be characterized by the salinity gradient which is codirected with the vector of acceleration of gravity. In the initial state the buoyancy frequency

$$N(y) = \sqrt { - \frac{g}{{\rho (y)}}\frac{{d\rho (y)}}{{dy}}} $$

is constant over the entire tank. The no-slip conditions for the velocity \({{v}_{n}} = 0\) and \({{v}_{\tau }} = 0\) are fulfilled on all the boundaries (except for the upper boundary). For the salinity the insulation condition \(\partial {{\rho }_{s}}{\text{/}}\partial n = 0\) must be satisfied on all the boundaries. The vertical velocity component is specified on the upper boundary of the tank. This simulates the wave generator (wave maker) which performs harmonic oscillations in accordance with the following expressions:

Fig. 1.
figure 1

Diagram of the computational domain (a); field of the vertical velocity component under the monochromatic impact at the amplitude a = 0.02 cm (b); corresponding graph of the vertical velocity in the middle of the first attractor ray as a function of time (c); kinetic energy as a function of time (d).

Full size image

for monochromatic oscillations \(y(0,x,t) = a\sin(2\pi x{\text{/}}{{L}_{1}})\cos({{\omega }_{0}}t)\);

for biharmonic oscillations \(y(0,x,t) = {{a}_{1}}\sin(2\pi x{\text{/}}{{L}_{1}})\cos({{\omega }_{1}}t)\) + \({{a}_{2}}\sin(2\pi x{\text{/}}{{L}_{1}})\cos({{\omega }_{2}}t)\) .

In what follows, we will consider the case \({{a}_{1}} = {{a}_{2}} = a\) as a model example. The small density drop over the tank depth \(H\) makes it possible to use the Navier–Stokes equation in the Boussinesq approximation

$$\frac{{\partial \vec {V}}}{{\partial t}} + (\vec {V} \cdot \nabla )\vec {V} = - \frac{1}{{{{\rho }_{m}}}}\nabla \hat {p} + \nu \Delta \vec {V} + \vec {f},$$
(1.1)
$$\frac{{\partial {{\rho }_{s}}}}{{\partial t}} + \vec {V} \cdot \nabla {{\rho }_{s}} = \nabla \cdot \frac{\nu }{{{\text{Sc}}}}(\nabla {{\rho }_{s}}),$$
(1.2)
$$\nabla \cdot \vec {V} = 0.$$
(1.3)

Here, \(\vec {V}\) is the velocity vector with the components \(\{ {{v}_{x}},{{v}_{y}}\} \), \(\nu \) is the kinematic viscosity, \({{\rho }_{m}}\) is the density on the upper boundary, \({{\rho }_{s}}\) is a contribution to the density due to salinity, the reduced pressure \(\hat {p}\) is the difference between the total and hydrostatic pressures for fresh water (at the density \({{\rho }_{m}}\)), \(\vec {f} = {{\rho }_{s}}{\text{/}}{{\rho }_{m}}\vec {g}\), \({\text{Sc}} = \nu {\text{/}}\kappa \) is the Shmidt number, where κ is the diffusion coefficient.

In the linear regime the characteristic wave beam width is developed by the system as a result of balance between the effects of geometric focusing in reflection of the wave beam from the inclined wall and viscosity [18, 23]. As the parameters characterizing the forcing, we can introduce the dimensionless frequency \(\omega {\text{/}}N\) and the dimensionless amplitude \(a{\text{/}}H\). In the numerical experiments the characteristic dimensionless amplitude is of the order of O(10–3), while the frequency \(\omega {\text{/}}N\) is specified over the range of existence of the attractor whose beam “skeleton” has the shape of a parallelogram with a single point of reflection on each boundary. Following [5], the attractor of such a form will be called by the attractor of (1, 1) type. The degenerate forms of the (1, 1) attractor which correspond to the boundaries of its existence range are the trapezium diagonals. As a dimensionless parameter that characterizes the perturbation, we can introduce the Reynolds number determined from the maximum wave-maker velocity and the depth: \({\text{Re}} = a{{\omega }_{0}}H{\text{/}}\nu \). The present study is the continuation of the series of studies [20, 23–25] on the experimental and numerical investigation of wave attractors in the laboratory-scale setups. For the sake of convenience in carrying out the comparisons, the numerical values of the parameters are given in the dimensionless and dimensional quantities.

In the two-dimensional formulation the numerical simulation of the system (1.1)–(1.3) was implemented in [18, 26] using the control volume approach and the qualitative agreement between the form of attractor in the calculation and linear theory was obtained. Approximately twofold difference between the amplitudes of wave motions was observed in quantitative comparison with the experiment. In [20, 24] the first three-dimensional simulation was carried out with the use of the spectral element method [27] and both the quantitative and qualitative agreement with the experiment was obtained.

All the following results were obtained from the numerical simulation of the nonlinear system of the Navier–Stokes equations in the Boussinesq approximation.

2 MONOCHROMATIC PERTURBATION

In the first stage of investigation, the monochromatic regimes were simulated with the aim to estimate the effectiveness of attractor generation under a given wave-maker amplitude and variation in the perturbation amplitude over the existence range of the attractor of (1, 1) type [5]. For the geometry shown in Fig. 1 the lower and upper boundaries of the existence range of the attractor of (1, 1) type correspond to \({{\omega }_{{{\text{cr1}}}}}{\text{/}}N = 0.55\) and \({{\omega }_{{{\text{cr2}}}}}{\text{/}}N = 0.74\), respectively. In reaching these critical frequencies, the parallelogram degenerates into a trapezium diagonal. As the integral dimensional measure of the effectiveness of attractor generation under constant wave-maker amplitude and invariable shape of the tank, we took the kinetic energy of liquid integrated over the trapezium area S

$${{E}_{k}}(t) = \int\limits_S {\frac{{{{\rho }_{m}}}}{2}} [v_{y}^{2}(t) + v_{x}^{2}(t)]dS.$$

For this measure we can introduce a value averaged over a sliding time window with the use of a fairly large number of the oscillation periods \(\left\langle {{{E}_{k}}(t)} \right\rangle \) and variation with respect to the mean value calculated as follows: \(r = D\left( {{{E}_{k}}(t) - \langle {{E}_{k}}(t)\rangle } \right){\text{/}}\langle {{E}_{k}}(t)\rangle \), where \(D\left( {{{E}_{k}}(t) - \left\langle {{{E}_{k}}(t)} \right\rangle } \right)\) is the variance with respect to the mean value. The dimensionless quantities \({{\bar {E}}_{k}}(t)\) and \(\left\langle {{{{\bar {E}}}_{k}}(t)} \right\rangle \) are determined by dividing by \({{\rho }_{m}}S{{(a\omega )}^{2}}{\text{/}}2\). It is well known that the regimes of motion in attractors can be similar to both the progressive and standing waves [25]. The quantity \(r\) makes it possible to give the qualitative estimate of proximity of the regime observed to one of the limiting cases [25].

2.1 Linear Regime

In Fig. 1 we have reproduced the characteristic form of the dependences observed in the monochromatic regime at the low amplitude of the forcing for a = 0.02 cm (\(a{\text{/}}H = 5 \times {{10}^{{ - 4}}}\)) and \(\omega {\text{/}}N = 0.63\). The characteristic time of a transient to the stationary regime is of the order of 30 oscillation periods, the signal spectrum is monochromatic with high accuracy, and the oscillations of the kinetic energy with respect to the mean value have a low amplitude (\(r = 0.103\)). The beam with the maximum energy density appearing after focusing reflection from the inclined wall is notated as the first ray of attractor. In Table 1 we have given the integral parameters that characterize the linear monochromatic regimes at a given value of a/H = \(5 \times {{10}^{{ - 4}}}\) over the frequency range from \({{\omega }_{{{\text{cr1}}}}}{\text{/}}N = 0.55\) to \({{\omega }_{{{\text{cr2}}}}}{\text{/}}N = 0.74\). It can be seen that the kinetic energy of attractor is maximized at \(\omega {\text{/}}N = 0.63\) when the oscillation amplitude is kept constant. Obviously, at this value of the forcing frequency strong nonlinear effects should be expected with increase in the wave-maker oscillation amplitude. At \(\omega {\text{/}}N = 0.63\) the quantity r reaches a minimum: in the attractor the motion is presented by the progressive wave.

Table 1. Kinetic energy under the monochromatic impacts at the amplitude a = 0.02 cm
Full size table

2.2 Nonlinear Regime

In Fig. 2 we have reproduced the characteristic flow patterns and the dependences observed in the case of weakly nonlinear regime at \(\omega {\text{/}}N = 0.63\) for a = 0.05 cm (\(a{\text{/}}H = 1.25 \times {{10}^{{ - 3}}}\)). In the weakly nonlinear regime the triad resonance [15] takes place, two subharmonic daughter waves of low amplitude being generated. The time-frequency diagram shown in Fig. 2 represents the signal spectrum calculated in the sliding window and averaged over the neighborhood of a point lying in the middle of the first attractor branch. In this regime the frequency spectrum of internal waves is discrete with the predominant contribution that corresponds to the perturbation frequency \({{\omega }_{0}}\), two subharmonic frequencies \(\omega _{1}^{*} + \omega _{2}^{*} = {{\omega }_{0}}\), two superharmonic frequencies \(\omega _{1}^{{**}} = \omega _{1}^{*} + {{\omega }_{0}}\) and \(\omega _{2}^{{**}} = \omega _{2}^{*} + {{\omega }_{0}}\), and the double frequency \(2{{\omega }_{0}}\).

Fig. 2.
figure 2

Field of the vertical velocity component under the external monochromatic impact at the amplitude a/H = \(5 \times {{10}^{{ - 4}}}\) (a) and the corresponding pressure field (b), dimensionless velocity in the middle of the first attractor ray (c), dimensionless kinetic energy (d), spectrum (e), and time-and-frequency diagram (f).

Full size image

The cascade of triad interactions develops with further increase in the perturbation amplitude to a = 0.1 cm (\(a{\text{/}}H = 2.5 \times {{10}^{{ - 3}}}\)). In Fig. 3 we have reproduced the characteristic patterns of the wave fields, the spectra, and the time development of the oscillation process and the kinetic energy of the system. The discrete components that correspond to the frequencies of daughter waves appearing in the triad resonance, which are similar to the spectrum components appearing in the weakly nonlinear case (a/H = \(1.25 \times {{10}^{{ - 3}}}\)), predominate in the frequency signal spectrum. In this case the complete signal spectrum represents the superposition of the discrete and continuous spectra. The presence of the continuous spectrum testifies that the developed wave turbulence regime appears [24, 25]. In Table 2 we have given the corresponding characteristics for the kinetic energy of system in the strongly nonlinear regime. From a comparison of Tables 1 and 2 it can be seen that in the case of the developed wave turbulence regime the global dimensionless energy characteristics of the system (the average energy \(\left\langle {{{{\bar {E}}}_{k}}} \right\rangle \) and variations with respect to the mean value r) differ only slightly from the dimensionless quantities characteristic of the linear regime. A comparison of the wave patterns in the linear and nonlinear cases shows that in the second case the energy is distributed more uniformly over the domain considered, namely, the attractor branches are wider and the daughter waves occupy the entire space.

Fig. 3.
figure 3

Field of the vertical velocity component in the monochromatic regime at the amplitude a = 0.1 cm/s during formation of an attractor (a) and after a cascade of instabilities (b), velocity in the middle of the first attractor ray as a function of time (c), kinetic energy as a function of time (d), frequency spectrum of the velocity in the wave turbulence regime (e), and time-and-frequency diagram (f).

Full size image
Table 2. Kinetic energy under the monochromatic impacts at the amplitude a = 0.1 cm
Full size table

3 REGIMES WITH BIHARMONIC PERTURBATION

3.1 Linear Regime

In Fig. 4 we have reproduced the characteristic example of the wave pattern and the main qualitative and quantitative characteristics of the system in the linear case under the biharmonic external forcing for following parameters: \({{\omega }_{1}}{\text{/}}N = 0.58\), \({{\omega }_{2}}{\text{/}}N = 0.66\), and \(a = 0.02\) cm. It can be seen that the system tends to the quasi-steady-state beating regime at time-scale of the order of 40 oscillation periods. This is close to the characteristic time of transient to the process of steady-state oscillations in the monochromatic case. The peaks corresponding to the frequencies of external perturbation predominate over the frequency spectrum; there are also peaks corresponding to the frequency \(2{{\omega }_{1}}{\text{/}}N\) and the difference frequency \(({{\omega }_{2}} - {{\omega }_{1}}){\text{/}}N\); however, their value is less than the main peak by more than two orders of magnitude. The time instants corresponding to the maximum kinetic energy lag significantly behind the time instants corresponding to the maximum wave-maker oscillation amplitudes. Of importance also to note that when the system reaches the steady-state beating regime the average kinetic energy of the system excited by a biharmonic perturbation is equal to the sum of the energies of attractors exited by the monochromatic perturbations in isolation from one another \({{\bar {E}}_{k}} = 21.7 \times {{10}^{{ - 4}}} \approx {{\bar {E}}_{{k1}}} + {{\bar {E}}_{{k2}}} = \) \((8.45 + 13.2) \times {{10}^{{ - 4}}}\) = 21.65 × 10–4 erg/cm2. Thus, in the linear regime the linear superposition principle holds with high accuracy. This is also fulfilled at small frequency difference \(({{\omega }_{1}} - {{\omega }_{2}}){\text{/}}N\).

Fig. 4.
figure 4

Vertical velocity component \({{v}_{y}}\) (a) and pressure (b) in the biharmonic regime at the amplitudes a = 0.02 cm/s and the corresponding vertical velocity component in the middle of the first attractor ray (c), kinetic energy (d), spectrum (е), and time-and-frequency diagram (f). Black curves in (c, d) show the envelopes of the wave-maker oscillation amplitude.

Full size image

3.2 Nonlinear Regime

In Figs. 5 and 6 we have reproduced the examples of nonlinear dynamics of wave attractors generated by wave-maker biharmonic oscillations for \({{\omega }_{1}}{\text{/}}N = 0.66\), \({{\omega }_{2}}{\text{/}}N = 0.628\), and \(\delta \omega {\text{/}}N = 0.031\) (Fig. 5) and \({{\omega }_{1}}{\text{/}}N = 0.628\), \({{\omega }_{2}}{\text{/}}N = 0.641\), and \(\delta \omega {\text{/}}N = 0.013\) (Fig. 6). In all the cases the wave-maker oscillation amplitudes were equal to \({{a}_{1}} = {{a}_{2}} = 0.05\) cm. It can be seen that the wave motion with a complex frequency spectrum is formed in both cases, the trend to the more dense spectrum “population” being observed with decrease in the frequency mismatch \(\delta \omega \). The distinctive beating process can be seen on the graphs of the vertical velocity as a function of time, but, in addition to mean energy oscillations, nontrivial dynamics of high-frequency energy fluctuations takes place, namely, the fluctuation amplitudes can differ by an order of magnitude in the phases of growth and decrease of the envelope of the wave-maker oscillation amplitude. Thus, the periodic “bursts” of wave turbulence are distinctive in the nonlinear biharmonic regime. Such “bursts” can clearly be seen on the time-and-frequency diagrams shown in Figs. 5 and 6. In particular, on the time-and-frequency diagram given in Fig. 5 we can see that at the frequency close to the frequency of the external forcing the signal amplitude “beatings” are shifted in time with respect to the “beatings” of daughter waves. Thus, the “beatings” of the wave-maker oscillation envelope, the “beatings” of the mean kinetic energy, and the “bursts” of wave turbulence are mismatched in time. We can assume that in the natural systems there is a time mismatch between the envelope of the internal tide amplitude and the enhancement of internal wave turbulence and mixing. A preliminary investigation of the energy of attractors generated by the biharmonic perturbation shows that in the nonlinear case the mean energy of the system differs significantly from the sum of energies of the components.

Fig. 5.
figure 5

Vertical velocity component in forming (a) and establishing (b) the biharmonic attractor with the amplitudes a = 0.05 cm/s and the relative frequency difference 0.05 and the corresponding vertical velocity (c), kinetic energy (d), spectrum (e), and time-and-frequency diagram (f). Black curves in (c, d) show the envelopes of the wave-maker oscillation amplitude.

Full size image
Fig. 6.
figure 6

Vertical velocity component (a), kinetic energy (b), spectrum (c), and time-and-frequency diagram (d) in the biharmonic regime with the amplitudes a = 0.05 cm/s and the relative frequency difference 0.02. Black curves in (a) and (b) show the envelopes of the wave-maker oscillation amplitude.

Full size image

SUMMARY

Dynamics of internal wave attractors under the external biharmonic forcing are first investigated. To get the intervals of parameters of greatest interest, the regimes under monochromatic forcing were first studied, and the frequency range was determined, where the attractor generation is the most effective. The investigation of the behavior of attractors under external biharmonic forcing showed that in the linear case the superposition principle is valid, namely, the attractors generated by each of the components of the biharmonic perturbation almost do not interact between themselves. In the nonlinear case, the beating regime accompanied by the bursts of wave turbulence developed as a result of a cascade of triad interactions can be observed under the external biharmonic forcing. In this case the level of kinetic energy fluctuations in the phase of growth of the envelope of the wave-maker oscillation amplitude can be by an order of magnitude higher than the level corresponding to decrease in the wave-maker oscillation amplitude.