Asymptotes of the nonlinear transfer and wave spectrum in the frame of the kinetic equation solution

The kinetic equation for a gravity wave spectrum is solved numerically to study the high frequencies asymptotes 10 for the one-dimensional nonlinear energy transfer and the variability of spectrum parameters that accompany the long-term evolution of nonlinear waves. The cases of initial two-dimensional spectra S(ω,θ) of modified JONSWAP type with the frequency decay-law S(ω) ~ ω -n (for n = 6, 5, 4 and 3.5) and various initial functions of the angular distribution are considered. It is shown that at the first step of the kinetic equation solution, the nonlinear energy transfer asymptote has the power-like decay-law, Nl(ω) ~ ω -p , with values p ≤ n-1, valid in cases when n ≥ 5, and the difference, n-p, changes 15 significantly when n approaches 4. On time scales of evolution greater than several thousands of initial wave periods, in every case, a self-similar spectrum Ssf(ω,θ) is established with the frequency decay-law of form S(ω) ~ ω -4 . Herein, the asymptote of nonlinear energy transfer becomes negative in value and decreases according to the same law (i.e., Nl(ω)~-ω -4 ). The peak frequency of the spectrum, ωp(t), migrates in time t to the low-frequency region such that the angular and frequency characteristics of the two-dimensional spectrum Ssf(ω,θ) remain constant. However, these characteristics depend 20 on the degree of angular anisotropy of the initial spectrum. The solutions obtained are interpreted, and their connection with the analytical solutions of the kinetic equation by Zakharov and co-authors for gravity waves in water is discussed.


Introduction
In the study of random field features for the nonlinear gravity waves in water, the four-wave kinetic equation (KE) plays a very important role.In (1), ()

M ( , , , )
k k k k is the second power of the matrix elements corresponding to the four-wave nonlinear interactions, 0 1 2 3 delta-function responsible for the resonant feature of the four-wave interactions, and ω i = ω(k i ) is the radian frequency of the wave component with wave vector k i .The relationship between wave vector k = (k x , k y ) = (k, θ) and wave frequency ω is given by the dispersion relation.For gravity waves in deep water, considered here, it is;  , and related linearly to N( ) k (e.g., Badulin et al., 2005) 42 2 S( , ) N( , ) ( / g )N( ) The expressions for matrix elements 0 1 2 3 , , , M are well known (e.g., Hasselmann, 1962;Webb, 1968;Zakharov, 1968; see also Badulin et al., 2005), but we use the elements from Crawford et al. (1980).In any representation, formulas for elements 0 1 2 3 , , , M are very cumbersome, so they are not given here.It is important to note only that the KI formally preserves the total wave energy, the total wave action,
First, for the case of an angular isotropic spectrum spread throughout the infinite frequency band (0 < ω <∞), Zakharov and Filonenko (1966)  .They have proposed to interpret these solutions as the Kolmogorov-type spectra of the constant energy-flux, P E , directed upward in frequencies, and constant wave-action flux, P N , downward in frequencies, respectively.
In this interpretation, the mentioned spectra acquire the proper Kolmogorov-type representation of the form (Zakharov and Zaslavskii, 1982; see also Pushkarev et al., 2003;Badulin et al., 2005) and 1/3 4/3 11/3 2 ( ) g where c 1 and c 2 are the dimensionless Kolmogorov's constants.The proper sources and sinks of these fluxes, obligatory in the Kolmogorov's theory (Monin. and Yaglom, 1971), were assumed to be located at the zero and infinity frequency-points, in dependence on the flux respectively (Zakharov and Zaslavskii, 1982).It should be noted, however, that this interpretation was a hypothesis, at that time, because it was based only on dimensional considerations and similarity to the Kolmogorov's theory, but does not follow from the Euler equations without external forces, having no sources and sinks.Later, in a large series of papers by Zakharov and co-authors (see below), the mentioned analytical results were sophisticated in different directions: generalisation for a case of anisotropic spectra in (Zakharov and Pushkarev, 1999), basing on paper by Katz and Kontorovich (1974); description of the spatial and temporal asymptotes of the self-similar spectra of form (6) (Pushkarev et al., 2003;Badulin et al., 2005); and searching for convergence conditions for the KI (Geodjaev and Zakharov 2017;Zakharov 2017).These analytical results are represented in details in the cited papers and in numerous references therein.
In parallel, a study of KE was carried out numerically, beginning with the works by Webb (1978), Masuda (1980), andHasselmann andHasselmann (1981), devoted to the development of methods for calculation of the KI.Finally, only two methods have survived: by Webb (1978) and Masuda (1980).The first of them was later elaborated by Resio and Perrie (1991) and Van Vledder (2006).It is used by these authors and the Zakharov's group.The second one elaborated by Polnikov (1989), Lavrenov (2003), Komatsu and Masuda (1996), Gagnaire-Renou, et al, (2010), who use it.Evaluation of these methods was done in (Benoit, 2005).Here it should be noted that the most important part of the KI-calculation algorithm is the method of estimating the contribution of integrable singularities in the integrand (the contribution of 'singular points'), resulting from integration of the frequency delta-function (Masuda 1980;Polnikov 1989).The different approaches to this point represent the main differences in the algorithms mentioned (e.g., Masuda, 1980;Van Vledder, 2006;Polnikov, 1989;Lavrenov, 2003).Herewith, Polnikov (1989) have showed numerically that the contribution of the singular points determines crucially the conservation balance for the integral values E, N, M (Eqs, 3-5), and the accuracy of the KIestimate as a whole.
For the first time, a detailed numerical study of the KI properties was carried out by Masuda (1980).It was expanded by Polnikov (1989) via introducing classification for the two-dimensional spectral shapes, S( , )  .After that, Polnikov (1990)   has performed a numerical solution of KE (1), where it was shown the fact of establishing the self-similar shape of wave  (ω,θ), in the course of long-term evolution of nonlinear waves.This result was later confirmed in a series of papers by the Zakharov's group (see references in Pushkarev et al., 2003;Badulin et al., 2005) and others (e.g., Lavrenov, 2004;Gagnaire-Renou and Benoit, 2007).
In addition to Polnikov (1990), Pushkarev et al. (2003) and Badulin et al. (2005) have established that the tail of the selfsimilar spectrum S sf (ω) falls according to the law S sf (ω,θ)~ ω -4 , and the shape of spectral peak for S sf (ω,θ) is close to one for the JONSWAP spectrum with γ = 3.3 (Badulin et al., 2005).They first noted that the only condition of conservation for total wave action N takes place, during wave spectrum evolution according to KE (1), and this evolution is accompanied by a leak of the total wave energy E and appearance of a certain energy flux P E forwarded to the higher frequencies.This result was treated as the formation of the Kolmogorov-like spectrum of form (6) in the frame of solving KE (1).
The mentioned effects represent the most interesting points for further investigations, taking in mind that all three integral values: E, N, and M, should be formally preserved simultaneously, as far as KE ( 1) is derived from the conservative Euler equations without external forces (Hasselmann, 1963;Zakharov et al, 1992).The present work is devoted to clarification this point in more details.In addition to this, other essential features of the spectral peak shape and the process of forming twodimensional self-similar spectrum S sf (ω,θ) are still not sufficiently described, i.e., the high-frequency asymptotes of the nonlinear energy transfer at short-and long-time evolution, integral parameters of the spectral peak and angular distribution for S sf (ω,θ).The mentioned points will be elaborated in the present paper also.
Besides of studying the shape of self-similar spectrum S sf (ω,θ), a lot of papers are devoted to a searching for establishing the Kolmogorov-spectra of forms ( 6) and ( 7).First, under advice of Zakharov, Polnikov (1994) has showed numerically that these spectra are actually formed as the result of numerical solution for the extended KE having the form In (8),

Dis  
are the source and sink functions, respectively, which are separated in the frequency space, to ensure the presence of the inertial interval, as the main condition for applicability of the Kolmogorov's theory (Monin and Yaglom, 1971).Despite of technical limitations in performing calculations, these results (Polnikov, 1994) gave grounds for the validity of the hypotheses by Zakharov and Zaslasvskii (1982).Similar, and more detailed results were later obtained in numerous studies (e.g., Pushkarev et al., 2003;Lavrenov, 2004;Badulin et al., 2005;Badulin and Zakharov, 2017;and references therein) in which other algorithms and thicker numerical grids were used to compute the KI.They confirmed the fact of establishing spectra of forms ( 6) and ( 7), depending on the configuration of the source and sink locations on the calculating frequency band.In addition, Pushkarev et al. (2003) and Badulin et al. (2005) have analysed the long-term asymptote of the peak-frequency downshift, ω p (t), estimated the Kolmogorov's constant values, and determined the rate of wave-energy leakage E(t) in the case of preservation of total wave action value N.These results were applied in the recent paper by Badulin and Zakharov (2017), devoted to a modelling numerous details of the spectral shape for self-similar spectrum of swell, including the well-known effect of the swell angular bimodality at high frequencies.S sf (ω,θ), as well as their dependence on the initial conditions, are not sufficiently specified (despite of the recent paper by Badulin and Zakharov, 2017).There is also a need to analyze the impact of the conservation condition for total wave energy E (besides of conservation of wave action N) on the shape of the self-similar spectrum.This task could be especially realized as a theoretical scenario for the numerical solution of KE (1) without sources and sinks.The comparison of these two cases of the numerical solution scenario for KE is very interesting for understanding the nature of establishing self-similar spectrum S sf (ω,θ) mentioned above.
This paper is devoted to studying these details of numerical solution for KE (1).The methods of research and computations are presented in Sect. 2. Calculation results and analysis are given in Sect.3. Section 4 is devoted to the interpretation of the results, and Sect. 5 contains conclusions.Appendix (A) contains the analytical derivations of asymptotes for the peak frequency downshifting, ω p (t), accompanying the self-similar spectrum evolution.
2 Methods of research

Scenarios of the KE solution
First of all, let us dwell shortly on the point of the laws of conservation in the frame of KE (1).From the very beginning, it was proved that the Euler equations without external forces are the conservative system (Zakharov, 1968), i.e. the energy of the system does not change in time.Thus, the KE, derived from the mentioned equations, should also preserve wave energy, at least.Eventually, it was proved that KE (1) formally preserve three integral wave values: energy E, action N, and momentum M, given by Eqs.(3-5), as was already mentioned.Though, a proper special testing of the known and actively used algorithms for the KI-calculations is not presented in literature in all details, being, usually, the intrinsic matter of the authors (e.g., Masuda, 1980;Polnikov, 1989;Resio and Perrie,1991;Van Vledder, 2006).
Herewith, it is important to note that all the mentioned algorithms provide establishing the self-similar spectrum, S sf (ω,θ) (see numerous references above), what could contradicts to the preserving three integral values: values: E, N, and M, as was many times noted (e.g., Pushkarev and Zakharov, 2000;Pushkarev et al, 2003;Badulin et al, 2005;Geodjaev and Zakharov, 2017;Zakharov, 2017).Indeed, let us introduce representation for the one-dimensional self-similar spectrum, Nl  , should also be self-similar, as it is the integral from self-similar functions.Thus, one may write where c=π/16 is the theoretical coefficient, R Nl = 4 3 11 pp gS   is the dimension of non-linear energy transfer Nl(ω).(Polnikov, 1989), and ( / ) sfN p F  is the n.d.self-similar spectral shape of Nl sf (ω).Then, conditions (3-4) lead to the system of the following equations (here consider only two values E and N, for simplicity) where / p     is the n.d.frequency.System (11) gives restrictions for the shape of () sfN F  (and, implicitly, to   ()   sf F  ).The question under consideration is: does a solution of system (11) exist?As function () sfN F  is prescribed (it is defined by the KI) and has one positive lobe in the range σ ≤ 1 and one negative lobe for σ >1 (as well known, e.g., Pushkarev et al., 2003;Badulin et al., 2005; and evidently shown below), it follows that system (11) is inconsistent, and simultaneous preservation of two values, E and N , is impossible.That is the paradox of the KE, which needs its own explanation.At some extent it is done in (Zakharov, 2017) via the bad convergence of the KI for slow-falling spectra, though it needs more details.
Many years ego, on the basis of specific numerical simulations with KE, it was pointed out that the wave energy should leak to the higher frequencies due to the four-wave nonlinear interactions (Pushkarev and Zakharov, 2000).Since that time this idea is applied in any case of considering KE (1) by the Zkharov's group (e.g.Pushkarev et al., 2003;Badulin et al., 2005;Geodjaev and Zakharov, 2017;Zakharov, 2017), though, in our mind, there is no convincing proof of it.Usually, these authors say about conservation of N and leakage of E and M, do not dwelling on physical reasons of zero-balance mismatch for them in the frame of KE solution.As we have no convincing (say, categorical) prohibition for saving E during the solving of KE (1), in our mind, two scenarios of preservation for both of these quantities may be chosen, at least, as an interesting theoretical alternative.It occurs that this way gives another paradox result (see below).

Numerical details
In the numerical study of the KI-properties and features of the KE-solution, the choice of a certain frequency-angular numerical grid (ω,θ) provides (at a certain extent) the accuracy of calculations (Polnikov, 1999;van Vledder, 2006).In our study, to examine the asymptote of the non-linear transfer of energy (NLT), the extended frequency grid and enhanced angular resolution are chosen.They are given by the ratios: and the previous values for ω 1 and q.This was done for the reasons of arising technical difficulties during numerous and long-term calculations on grid (12a,b) (high time consuming and appearing numerical instability).
The initial wave spectrum was given in the slightly modified form for the well-known JONSWAP spectrum (J-spectrum) where  = 0.07 to 0.09 is the n.d.peak-width parameter of the J-spectrum, γ = 1 to 7 is the n.d.'peakness' parameter, ( , )   is the n.d.frequency-angular form, and is the Pearson-Moskowitz spectrum, modified to the case of an arbitrary degree of the spectrum-tail decay, n (following to Polnikov 1989).The initial peak frequency, ω p (0), was assumed to be 2 rad/s, whilst the initial form of function ) , (    was assumed to be independent of the frequency and given in the form Variation of the parameters a and m ensure the assignment of any degree of angular anisotropy of the spectrum used for calculation of the KI at the initial stage of spectrum evolution. The algorithm for calculating the KI was used according to work by Polnikov (1989), where it was shown a high smoothness and numerical stability of the proper NLT-estimates.Comparison of the KI-estimates on grids (12a,b) and ( 13) has showed that the one-step KI-calculations have relative deviations within 3% to 5% (for peak values of the NLT), what gives the upper limits for the errors of calculations.
To check the conservation balance Q, we have used the ratio where the sub-indexes E and N mean the cases of considering NLT for the energy spectrum or wave action one, respectively.
Our experience shows that for the fast falling spectra (faster ω -5 ), balance Q has the value about of 1-2% for both values: E and N (M was out of our considerations), For the slower falling spectra, the balance for N retains the same, whilst for E it becomes negative and has the order of 5-10% (depending on the shape of spectrum and stage of evolution, see below).The more detailed consideration of this point is out of this paper.It needs a special research, as mentioned already.Nevertheless, the presented estimates are quite acceptable for our purposes to calculate the wave evolution separately: with the exact preservation of E, and the same for N.In the KE-solution, the exact balance, Q E = 0 or Q N = 0, is made by correcting the whole negative part of the relevant NLT , , ( , ) EN Nl  , with a proper coefficient, at each time-step.It is done for the "purity" of the scenarios mentioned.This correction does not influence on the high frequency asymptotes of NLT, absolutely.Herewith, it saves the physics of spectrum evolution, because the mentioned corrections are weaker than one-step changes of the spectra in the course of the KE-solution.
The behaviour of the high-frequency asymptote for the one-dimensional NLT, NL(ω), was determined by the leastsquares estimation of the decay parameter p in the power-law dependence for NL(ω) of the form The estimation of p was executed in frequency band 5 ω p < ω < (15-20)ω p for numerical grid (12).The quantitative characteristics of the two-dimensional spectrum-shape were determined in the integral representation form by introducing the spectral frequency-width parameter, B, and the angular-narrowness function of spectrum, A(ω) ( Polnikov, 1989( Polnikov, , 1990)): Here, is the value of the one-dimensional spectrum at the peak frequency, and p  is the direction of propagation of the peak-component for a two-dimensional wave spectrum.Hereafter, p  = 0.
The numerical solution of KE (1) was performed according to the algorithm by Polnikov (1990), which includes an explicit numerical scheme of the first order of accuracy, linear interpolation of the spectrum from the nodes of the computational grid to the resonance points of four-wave quadruplets, and the choice of the time-step, Δt, according to the condition Ratio (20) indicates restriction of the spectrum-intensity change to 3% or less of the current value, S( , ,t )  , at each time- step.Although this choice of Δt led to a significant increase of computing time, it was adopted, because it ensures smoothness and numerical stability of the solution, resulting in good final results.
In this paper, all calculations were performed in the dimensional units.First, the scenario of preservation for E was realized, in addition to considerations by the Zakhariv's group (see references above), assuming the wave system is the conservative one.Then, some KE-solutions were done with the scenario of preservation for N, to demonstrate the difference between them.Below, these scenarios are named as the algorithms versions for the KE-solution.Analysis of the spectral forms of numerical solutions for KE (1) was carried out at evolution times exceeding 10 5 /ω p (0). 3 Calculation results and analysis

Asymptotes of the NLT
The one-step calculations of KI on grid (12) were performed for several spectral forms, the representative set of which is given in Table 1, where the estimates of parameter p for the power asymptote of the NLT of form ( 18) are also shown.The typical forms of the calculated one-dimensional functions NL(ω) and their asymptotes are shown in Figs 1a,b.  1. Line 1 corresponds to the whole tail of NLT, line 2 to the part of tail (10 <ω <30) rad/s (with a weight of 3, to separate the lines).Line 1' is the root-mean-square trend of the tail-part 1 (equation, y = 0.0004x -4.4 ), and line 2' is the same for the tail-part 2 (equation, y = 0.002x -4.9 ).The main features of NL(ω)an alternating two-lobe shape near the spectral peak and an infinite positive high-frequency tail (Fig. 1a)are well known (Hasselmann and Hasselmann 1985;Polnikov 1989;van Vledder 2006).Here, it was found that the high-frequency tail of the NLT for all kinds of J-spectrum (14-16) is positive only when n ≥ 4; when n < 4, it is negative and poorly defined due to the weak convergence of the KI for such spectra (Zakharov, 2017).In general, the features of the high-frequency asymptotes for NL(ω) at the first-step of solving the KE (1) can be formulated as follows: 1) The high-frequency tail of the NLT (for J-spectra) is positive when n ≥ 4 and negative when n <4; 2) Representation of NL(ω) in form ( 18) has a varying power-law of decay; 3) For values n ≥ 5 , parameter p for NL(ω) is close to value n-1 (i.e., Nl (ω) ~ ωS (ω)), in the intermediate frequency region, 5 < ω/ω p < 15; but in the entire frequency band, 5 < ω/ω p < 25, the decay of Nl(ω) is somewhat weaker than ωS (ω); 4) For fixed n > 5, parameter p increases with a decrease of peakness parameter γ and increase of angular anisotropy parameter a.
5) When spectrum-decay parameter n approaches to 4, the relative intensity of the NLT-tail decreases radically (Fig. 2a); and when n ≤ 4, the decay-features of NLT begin to depend significantly on the spectrum parameters: γ, a, relative frequency ω/ω p (Fig. 2b), and on the limits of the computing grid in units of ω/ω p .This feature is due to the slow convergence of the KI in such a case (see other numerical details in Polnikov and Uma, 2014).
(a) (b) Figure 2: (a) Form of NL(ω) at low frequencies for the initial spectrum of run 5 from table 1.(b) Shape of NL(ω) at high frequencies for the same spectrum.Line 1 is NLT, line 1' is the root-mean-square trend for the part of tail in the band: (10 <ω <30) rad/s (equation, y = 0.0024x -0.9 ).
From the behaviour of the asymptote for NLT, it follows that the nonlinear interactions at high frequencies are sensitive to both the frequency-shape of the spectral peak (i.e., variations of parameter γ) and the angular shape of the spectrum (i.e.,

Kinetic equation solutions on large time scales
To achieve the goals posed, we performed a large series of numerical solutions for KE of type (1) on grid (13) with the algorithm versions both with conservation of total wave action N and with retaining total wave energy E. The difference between these versions was realised by making the exact balance of positive and negative NLT for either N or E, respectively, at each time step of the KE numerical solution.For the first time, our calculations have shown that on large scales of evolution time, this difference impacts slightly on the features of the self-similar shapes for NL(ω) and S(ω), in the peak-domain (0.5 ω p < ω < 2 ω p ) , and does not impact on the asymptotes for their high-frequency tails (ω > 2 ω p ).To show these results, we dwell, mainly, on the features for one-dimensional functions.
For each version of the KE-solution, the following results were established: 1) The self-similar one-dimensional spectra, S sf (ω), in different cases, have slightly different frequency widths B and decay with the law 2) The self-similar shapes of Nl(ω) are slightly different, in different cases, and have the high-frequency asymptotes as Nl(ω)~ -ω -4.15 ± 0.05 ; 3) The shape parameters for the two-dimensional spectrum S sf (ω, θ): A p = A(ω p ) and B, vary slightly in time within ±5% (as spectra are evaluating), and their average values depend on the degree of anisotropy for the initial spectrum S(ω,θ). 1 The summarised results of calculations for a representative series of initial spectrum shapes S(ω,θ) are presented in digital form, for case E=const, in Table 2, together with the asymptote parameters for Nl (ω) and parameters A p and B, for self-similar shape S sf (ω, θ).In case N=const, values of B are on 1-2% smaller, whilst values of А p are the same (not presented due to similarity).Table 2 shows that the values of А p and B differ significantly in the cases of isotropic and anisotropic initial spectrum S(ω,θ), though they have intermediate values (run 7) for weakly anisotropic initial spectra.
1 Hereafter, all estimates of the powers are obtained in the EXCEL-shell with the installed trend-formulas (obtained by the method of least squares).Scattering the exponent powers, in ( 21), ( 22) and further, are due to the scattering of the obtained set of estimates, exceeding the statistical errors for the exponents.Note.Different degrees of shadowing show differences in angular forms for self-similar spectra.Value В=0.33 corresponds to γ=3.3, В=0.63 to γ=1.0 and А p = 0.16 to Ψ= const, А p = 0.64 to Ψ= сos 2 (θ).
For anisotropic spectra, the function of angular directivity, A(ω), includes the sharp peak slightly below peak frequency ω p (Fig. 3a).Shape of A(ω) does not depend of the algorithm version.With the frequency increasing, A(ω) quickly becomes constant, A ≈ 0.16 (when ω ≥ 2 ω p ), corresponding to the isotropic distribution of the spectrum (Fig. 3a).In turn, the frequency shape of one-dimensional spectrum S sf (ω) has the sharp peak too.For isotropic initial spectra, its shape is very similar to one for the J-spectrum with parameters γ = 3.3 and n = 5 (Fig. 3b).The slight difference of shapes for S sf (ω), in cases E=const and N=const, is shown in Fig. 3c, in terms of normalised spectra S nor (ω/ ω p ).For the anisotropic initial spectra, the shape of peak for S sf (ω) is similar to that above, although the peak is 1.5 times wider (in the value of B).Such a detailed description of the integral parameters for two-dimensional self-similar spectrum S sf (ω, θ) essentially supplements the results of earlier works (Polnikov 1990;Pushkarev et al. 2003;Badulin et al. 2005).It has a reasonable academic interest.
The shapes of self-similar one-dimensional functions Nl(ω).are shown in Figs 4a,b in the normalized form, Nl nor (ω/ ω p ) for several time moments, It is seen that Nl(ω) has only two lobes, which differ slightly for cases of E=const.and N=const, providing the proper balances (11).In terms of Q (Eq.17), they are about 0.05-0.1%.More details about features of Nl nor (ω/ ω p ) are out of this paper aims.The most important fact about shapes of Nl(ω) for the below analysis is that, in any case, their tail asymptotes have form (22). Line 1 is the mean for all anisotropic runs (7 through 11) from Table 2; line 2 is run 10; line 3 is run 11.(b) Self-similar shapes S sf (ω/ω p ) for run 3 from Table 2 (isotropic spectrum).Line 1 is the shape of J-spectrum according to (14-16) (the same run); line 2 is S sf (ω/ω p ) at time t ≈ 1·10 4 s; line 3 is S sf (ω/ω p ) at the time t =1·10 5 s.(c) The difference of shapes for normalised self-similar spectra S nor (ω/ω p ) in two cases.Run 4 from Table 2.
. Detailed information regarding the features for spectrum evolution is shown in Figs.5-8 (using the versions of KEsolution with constant wave action N are marked in the proper captions).The time-history of a one-dimensional spectrum S(ω) for run 1 from Table 2 is presented in Fig. 5a.The corresponding time-history of one-dimensional NLP Nl(ω) is shown in Fig. 5b.It can be seen that when the self-similar shape of the spectrum is practically formed (at a time of the order of 4•10 4 s, or more), the NLT changes principally in such a way that the intensity of the tail for Nl(ω) becomes fast going down.
Pay attention that for ω>2ω p , Nl(ω) is never equal to 0, and its asymptotic behaviour is defined by formula ( 22).The same situation is realised for all the cases considered (with small variations for parameter p, as shown in Table 2).2. Line 0 corresponds to t = 0 s; line 1 to t = 1.4•10 5 s; line 2 to t = 3.0•10 5 s; line 3 to t = 1.2•10 6 s; line 3' is the tail part of line 3 (with a weight of 0.3); 3" is the trend of tail 3' (equation y = 0.0792x -4.1 ).(b) Time-history of one-dimensional NLT Nl (ω, t) (the same run).Line 0 corresponds to t = 0 s; line 1 to t = 470 s; line 2 to t = 4.5•10 4 s; line 3 to t = 3.7•10 4 s.In this respect, run 6 from Table 2 is very instructive (Figs.6a, b). Figure 6a demonstrates the evolution of the initially isotropic spectrum (a = 0) with parameters γ = 1 and n = 4.In this case, for the initial J-spectrum decaying with the law: S(ω) ~ ω -4 , Figure 6b shows that, at the first time-step, the intensity of the tail for nonlinear transfer Nl(ω) is not small, as it might be expected in accordance with the conclusions by Zakharov and Filonenko (1966).Later, on the relatively great evolution scales, t > 10 4 •1/ω p (0), the self-similar spectrum has already formed, having the same decay-law as the initial spectrum: i.e., S(ω) ~ ω −4 .Though, on these scales, the nonlinear transfer gets the rapidly decaying tail of form ( 22), and the intensity of this NLT-tail is already very small, corresponding to the theory of Zakharov and Filonenko (1966).The same results we have got using initial spectra with n < 4 (e.g., the case with n = 3.5, which is not shown, as it is very similar to run 6, and Figs.6a,b).The said is very interesting and important result shown here for the first time.
Undoubtedly, the main reason for the small values of the tail-intensity for function Nl(ω) at a large t is the formation of the specific, self-similar shape S sf (ω,θ) for the entire spectrum.In this case (run 6), the main difference between initial form S(ω,θ) and the final one, S sf (ω,θ), is determined only by the self-similar frequency-shape of the peak-domain for S sf (ω,θ), because the angular distribution of the spectrum is always isotropic.This fact testifies, first, the nonlocality of the four-wave nonlinear interactions (already mentioned above) and, second, the fundamental role of the spectrum-peak shape in forming the tail of NLT-function Nl(ω) with the frequency-asymptote of form ( 22).(b) The version with condition N = const.Line 0 corresponds to t = 0 s; line 1 to t = 2500 s; line 2 to t = 6.7•10 4 s, line 3 to t = 1.7•10 5 s, line 4 to t = 6.1•10 5 s, line 4' to the tail part of line 4 (with a weight of 0.5), line 4'' to the trend for line 4' ( equation y= 0,09•x -4.06 ) At the same time, some numerical values are different for the various versions of the algorithms, as seen from the temporal asymptote for certain spectral parameters.It is evidently true for total wave action N(t) and total wave energy E(t).
The same can be said about downshifting the peak frequency, ω p (t) (see below).But more interesting is the behaviour of the fluxes for wave energy P E and wave action P N , taking place in two different versions for the KE-solution algorithms.
To see this, the frequency functions of fluxes for wave energy P E and wave action P N were estimated with the formulas similar to ones from (Pushkarev et al. 2003;Badulin et al. 2005): where ω 1 is the lower limit of the numerical frequency band.According to (Pushkarev et al. 2003;Badulin et al. 2005), the positive value of P E (or P N ) indicates the flux upward in frequency and the negative value does downward.
For run 2 from Table 2, the time-history of fluxes P Е (ω,t) and P N (ω,t), under the condition of constant wave action N, is  2 (when N is constant).Line 0 corresponds to t = 0 s; line 1 to t = 700 s; line 2 to t = 2500 s.(b) Time-history of flux P N (ω,t) for the same run 2. Lines correspond to the same time moments.action N is constant, a quasi-constant tail for upward energy flux P Е (ω) is established in frequencies range ω > 2ω p (t), on long evolution scales: t > 10 4 /ω p (0).The latter means the loss of wave energy in the system.The flux at the upper edge of frequency band, P N (ω max ,t), is always zero, what means preservation of N.Under the condition of constant wave energy E, the same is true for downward flux P N (ω) (Figs.9a, b).Pay attention that, in proper cases, the tails for both fluxes are not constant in ω but to quasi-constant, because of Nl(ω) is not equal to 0 at higher frequencies.The fluxes asymptotes are defined by formulas (22,23).Due to very week dependences of these tails for P Е (ω) and P N (ω) on ω, we could call them as the constant ones, following to (Pushkarev et al. 2003;Badulin et al. 2005) for simplicity of the further analysis.Herewith, the shape of the self-similar spectra obtained in both versions of the KE-solution is nearly the same, though, the spectra differ in their intensities (Figs.7a,b) and in the rate of formation, defined by different versions of the KE-solution algorithms (see below).This difference is manifested via the different temporal asymptotes for peak frequency ω p (t), values of spectral peak, S p (t) = S sf (ω p ), and, evidently, temporal asymptotes for integral characteristics: E(t) or N(t) (in proper case).
All of these asymptotes could be derived analytically, among them the basic one is the temporal downshifting of peak frequency, ω p (t).

Correspondence between numerical and analytical results
Basing on the self-similarity for S sf (ω,θ), one can analytically show that for constant N, it should takes place the asymptote (Pushkarev et al. 2003;Badulin et al. 2005) whilst, for constant E, it is true (not presented in literature) Details of derivations for ratios (24, 25) can be found in the Appendix (A) of this paper.But here we dwell on manifestations of differences in asymptotes for E(t), N(t), and S p (t) in different version of the KE-solution algorithms, attracting proper our results for isotropic initial spectra.
First, let us consider the case of total wave energy E preservation (as the less studied).If the self-similar one-dimensional spectrum, S sf (ω), is represented by ratio (9), then we have Consequently, the peak of spectrum, S p (t) , should grow as and the total wave action should grow at the same rate: (This growth of N will be discussed below).Our calculations show that under the condition of constant wave energy, the peak frequency of the self-similar spectrum decreases according to the law Within the limits of the numerical estimates scattering, ratio (29) corresponds to (25) rather well.As seen, dependence (34) is remarkably weaker then ( 29), what is in full accordance with the analytics.Despite of ratio ( 34) is slightly weaker than the theoretical dependence, within the limits of estimates scattering, it is rather close to (24).
Herewith, in this case, our calculations show that and what means that the expected theoretical dependences (32, 33) are very well fulfilled.
Thus, in general, the correspondence between numerical and analytical results is rather well.The mentioned deviations of numerical ratios (34) and (29) for ω p (t) from the analytical ones,( 24) and ( 25), are apparently due to the limited range for variability of values ω p in our calculations, determined by grid (13), and to relatively small scales of the considered evolution time.The sophisticated calculation could be executed in a separate study with using the properly chosen grid.
According to theoretical estimates ( 27) and ( 33), the only visually observable effect of the differences in these versions of the KE-calculation algorithms is the different asymptotes of the spectrum peak S p (t).This effect is really observed in our calculations (Figs.7a,b).
In addition to the said, it is interesting to give estimates of the self-similar spectrum values at any fixed frequency where the tail for self-similar shape-function sf p

F ( / )
 is explicitly written.In the case of constant E, from (37) and ( 27), the following dependence should take place For constant N, according to (36, 37), it must be valid  34), if one accounts for the above scattering in the estimates of powers.Pay attention that in the case of constant N, the intensity of spectral tail decreases faster, to support the proper self-similar spectral shape.These estimates complete our description of new peculiarities for the numerical solution of KE, established in this work.

Interpretation of results
The peculiarities of the frequency asymptote of NLT, Nl(ω), and the features of the self-similar spectrum formation (Sects. 3.1, 3.2), which are accompanied by the peak-frequency asymptotes ( 24) and ( 25) and related asymptotes (32,33) and (27,28) (Sect.3.3), do not require any special treatment.The only principal point for treatment is the search for reasons of establishing self-similar spectrum S sf (ω) with the tail of form (21), because it is realised for any version of the algorithm for solution of KE (1).Although this numerical result, S sf (ω) ~ ω -4 , was obtained earlier under the condition of wave-action constancy (Pushkarev et al. 2003;Badulin et al. 2005), here the question arises: why the one-dimensional self-similar spectrum, S sf (ω), always has the tail of form (21), both in cases of preservation of the total wave action and in cases of conservation of the total wave energy.An interpretation of this result is given below.
First, we note that the possibility of separating these algorithms is stipulated by the fact that two integral parameters of the wave-system, N and E, are not simultaneously conserved in the numerical solution of KE (1).This statement was repeatedly noted by Zakharov and co-authors (e.g., Pushkarev and Zakharov, 2000;Pushkarev et al., 2003;Badulin et al., 2005;Zakharov 2017, and references therein).The absence of the simultaneous preservation of N and E for a self-similar shape solution of KE (1) was discussed in Sect. 2.Here we may note some technical and mathematical reasons for nonpresentation of values for N or E in the course of the KE numerical solution.The former include: using restricted numerical grid, systematic errors in the calculation of KI, provided by approximate estimations of contributions from the singular points in the integrand (Polnikov 1989, van Vledder 2006), and insufficient resolution.The mathematical reasons are: the weak convergence of the KI for slow-falling spectra (Zakharov 2017) and typical shortages of numerical solution for the KE.These reasons are not due to the physics of the phenomenon, because KE of form (1) has no external sources or sinks, which could violate the constancy of N or E. Therefore, from the viewpoint of mathematics and physics for the system under study, both versions of the computing algorithm for KE (1) have a reason to be used, at least as a compensating the mentioned errors.As shown above, the differences among these versions of the KE solutions lead to fundamental differences in the formation of fluxes for wave action N and wave energy E. When N is constant, on the scales of long-term evolution, the (quasi-)constant upward energy flux P E (ω) is realised in the tail frequency-range.In turn, if E is constant, the similar downward action flux P N (ω) occurs (Figs.8a and 9b).Thus, two fundamentally different physical situations arise, provided by the choosing the algorithm version.Nevertheless, the tail of one-dimensional self-similar spectrum, S sf (ω,), always takes form (21) (Figs.7a,b).The interpretation of this result requires a separate discussion, and, may be, a separate research.Indeed, according to Zakharov and Zaslavskii (1982), the spectra of forms ( 6) and ( 21) are treated as the Kolmogorovtype spectra of constant energy flux P E (ω) directed upward in frequencies.Therefore, they must be formed only in the version of the KE solution with the condition of total wave action preservation, when proper flux P E (ω) occurs.According to this logic, in the case of the KE-solution with the condition of constant wave energy, the spectra with the tail of form ( 7) should be formed, that is, the tail should take the form S sf (ω) ~ ω -11/3 .If it is so, there are not any problems, and the treatment by Zakharov and Zaslavskii (1982) could be accepted totally.
However, the spectrum tail of form S sf (ω) ~ ω -11/3 is not observed in our KE-solutions with the version of constant E. This result is not related to the accuracy of the calculations, because the estimated numerical errors of our calculations are only 3% to 5% (Sect.2.1), and the errors of estimating the exponents in the asymptotes are not more than 1% to 2%.That gives grounds for distinguishing the '-4' and '-11/3' decay laws.Therefore, the numerical result of invariability of the solution of form (21) should be considered as the quite reliable one.
Our treatment of the matter under discussion is based on the results of estimating the high-frequency asymptote for nonlinear energy transfer Nl(ω), obtained above (Sect.3.1).Indeed, by setting up a standard initial spectrum of waves in form (14-16) with the decay index n > 4, we found that, for such n, the tail of Nl(ω) is always positive and decays rapidly (see Figs. 1a, 2a and the asymptotes in Table 1).Consequently, in the course of spectrum evolution in terms of KE (1), at the beginning stage, the intensity of the spectrum tail increases, resulting in decreasing the decay-power of current spectrum.On the time scales when the tail of the spectrum approaches to form S(ω)~ω -4 , the asymptote of Nl(ω) has already changed its sign and took form (23), i.e., close to Nl(ω) ~ -ω -4 .Mathematically, it is evident that this change of Nl(ω) stabilises the spectral shape, providing the main reason of establishing self-similarity of the spectral shape with the tail of form S sf (ω)~ω -4 .
Pay attention that a choosing the algorithm version does not impact to this fact, as is dose not impact of the tail of Nl(ω).
It should be noted here that, due to the nonlocality of four-wave nonlinear interactions, the asymptote of Nl(ω) is determined not only by the asymptotic behaviour of the spectrum tail but does also by the shape of its peak.The latter is clearly seen from the results of the calculations for run 6 from Table 2, when the initial spectrum is falling with the '-4' law (Figs. 6a,b).This manifests the role of the spectrum-peak domain in formation of the entire shape of self-similar spectrum S sf (ω,θ).Thus, due to reforming both the shape of spectral peak domain and the spectrum tail, in the course of the KEsolution, the coincidence of power-forms for asymptotes for S(ω) and Nl(ω) is realised, which stabilises the shape of the entire spectrum, providing the tail of form (21).Moreover, our numerical experiments for initial spectra with n < 4 showed that, although the KI is not well definite for such spectra, the stabilising process is even faster, because of the tail of Nl(ω), in In turn, according to the theoretical point of view (Zakharov at al., 1992), the wave action has a formal meaning of virtual "wave particles".Thus, the fact of increasing N, in the case of constant energy of the system, could be interpreted as the formal effect of 'condensation' of wave-particles from the tail of the spectrum during the nonlinear wave evolution in the frame of KE (1).
Another situation, relevant to the Kolmogorov-type turbulence, is realised in the case of presence of source and sink external to the wave system, separated in the frequency band.It was modelled in numerous works devoted to numerical solution of the KE of form (8) (e.g., Polnikov 1990Polnikov , 1994;;Pushkarev et al. 2003;Badulin et al. 2005; among others).
However, the detail discussion of the impact of the mentioned algorithms on process for Kolmogorov-like spectra formation in forms ( 6) and ( 7), detailed analysis of these processes, and their reconciliation with the conclusions of this paper requires a further study and separate presentation in future.

Conclusions
Summarising the results above, we can draw the following conclusions.
First, the features of the high-frequency asymptotes for nonlinear energy transfer Nl(ω), which are close to a power-law form, play a fundamentally important role in understanding the evolution features of nonlinear gravity waves in water, governed by the four-wave KE of form (1).In particular, the behaviour of these asymptotes indicates the fact of nonlocality of nonlinear interactions under consideration.The change of the sign for Nl(ω) and its asymptote to form Nl(ω) ~ -ω -4 , associated with a change in the decay-law of spectrum tail to form S(ω) ~ ω -4 , is responsible for maintaining the tail of spectral shape in form S sf (ω)~ ω -4 .Second, the shape of the peak-domain for two-dimensional spectrum S(ω,θ), the characteristics of which are presented in Table 2 and in Figs. 3a,b,c, is responsible for the formation of self-similar two-dimensional spectrum S sf (ω,θ) as a whole, what is clearly visible on the example of calculations for run 6 from Table 2, for the initial spectrum with n = 4 (Figs.6a, b).Third, self-similar two-dimensional spectrum S sf (ω,θ) with the tail of form S sf (ω) ~ ω -4 , supported by a nonlinear energy transfer with the tail of form Nl(ω) ~ ω -4.15 ± 0.05 , is established in solutions of KE (1) for any initial spectral forms (14-16), in the versions of both preservation of total wave energy E and constancy for total wave action N (Figs.7a,b).As a result, the establishing self-similar spectra with the tail S sf (ω) ~ ω -4 is the purely mathematical property of the four-wave KI, that was analytically found and proved by Zakharov and Filonenko (1996).
The (quasi-) constant downward fluxes of wave action P N (ω) and the upward energy fluxes P E (ω), taking place in the spectrum-tail range, are not the reasons for self-similar spectrum formation in the course of solving KE (1), because of the deviation of this situation from the Kolmogorov's one.The differences between the frame of KE (1) and Kolmogorov-type turbulence are: a) the presence of a distinguished frequency scale in the system, provided by the peak of spectrum, b) the nonlocality of four-wave nonlinear interactions formatting the entire spectrum, and c) the absence of sources and sinks external to the wave system.
. Later we use KE (1) in the form written for the wave-energy spectrum, S( , )  , given in the frequency-angular space, ( , ) have found the analytic solution of KE (1) of form the KI identically to zero (i.e.,

Nonlin.
Processes Geophys.Discuss., https://doi.org/10.5194/npg-2018-35Manuscript under review for journal Nonlin.Processes Geophys.Discussion started: 9 August 2018 c Author(s) 2018.CC BY 4.0 License.Nevertheless, studying the features of the KE of form (1) requires its continuation.Namely, asymptotes of the nonlinear energy transfer, frequencies, both at the initial stage and on long-term scales of the spectrum evolution, is not yet described.The quantitative characteristics of the two-dimensional self-similar spectrum shape, ω 1 = 0.64 rad/s, q = 1.05, 1 ≤ i ≤ I = 90, and Δθ = 5 o .(12b) Nonlin.Processes Geophys.Discuss., https://doi.org/10.5194/npg-2018-35Manuscript under review for journal Nonlin.Processes Geophys.Discussion started: 9 August 2018 c Author(s) 2018.CC BY 4.0 License.The initial high-frequency asymptote of NLT is determined from the results of calculating KI at the first time-step of the KE solution.To study the asymptote of NLT and spectral shape resulting from numerical solutions of KE on long-term scales for spectrum evolution, we have used the truncated frequency grid and lower angular resolution: 0.64 ≤ ω ≤ 7 rad/s with I = 50, Δθ = 10 0 , (13) Nonlin.Processes Geophys.Discuss., https://doi.org/10.5194/npg-2018-35Manuscript under review for journal Nonlin.Processes Geophys.Discussion started: 9 August 2018 c Author(s) 2018.CC BY 4.0 License.
Geophys.Discuss., https://doi.org/10.5194/npg-2018-35Manuscript under review for journal Nonlin.Processes Geophys.Discussion started: 9 August 2018 c Author(s) 2018.CC BY 4.0 License.variations of parameter a) for fixed n.Indeed, according to Table1, despite of the tail of spectrum S (ω) is the same, the tail of Nl(ω) is different for different shapes of the spectral peak.It means: the nonlinear interactions have an explicit non-local feature.This conclusion has already been noted in literature(Hasselmann and Hasselmann 1985;Polnikov 1989;van Vledder 2006); however, the non-locality of the four-wave nonlinear interactions is expressed in the asymptotes of NLT in the clearest manner.As shown below, the non-locality of NLT is very important in explaining the reasons of the self-similar spectral shape formation.

Figure 3 :
Figure3: (a) Typical shape of the self-similar function for the angular directivity, А(ω).Line 1 is the mean for all anisotropic runs (7 through 11) from Table2; line 2 is run 10; line 3 is run 11.(b) Self-similar shapes S sf (ω/ω p ) for run 3 from Table2(isotropic spectrum).Line 1 is the shape of J-spectrum according to (14-16) (the same run); line 2 is S sf (ω/ω p ) at time t ≈ 1·10 4 s; line 3 is S sf (ω/ω p ) at the time t =1·10 5 s.(c) The difference of shapes for normalised self-similar spectra S nor (ω/ω p ) in two cases.Run 4 from Table2.. Detailed information regarding the features for spectrum evolution is shown in Figs.5-8 (using the versions of KE-
Geophys.Discuss., https://doi.org/10.5194/npg-2018-35Manuscript under review for journal Nonlin.Processes Geophys.Discussion started: 9 August 2018 c Author(s) 2018.CC BY 4.0 License.16Here, it is important to emphasise that the results presented above were repeated for versions of solutions for KE (1) both with preservation of total wave energy E and total wave action N too (e.g., Figs.7a,b below).Moreover, the same results were found in the algorithms without the exact preservation of N or E, when the proper control was not carried out.(a)(b) Figure7: Time-history S(ω,t) for two versions of the KE-solution (run 4): (a) The version with condition E = const.Line 0 corresponds to t = 0 s; line 1 to t = 1830 s; line 2 to t = 0.5•10 5 s, line 3 to t = 3.5•10 5 s, line 3' to the tail part of line 3 (with a weight of 0.5), line 3'' to the trend for line 3' (equation y= 0,135•x -4.06 ).

Figure 8 :
Figure8: (a) Time-history of flux P Е (ω,t) for run 2 from table 2 (when N is constant).Line 0 corresponds to t = 0 s; line 1 to t = 700 s; line 2 to t = 2500 s.(b) Time-history of flux P N (ω,t) for the same run 2. Lines correspond to the same time moments.action N is constant, a quasi-constant tail for upward energy flux P Е (ω) is established in frequencies range ω > 2ω p (t), on

Figure 9 :
Figure9: (a) Time-history of flux P Е (ω,t) for run 2 from table 2 (when E is constant).Line 0 corresponds to t = 0 s; line 1 to t = 250 s; line 2 to t = 1830 s.(b) Time-history of flux P N (ω,t) for the same run 2. Lines correspond to the same time moments.
. They are determined by the ratios (for any version of the algorithm for the KE-solution)

Table 1 .
Asymptotes of the NLT at the first time step of KE solution.
(b) Form of NL (ω) at high frequencies, run 1 from table

Table 2 .
Asymptote of Nl(ω) and shape parameters of S sf (ω,θ) in the long-term KE solution.