1. Introduction
The study of ocean waves is significant in many ways. First, waves are one of the main mechanisms of interaction between air and water. As such, they determine the effect of the ocean on the weather (both locally and globally) as well as the weather effect on the ocean motion. The impact of the ocean on the Earth’s climate cannot be overestimated. Another important field is the study of the waves themselves. This is of great value for ship navigation, coastal construction, and other fields.
It is widely accepted nowadays that [
1,
2] the wave kinetic equation based on slowly changing wave amplitude (WKB approximation), also known in physical oceanography as the radiative transfer equation [
3], statistically describes the evolution of nonlinear ocean surface waves:
where
is the wave energy spectral density,
and
are the real space vector coordinate and Fourier space wave vector,
k and
are the modulus of the wave number and the angle of wave propagation,
is deep water dispersion relation, and
g is the gravity acceleration.
,
, and
are the nonlinear wave interaction, the wind energy input, and the dissipation source terms, respectively.
For more than 40 years, significant efforts have been devoted to the construction of different forms of
and
. In the current research, we are focusing on the assessment of recently developed ST6 wind energy input
and wave energy dissipation
source terms, recently included as a package in the operational wave forecasting model WAVEWATCH III [
4].
The ST6 model history originates from 3-year experimental observations, performed in Lake George, Australia, from 1997 through 2000 by [
5,
6]. These measurements resulted in several conclusions, including that
is the nonlinear function of the spectrum. For the wave-breaking dissipation, the conclusion of the threshold existence for “inherent” wave-breaking has been made in terms of “significant” wave steepness. Additionally, the opinion regarding the existence of the “two-phase” behavior of wave dissipation has been expressed: the wave-breaking can be realized at any scale, but for the above spectral peak waves, the breaking is enhanced due to the compression of the shorter waves by the longer ones.
Conceptually, a different approach to the formulation of wind energy input and wave breaking dissipation source terms, based on self-similar properties of Equation (
1) and known as ZRP balanced source terms model, has been recently realized in [
7,
8]. It relied on the self-similar solutions of Equation (
1) with the “footprints” in the form of specific indices of power-like total energy and mean frequency dependencies on the fetch, predicted theoretically by Zakharov with co-authors [
9,
10,
11] and observed in multiple field experiments analyzed in [
12]. It was also shown in [
7] that out of five numerically tested historically known wind energy input source terms with similar wave dissipation, four exhibited complete or quasi self-similar behavior, and one failed the self-similarity test. Thus, self-similarity is a well-known ocean surface wave property, confirmed experimentally [
13] as well as numerically [
14] but generally overlooked in applications.
In the current paper, we are applying the above-mentioned self-similarity testing technology [
7] to ST6 source terms model. The comparison, performed with respect to the numerical results from [
15] and with earlier work by Gagnaire-Renou et al. [
16] detects the fetch span, where quasi-self-similar behavior is realized. Motivated by these observations, we propose an alternative source terms model, abbreviated to PGZ-2, which improves the ST6 model properties demonstrated in [
15].
The kinetic equation of type (
1) is a well-establish basis for the study of wave spectra. Nevertheless, one should remember that such an approach has its restrictions which arise from the fact that the kinetic equation is an approximation and neglects some correlation properties of the nonlinear waves [
17]. In our article, we deal with the conventional Hasselmann equation as the base of our model; the real ocean may be more complex.
As the PGZ-2 model is still in development, in this paper we present the application of the model to the steady wave solution and perform the comparison of the model to the ST6 one.
2. Materials and Methods
2.1. Numerical Models
Recently, the full version of Equation (
1) has been modeled in the limited fetch statement [
15] with constant offshore wind of 10 m/s orthogonally to the fetch (
Figure 1) in the framework of the WAVEWATCH III operational model. This was an extremely difficult computational problem, first and foremost, due to the necessity of the calculation of the exact nonlinear wave interaction in every real-space 2D grid point on every time step. See [
15] for the details.
In the presented research, the simplified and much more computationally easy approach has been used. We considered the stationary version of Equation (
1) in limited fetch formulation:
used for numerical simulations by [
7,
8,
18]. These works utilized the full four-wave interaction term
, calculated with the help of WRT algorithm by [
19,
20]. The first order forward spatial derivative numerical scheme in the fetch coordinate has been used for the integration of the Equation (
2) but led to the appearance of quasi-monochromatic waves in the direction orthogonal to the wind observed in the numerical experiments of [
21,
22,
23,
24] as well as the field observations by [
25,
26,
27,
28,
29,
30,
31]. The appearance of quasi-monochromatic waves orthogonal to the wind and leading eventually to the numerical instability was circumvented in the publications of [
7,
8,
18] through zeroing of one-half of the Fourier space corresponding to the waves running in the wind direction. The absence of such truncation of phase space for waves, running against and orthogonal to the wind, leads to the eventual blow-up in the numerical simulation of Equation (
2).
The advantage of the described approach consists in the fact that the described model with truncated nonlinear wave interactions uses significantly finer spatial resolution step of just 1 meter for the the whole fetch span. At the same time, the disadvantage of this approach is neglect of wind-opposite and wind-orthogonal wave generation. Nevertheless, this type of approach has, apparently, the physical sense. The first reason to believe in such approximation is the initial development of the wave propagation in the form of combined self-similar solutions for fetch-limited and duration-limited situations [
22,
23]. The second reason consists in the fact that in the asymptotic state, observed in Pushkarev [
21], Pushkarev and Zakharov [
22], Pushkarev [
23], the self-similar structure is preserved for most of the waves running in the wind direction, and finally, the third reason consists in the tendency of the appearance of waves in the direction “almost orthogonal” to the wind, even in the truncated version of the wave model Equation (
2).
2.2. Theoretical Background: Self-Similar Solutions
For the reader’s convenience, we repeat self-similar aspects of wave ensemble behavior studied in [
8,
32]. The equation
allows self-similar substitution
subject to the relations
Equation (
5) is known as the “magic relation”. The self-similar behavior of the wave ensemble can be identified through the simple test: if the integral characteristics of the wave ensemble have power-like dependencies on the fetch, i.e., total wave energy
and mean frequency
, and the indices
p and
q satisfy the relation (
5), then the system behavior is self-similar.
The spectrum shape is mostly defined by the nonlinear interactions which are governed by the Hasselmann equation. Such a situation allows for self-similar solutions which scale both in frequency scale and in energy value. Those scalings are connected to each other.
The system of Equations (5) and (6) do not define values of
p,
q, and
s, since it needs the third relation, which can be taken from the experimental field observations. This procedure has been performed in [
8], where it was taken in the form of the regression line for the specific spectral characteristic of the experimental field measurements [
8], allowing the finding of the set of values
which define the so-called ZRP wind input term.
Below, we will perform a different task: knowing target field experimental relation index p, we will demonstrate that it is possible to come up with the specific value of s, allowing us to build the new model, reproducing and improving partial quasi-self-similar properties of the ST6 model.
Still, this task is not straightforward. The reason is that there is an initial stage of wave turbulence development which is not self-similar. Therefore, the value of s obtained with the help of Equation (6) would be an estimate. Nevertheless, it can be used as the first approximation. Further tuning is required to obtain the newly constructed source term model.
2.3. Theoretical Background: Nonlinear Interactions
The nonlinear interactions are governed by Hasselmann equation [
33], which is the kinetic equation for wave modes of different frequencies and directions exchanging energy by nonlinear interactions.
Here, is the wave number spectrum (the wave number is the wave motion invariant connected to energy spectrum as ), is the angular frequency (for deep water waves ), g is the gravity acceleration, is the wave vector. The nonlinear interaction rate is defined by the coupling coefficient T. The equation includes two -functions which describe the resonances required for the waves to interact. The resonances are both in space (wave vectors ) and in time (frequencies ). Note that the 3-wave resonances are impossible for gravity surface waves because of their dispersion condition. This is why the interaction is through 4-wave interaction. The relative smallness of this interaction provides for the ability to deal with wave modes independently (so called “weak turbulence”).
To calculate the spectra evolution, the new algorithm was used. This algorithm is based on the pre-configured grid of the interacting 4-wave quadruplets. The grid was chosen based on different quadruplet significance with the emphasis on the quadruplets near the center of the so-called Phillips curve. Those quadruplets carry the majority of the energy transfer.
The quadruplet interaction uses quadruplets satisfying the resonance condition:
The exact obedience of the resonance condition for each quadruplet ensures energy conservation, preventing error accumulation. As the quadruplet vectors do not coincide exactly with the grid points of the energy spectrum, we have to use the bilinear interpolation to determine the spectrum value at each quadruplet vector. Still, energy conservation is ensured by the calculation algorithm which specifically maintains that nonlinear interactions result only in redistribution of the energy.
The use of the quadruplet grid allows for easy parallelization of the calculations. The interaction for each quadruplet may be calculated independently, and therefore, each processor core may be used for the calculation of some part of the quadruplet grid without any loss of the processor time.
The backward scattering presents another problem in calculations. We use a small initial wave spectrum near the coast as the boundary condition. The wind pumping produces mostly the forward waves (which propagate from the coast to the open ocean). The nonlinear interaction mostly keeps this outward direction. Therefore, to simplify the calculation, we ignored the waves which are directed towards the coast. When doing the calculations, we manually zeroed that part of the spectra at each numerical step.
The significant part of the model is the high-frequency tail addition for the energy spectrum. Both the ST6 and PGZ-2 models include, together with the pumping, some kind of energy dissipation, but the main part of the energy pumped stays in the energy spectrum. This energy is being transferred from the lower frequencies to the high frequencies by the Kolmogorov cascade transfer, which results in the formation of Zakharov–Filonenko spectrum [
34]. This spectrum has
with the energy transfer proportional to the energy scale in the 3rd power. The numerical calculations may cover only some part of the frequency scale, and the the energy transfer outside this part is unavoidable. Thus, one needs to provide for some kind of continuation.
O. Phillips in [
35] suggested another form of the spectral tail, namely
Here,
is the universal Phillips constant. This asymptotic is observed for the utmost high frequency area, while in the energy capacity spectral band, the ZF spectrum is realized.
We use the classic Phillips spectrum as the tail which we added artificially outside the numerical calculation zone.
2.4. ST6 Model Numerical Simulation
We reproduced the numerical experiments, performed in the frame of full model Equation (
1) in [
15], using truncated model Equation (
2). Numerical experiments of [
15] targeted the field experimental benchmark for the total dimensionless wave energy dependence on the fetch [
36], hereafter KC1992:
where
is significant wave height,
E is the total energy,
is dimensionless fetch,
is friction velocity, and
. Omitting the details of the ST6 model, we present the following tuning parameters of the ST6 source terms, described in detail in Liu et al. [
15] for WRT nonlinear interaction case:
The presented numerical results were obtained using this parameters set.
We present our set of parameters value below at the
Section 3 (see Equation (
16)).
Figure 2 shows the total wave energy behavior dependence on the fetch for various models and field experimental observations as well as different values of the total energy on the shoreline.
One can see that, for the ST6 model described by Equation (
2), the curve corresponding to the shore total energy value
m
2 is in a good agreement with the target dependence KC1992 for the 5–40 km fetch span as well as the numerical results of [
15]. The situation is different, however, for other curves corresponding to lower shore total wave energy levels
m
2 and
m
2, which exhibit noticeable deviation from the KC1992 target dependence. Additionally, in the simplified model Equation (
2), we were not able to reach the “steady sea state” observed in the numerical experiments by [
15] for the 200–600 km fetch range. One more difference (idealization) with WW3 and other wave models is disactivating the linear input term, which is important at the initial stage of wave development. The effect of this is seen as a difference between PGZ-2 and other models for short dimensionless fetches. We plan to perform more detailed tuning of parts of our model to minimize the effect of numerical problems [
37].
Figure 3 shows the local index
p of power-like total wave energy dependence, derived from the data, plotted on
Figure 2. It shows that the value of
p in the ST6 case, described by Equation (
2), asymptotically converges to
with fetch growth independently on the total wave energy value at the shoreline.
Similarly,
Figure 4 plots the local index
q of power-like mean frequency dependence in the ST6 case, described by Equation (
2), converging asymptotically to
.
Figure 5 shows the combination
for the same case, asymptotically converging to
. The observed value of this combination is not equal to 1 as it is supposed to be in Equation (
5), and we are observing, therefore, the above-mentioned partial (i.e., related to portion of the fetch) quasi-self-similarity case for the ST6 model described by Equation (
2).
Figure 6 shows 3D wave energy spectrum for the fetch
km as the function of frequency
f and angle
, localized in the angular spread of
around the spectral maximum frequency for the ST6 case described by Equation (
2).
Figure 7 shows the angle-averaged wave energy spectrum corresponding to the previous figure as the function of the logarithm of frequency.
One can see the area of the spectral maximum, Kolmogorov–Zakharov (KZ) spectrum
[
34], and Phillips dissipation area
. These observations confirm the already known properties of surface wave turbulence observed in many field observations as well as numerical simulations before.
The general spectra structure, as obtained in numeric experiments, is the following. The pumping at some intermediate frequency results in the initial growth of energy spectrum there. This growth is bounded by the nonlinear interactions. Those nonlinear interactions result in the formation of the Kolmogorov spectrum which provides for the energy flow into the high frequencies by the cascade transfer of energy. The high frequency part of the spectra stabilizes in this way while the lower frequency part of the spectra continues its evolution. The energy pumping into the higher frequency part of the spectra results in the formation of the inverse cascade (which is the cascade transmittion of the wave number into the lower frequency part of spectra). Some energy is also transferred there, but it is much slower than the energy transfer into the high frequency, and this transfer is negligible in the total balance of energy transfer. This general pattern is present in all numerically calculated spectra and is also found in experimental ones. The cascading energy transfer by Zakharov–Filonenko spectra is very stable with respect to any disturbances.
4. Discussion
This is the first article about our new model. The PGZ-2 model described here is a new development. The current version should be considered as a preliminary one. We plan to perform the various checks of it and compare it to the existing models and experimental data.
The development and the augmentation of the model may be performed in some different ways. The most obvious one is the further tuning of the coefficients. Other tuning may be a detailed studying of the transition between the and high-frequency spectrum tail. The transition point is defined by the balance between the nonlinear transfer of energy and damping and theoretically may be calculated from this balance.
The fetch-based calculations may be improved by including some kind of the back-scattering, i.e., the waves traveling (and transferring energy) to the coast from the open sea. Additionally, some energy may be transferred by the waves traveling along the coast, which also require consideration.
We believe that the numerical experiments using our new algorithm may be performed to improve the PGZ-2 model and to better tune it to the actual ocean wave conditions. We plan to continue the work on the PGZ-2 model and use it both for the fetch-based and duration-based calculations.
We performed the assessment of the ST6 model of wind energy input and wave dissipation source functions through numerical simulation of a stationary version of the Hasselmann equation with an exact four-wave interaction nonlinear term in the truncated Fourier space.
It is shown that the ST6 model exhibits partial quasi-self-similar properties, reproducing the reference field experimental total wave energy behavior for the fetch span from to km. For smaller values of the fetch, the ST6 model exhibits deviations from the target total wave energy reference curve KC1992, depending on the boundary value of the wave energy on the shore line.
We were neither able to reach the “steady sea state” state with the help of the truncated model nor to see any physical reason for its existence.
We formulated the new model of the balanced wind energy input and dissipation source terms with two unknown tunable parameters, the amplitude and index of the power-like wind source term, and call it PGZ-2. Those unknown values are obtained with the help of matching the PGZ-2 numerical results with the field experimentally observed energy growth curve KC1992 and theoretically predicted PGZ-2 self-similar relations.
The formulated PGZ-2 model significantly improves the correspondence of the total wave energy behavior for smaller fetches as well as larger ones with experimental curve KC1992.
Both ST6 and PGZ-2 models exhibit dependence of the wave energy growth along the fetch on the total wave energy value at the shore. The scatter of this dependence is lower for the PGZ-2 model.
One should note that the alternative PGZ-2 model should be construed as the demonstrational one and not for practical purposes since the KC1992 experiment analysis presents not quite typical behavior for the pool of the multiple experimental observations, studied in Badulin et al. [
12], where the characteristic power-like total wave energy index value was
. In this relation, in the authors’ opinion, though close in spirit, the different ZRP approach to the formulation of the balanced source terms presents the better alternative.
The tuning of the PGZ-2 model was performed to be in agreement with the experimental data. More fine tuning may be performed on the basis of more numeric experiments.
We believe that the PGZ-2 model is a basis for further study. In the near future, we plan to study its properties in more detail and continue the tuning of the model.