A Continuous Non-ergodic Theory for the Wave Set-up

Inhomogeneities in the wave field due to wave groups, currents, and shoaling among other ocean processes can affect the mean water level. In this work, the classical and unsolved problem of continuously computing the set-down and the following set-up induced by wave breaking on a shoal of constant finite slope is tackled. This is possible by using available theoretical knowledge on how to approximate the distribution of wave random phases in finite depth. Then, the non-homogeneous spectral analysis of the wave field allows the computation of the ensemble average by means of the phase distribution and the inversion of the integral of the second moment for the special case of a shoaling process with uniform phase distribution. In doing so, I am able to obtain a direct effect of the slope magnitude on the phases distribution. Therefore, an analytical and slope-dependent mean water level with continuity over the entire range of water depth is provided.


I. INTRODUCTION
Conservation principles in the physical sciences often have significant consequences for the solutions of governing equations [1,2].However, discussions on the energy partition in wave hydrodynamics have risen very tardily [3][4][5][6][7][8].Not surprisingly, a strong debate surfaced over the proper formulation for the momentum of water waves [9][10][11][12] to understand the phenomenon of longshore currents [13][14][15][16].Other nearshore phenomena were concurrently discovered [17,18], most prominently the change in mean water level outside the surf zone (set-down) followed by a set-up within this zone [19,20].Experimental evidence supported the radiation stress theory [21,22] as it could seemingly explain all these phenomena.Computation of the set-down/set-up has important consequences for adjacent water wave processes.For instance, while studying the design of breakwaters Hunt Jr. [23] showed that the run-up of waves over a beach depends on the set-up.Furthermore, rip [24][25][26] and longshore currents [27,28] are also influenced by the alongshore variations in setup.Additionally, rogue wave occurrence in the transition between deep, intermediate and shallow waters has recently been linked with the analytical effect of slope on non-homogeneous evolution of irregular waves [29].Remarkably, the change in mean water level can be achieved by any ocean process in which the energetic and momentum balances are disturbed, such that the classical view of the wave group-induced set-down over a flat bottom [22] is complemented by shoaling-induced [30,31] and current-induced counterparts [32].
The set-down is usually computed from the Bernoulli equation [33], thus the sub-harmonic is an a posteriori term due to the conservation principles applied to a second-order irregular wave field.However, as discussed in Mei et al. [33], the typical calculation neglects the depth gradient, and thus overlooks the effect of slope.In addition, there is no analytical way to obtain a transition between set-down and set-up, as the latter occurs following wave dissipation and breaking which can not be applied to the Bernoulli equation.Experiments assessing the continuous evolution from wave group-induced set-down transitioning to shoaling-induced set-down until dissipation is strong, and a wave set-up appears have been conducted since the 1960s [34].Despite the elapse of half a century, no continuous theoretical model has been developed to reproduce such an evolution other than the piecewise formulation based on the radiation stress [21,22].Indeed, Battjes [35] provided the best empirically-driven closed-form model for the set-up, but does not overcome the piecewise formulation between setdown and set-up.Likewise, although Hsu et al. [36] generalized the set-down and set-up calculations of McDougal and Hudspeth [28] for oblique waves under the effect of refraction and Massel and Gourlay [37] generalized the problem to include bottom friction, these two examples show that theoretical and numerical advances have still not been able to overcome the piecewise decomposition.
The field of wave statistics is often treated as a mere consequence of the solutions to the governing equations of hydrodynamics, and fundamental ocean processes are not believed to be affected by extreme wave occurrence or joint probability densities of surface elevation and random phases.In this work, I follow the footsteps of Mendes [38] and demonstrate that wave statistics can indeed affect the fundamental physical process leading to the set-down and set-up, whose dynamical evolution can not be treated by a unique approach due to wave breaking.
compute integral quantities and their conservation such as mass, momentum or energy flux [5,7].Solving the Bernoulli equation leads to a clear dependence of the setdown on the mathematical form of the solution for the velocity potential [33]: (1) The velocity potential up to second-order in steepness can be written as [39,40]: with notation θ = k(z + h), Λ = kh and ϕ = kx − ωt.If ∇h ≡ ∂h/∂x = 0 is assumed, the horizontal component of the velocity vector reads (∂ϕ/∂x = k, Likewise, noting that ∂h/∂z = 0 by definition and ∂θ/∂z = k, the vertical component of the velocity reads: For the next step of taking the time average of the square of the velocity components, the reader is reminded of periodic averaging of trigonometric functions (see eq. ( 12)): lim for all n ∈ N * .Through integration by parts, one has as a corollary for all (m, n) ∈ N * : ⟨sin 2n+1 ϕ cos 2m+1 ϕ⟩ = ⟨sin 2n ϕ cos 2m+1 ϕ⟩ = 0 , (6) where the operator ⟨•⟩ denotes time averaging.The square of the velocity components will have only two nonvanishing terms ⟨sin 2 ϕ⟩ = ⟨cos 2 ϕ⟩ = 1/2.Therefore, by means of eqs.(1,(3)(4) the set-down may be computed over a flat bottom: which taking into account the hyperbolic identity cosh 2 θ − sinh 2 θ = 1 and the leading order in the dispersion (ω 2 = gk tanh Λ) leads to the formula: As far as the author is aware, the term of second order in steepness inside the brackets does not appear in the literature, likely being neglected without further discussion.In the limit of second order theory through the Ursell number Ur ⩽ 8π 2 /3 [25], the term inside brackets increases the overall set-down by not more than 15%.However, beyond the limit of the second-order theory, the additional term will cause the set-down to diverge.Thus, this is a symptom of the unsuitability of this approach in computing the set-down/set-up over continuous range in relative water depth.Alternatively, one may compute the gradient of the set-down from the horizontal gradient of the radiation stress, as the momentum balance equations lead to the radiation stress formula [21,22]: where the cross-shore component of the radiation stress is a function of the ratio between phase and group velocities: where p = p 0 + ρgζ = ρg(ζ − z) is the pressure field underneath the waves.With a little algebra, limited to the first order in steepness this integration leads to: .
This procedure verifies that eq. ( 8) up to the first order is the solution of eq. ( 9).Hence, there is no advantage in integrating the latter equation as compared to perform derivatives in eq. ( 1).Furthermore, the set-down computed above was derived without any need for shoaling formulae, which suggests there is little difference between wave-group set-down driven by difference in wave heights and shoaling of a very mild slope ∇h ≲ 1/100.Noteworthy, waves may break in the region of validity of secondorder waves, and as soon as the waves reach the surf zone and are subject to strong dissipation the solution to eq. ( 9) is no longer the equivalent of solving Bernoulli's equation.Nonetheless, the set-down caused by the crossshore component of the radiation stress can not be properly formulated in the surf zone even in eq. ( 9), except if one assumes that the wave height stays constant within this region [22,34,35].The major issue tackled by this work is to find a continuous analytical way to express the mean water level in these two physically different zones.This can not be done either through eq. ( 9) nor eq. ( 1).Hence, waves of second-order in steepness follow an approximated piecewise formula for the set-down and set-up, see for instance Chen et al. [31].In the next sections, I attempt to describe both set-down and set-up continuously over a plane beach within a single physical principle applicable to both zones.

III. RANDOM PHASE DISTRIBUTION AND WAVE STATISTICS
Although wave dissipation typically separates physical theories of ocean wave mechanics between prior to breaking and after breaking, statistical measures of irregular waves do not suffer from this dynamical problem, at least empirically.For instance, Glukhovskii [41] showed that is possible to connect in a continuous range for relative water depth (k p h) in terms of the ratio H s /h the distributions of wave heights in deep and shallow water, otherwise known to be restricted to piecewise solutions (See Wu et al. [42], Mendes and Scotti [43] and Karmpadakis et al. [44] for distributions dependent on this ratio).The latter ratio is often found to be a good proxy for wave breaking [35,45].As such, in this section I invoke the practical knowledge of wave statistics and reveal how the non-ergodicity of irregular water waves plays a role in oscillating the mean water level from its initial condition (without waves) continuously as they approach a plane beach.Let the ensemble average of a random variable X(t) be defined: where f (X) is the probability density of the random variable X(t).Then, a stochastic process is said to be ergodic if E [X(t)] = ⟨X(t)⟩ holds.Now, it is of interest to compute ensemble averages of powers of the surface elevation ζ n for n ∈ N * .The mean water level can be easily computed in the case of linear waves [46], and it will always leads to an ergodic process over a flat bottom with E[ζ] = ⟨ζ⟩ = 0. Through a change of variables [47] and applying the law of the unconscious statistician [48] to eq. ( 12), I compute the ensemble average as follows: Eq. ( 13) delineates how the probability density of the surface elevation is transformed into the distribution of random phases in computing the very same ensemble average.In the case of linear waves, ergodicity is corollary of an uniform distribution of phases (or alternatively, a Gaussian distribution of the surface elevation): Likewise, the second moment (the spectral energy density) also features ergodicity: = lim which is well-known to be related to the Khintchine [49] theorem.Consequently, it can be demonstrated as in appendix A that any signal which can be described as a two-dimensional random walk, is zero mean with finite variance and independent random variables, will result in a Rayleigh distribution [50,51].For a review, see section 7.7-2 of Middleton [52] and 6.2.2 of Mendes [53].However, if the distribution of random phases is not uniform the envelope will deviate from the Rayleigh law.Indeed, when waves travel over a shoal the distribution of surface elevation, crests or crest-to-trough heights are strongly deviated from the normal law [54][55][56].The latter can be understood as random phases no longer being uniformly distributed, as demonstrated by Bitner [57].The above statistical formulation computes the mean water level from the time average as if the shoaling process was ergodic.Although the latter is not factual, for a variety of purposes the ergodic approximation is useful and greatly simplifies calculations.Compared to the deterministic approach of section II, the possible generalization of the statistical approach in eq. ( 14) has the advantage of continuously assessing the mean water level in the surf zone.

IV. NON-ERGODIC CONTINUOUS SET-UP
In this section I attempt to compute the ensemble average from the random phase approach.Consider the case of spatial inhomogeneity due to a shoal.In this case, the time series is approximately weakly stationary.The evolution of the non-homogeneous spectrum is better formulated as a homogeneous spectrum corrected by a term Γ = ⟨ζ 2 ⟩/E that absorbs such inhomogeneity, where E is the spectral energy density factoring out the linear energy (1/2)ρga 2 .Up to second order in steepness, one would find a surface elevation in its simplest form [40]: (16) In the notation of Mendes and Kasparian [29] I can rewrite it as The Γ correction simply reads [58]: where ε = H s /λ denotes the irregular wave measure of significant steepness, with H s being the significant wave height and λ the zero-crossing period.The effect of the Γ correction on the moments of ζ in a uniform distribution of phases up to second order follows (note that the Airy case has E[ζ 2 ] = 1): If linear wave theory is assumed (ζ = a cos ϕ) but with slightly non-uniform phase distribution affected by 2πf (ϕ) = 1 − πε 4 √ χ1 cos ϕ, µ 2 is recovered to second order.This suggests one can estimate the magnitude of the perturbation over the otherwise uniform distribution 2πf (ϕ) − 1 ∼ − πε 4 √ χ1 cos ϕ.Therefore, a generalized distribution of the form is sought: where (Ξ 1 , Ξ 2 ) are coefficients to be uniquely determined by wave theories.As a remark, from the point of view of the probability distribution [58] the second moment has been normalized by the energy, now modified by a second-order correction.In this case, the true second normalized moment is µ 2 /E = Γ.Likewise, the first moment reads

A. Computation through Gram-Charlier Series
To study non-uniform distributions of random phases one must obtain the joint probability density of both phases and envelope (the analytical equivalent of the wave height).The non-uniformity of random phases is therefore closely related to the non-Gaussianity of the surface elevation probability density.Historically, deviations from Gaussian law are dealt through asymptotic expansions of the central limit theorem, having been rigorously defined since Laplace [59] with further refinements introduced by Chebyshev [60], Berry [61] and Esseen [62].The approximation through Gram-Charlier or Edgeworth series of the Gaussian distribution is therefore not a refinement of the Gaussian distribution, but rather an approximated distribution when this limit has not been reached.Because the joint probability density of the envelope of the surface elevation ζ and its Hilbert transform ζ can be well approximated by the expansion of third-order cumulants in terms of Hermite polynomials, Tayfun [63] showed that it can lead to the skewness µ 3 correction in deep water: Integrating the joint density over the envelope R = (ζ 2 + ζ2 ) 1/2 with respective pair (ζ, ζ) = (R cos ϕ, R sin ϕ), one finds the phase distribution in deep water [63]: In deep water, we typically have µ 3 ⩽ 0.6 [64,65], such that the peak of the phase distribution does not exceed the uniform distribution 2πf (ϕ) in 10%.By comparing eq. ( 20) with eq. ( 22) and the inequality of eq. 14 of Tayfun [64], we find the bound: In intermediate and shallow waters, especially in unsteady conditions similar to the experiments of Trulsen et al. [54], this departure from the standard uniform distribution is expected to be much larger, see for instance figure 13 of Bitner [57].However, approximations of Gram-Charlier series are reduced to either skewness or kurtosis effects, see Mori and Janssen [66] for the latter.Nevertheless, the full computation of a Gram-Charlier series is dependent on both moments of the surface elevation [67].Together with the explicit effect of bandwidth ν [7] from eq. 20 of Tayfun [64], the impact of both skewness and kurtosis as computed in eq.50 of Tayfun and Lo [68] reads: Powered by TCPDF (www.tcpdf.org) FIG. 1: Contour of the normalized set-down/up as a function of pre-shoal steepness ε 0 .
For a broad-banded sea state with typical JONSWAP spectrum of peakedness parameter γ = 3.3 the bandwidth is of the order of ν ∼ 1/2.By means of the relation between skewness and steepness [68], I find: Furthermore, since eqs.19-22 of Mori and Kobayashi [69] show that the kurtosis is computed as µ 4 = 48(ka) 2 in deep water (see appendix A of Mendes and Kasparian [70]), I further simplify the phase distribution as follows: which adds another 0.4% deviation to the maximum of 10% of the first term containing cos ϕ.Hence, the approximation of eq. ( 22) is justified in deep water.In intermediate water, the importance of the second term can increase by fivefold due to µ 3 ∼ 1.In order to express the phase distribution in terms of relative water depth and wave steepness, I adopt the finite-depth formula of eq.22 of Mori and Kobayashi [69].Taking into account the ratio between the coefficient in eq. ( 16) and those of Mori and Kobayashi [69] amounts to √ χ1 /(D 1 +D 2 ) ∼ 2, the phase distribution fully reads: where the vertical asymmetry between crests and troughs 1 ⩽ S = 2a/H ⩽ 2 is added to correct an otherwise symmetrical approach from the beginning.In comparison with eq. ( 20), the above solution is equivalent of finding Ξ 1 = 3/2S ≈ 1.25, upholding the lower bound in eq. ( 23).With the exact distribution for random phases in hand, I compute the change in mean water level: 284,000 In figure 1 the theoretical evolution of the normalized mean water level is displayed, decreasing evermore towards shallow water until it reaches its global minimum at the surf zone.Then, the mean water level starts to increase and reaches the plunging point at µ 1 = 0 and quickly increases to a normalized set-up, as expected from observation [71].Figure 1 also points to the fact that although increasing the pre-shoal steepness leads to an earlier peak in set-down, the magnitude of this peak seems to vary weakly with the steepness.This is  28) and the classical approach of eq. ( 8) adapted to irregular waves in eq. ( 29 tive water depth as a function of curves of fixed steepness according to the present theory and steepness-limited wave breaking [72] with and without vertical wave asymmetry. in qualitative agreement with the linear term of the setdown formula of eq. ( 8), because the latter converges to −a 2 /4h in the surf zone.Therefore, the new formula recovers both sides of the piecewise theory of Longuet-Higgins and Stewart [21], and being continuous explains theoretically when and how fast the transition between set-down and set-up occurs.
In order to compare the set-down computed from the random phase distribution with the radiation stress theory, eq. ( 8) is rewritten for irregular waves by using the same transformation from eq. ( 16): (29) Thus, figure 2 compares the non-ergodic set-down from the ensemble average with the time average from eq. ( 29).The set-down computation based on the radiation stress of waves over a flat bottom underpredicts the result of the present model, and this is expected because the latter is based on the Γ model that is valid for relatively steep slopes.This is in agreement with the fact that Longuet-Higgins and Stewart [21] underpredicts the set-down in the experiments performed by Saville [20] for steep slopes [31,73,74].
In addition, it seems that the point in which wave breaking occurs, which in the present model first moment plot corresponds to the transition between set-down and set-up, is well predicted by the mean.To find the breaking point one uses ε ⩽ tanh (kh)/7 and thus find kh ⩾ tanh −1 (7ε) according to Miche [72].Factoring out √ E which is of the order of O(1), The breaking point can be obtained by setting dµ 1 /d χ1 = 0, finding to leading order: Following the definition of the trigonometric coefficient χ1 in eqs.(16)(17)(18)(19), I solve the cubic equation in tanh k p h and find the critical relative depth for S ≈ 1.2: FIG. 4: Amplification in the classical set-down formula due to the slope effect of eq. ( 38).
I can also obtain the location of the plunging/spilling point by solving µ 1 (ε, (k p h) c ) = 0 which implies χ1 = 20/π 2 ε 2 S 2 : Figure 3 compares steepness-limited wave breaking critical points with the breaking and plunging points as estimated by the present theory, showing qualitative agreement.However, the model shows that depth-limited wave breaking turns the transition between set-down and set-up earlier than inferred from Miche [72] at low wave steepness.

V. SLOPE EFFECT
To fully compute the effect of steep slopes on the mean water level, one must start with its effect on the derivative of the velocity potential.The algebra can be split into slope-dependent and slope-independent terms: Hence, the slope-dependent term of the horizontal velocity component will lead to a modification in eq. ( 1): with the term ⟨ζ⟩ having been already obtained in eq. ( 8).
The complementary part of the horizontal velocity is computed by taking derivatives ∂/∂x on the hyperbolic functions in eq. ( 3) since ∂θ/∂x = ∂Λ/∂x = k∇h: and thus one can obtain: Because of eq. ( 6), the time average of u 0 ∆u vanishes.Although the hyperbolic coefficients appearing in (∆u) 2 are bulky and can not be simplified I may take the limit at z = 0 for these coefficients due to ⟨u 2 0 − w 2 ⟩/2g ∼ ⟨(∆u) 2 ⟩/2g < 0.1.In that case, it is straightforward to show that lim z→0 H 1 = 1 and lim z→0 H 2 = −2 sinh 4 Λ/ tanh 3 Λ.Therefore, I may approximate: Accordingly, up to second order in steepness ka and under the effect of an arbitrary slope ∇h without curvature (∇ 2 h = 0), the wave-driven set-down is fully computed at last:  In figure 4 the disparity between mild slope and steep slope mean water levels is displayed.Taking into account the magnitude of the maximal correction due to the terms proportional to (ka) 2 /(kh) 6 approaching the Ursell limit in shallow water, I can further simplify the slope-dependent set-down: The above slope correction is the first of the kind to be found analytically through the classical textbook methodology, albeit the dependence on the slope has been found previously through numerical simulations [74] and experiments [30].Although Chen et al. [31] have also found a dependence on (∇h) 2 for the shoaling-induced set-down, this comes from an a priori expansion of the velocity potential in powers of ka • ∇h and the formulae for the computation contains dozens of hyperbolic and trigonometric coefficients that make such a result computationally cumbersome.Notably, the present derivation does not apply to realistic reflective beaches, so the above formula is limited to cases with a negligible reflection coefficient K R ≈ 0.1ξ 2 , hence with small surf similarity ξ ∼ ∇h/ √ ε ≲ 1 [75].This limitation inhibits an otherwise divergence in the set-down due to a step ∇h → ∞.On the other hand, one may examine the maximum set-down beyond second-order theory in the limit when sinh kh ≈ tanh kh ≈ kh: Therefore, the maximum set-down using linear theory is amplified by a percentage of (50∇h) 2 % due to steep slopes, which for negligible shoal reflection lies in the range 25% − 250% for 1/10 < |∇h| < 1/3.This amplification magnitude is at par with two to threefold larger set-down than predicted by Longuet-Higgins and Stewart [21] as compared to experiments in Saville [20].However, the set-up can not be computed from eq. ( 39) and no continuous formulation can be extended.In the next section, I compute and show that the same shape of the set-down zone can be obtained, although it is more complex algebraically but with the advantage of computing the mean water level continuously from set-down up to set-up.

A. Non-ergodic Approach
In order to extract the slope dependence from the set-down, first one must generalize the random phase approach for eq.( 19).As discussed in Mendes and Kasparian [70], the shoaling slope effect induces a correction √ ∇h to the existing excess kurtosis of steep slopes, thereby decreasing the non-Gaussianity as the slope tends to zero.Combination of the "uniform" slopedependent phase distribution with the inhomogeneous slope-independent phase distribution of eq. ( 28) leads to (see appendix B): Naturally, the generalized slope-dependent distribution recovers both inhomogeneous and uniform distributions over mild slopes.Such distribution allows us to estimate the effect of arbitrarily steep slopes on the set-down, as an alternative to the computation of the Bernoulli equation (or radiation stress) leading to eq. ( 39).As such, the setdown driven by the equivalent inhomogeneous random phase distribution of eq. ( 28) now reads: As observed in figure 5, as the slope becomes milder (|∇h| ⩽ 1/100) the non-ergodic set-down converges to FIG. 6: Amplification in the mean water level formula due to the slope effect of eq. ( 42) as compared to a fixed mild slope of |∇h| = 1/1000.
the adiabatic classical set-down of Longuet-Higgins and Stewart [21].Comparing it with the slope-dependent classical set-down computation introduced in eq. ( 38), the magnitude of the global minimum is similar, but the shape of the set-down evolution is different.Moreover, the slope-dependent classical set-down only deviates from the adiabatic one near wave breaking, thus a shortcoming.Naturally, the present formula is superior as it is also capable of computing the set-up, in particular also describing its dependence on the slope.Lastly, figure 6 describes how a steep slope amplifies the set-down, but this growth is significant only in a narrow range of relative water depth 0.6 ⩽ k p h ⩽ 0.9.

VI. CONCLUSION
In this work the radiation stress-led computation of the set-down has been extended to include the effect of slope.Furthermore, I have shown that statistical reasoning removes the difficulties in dealing with wave dissipation for the computation of the mean water level when calculated with a non-uniform distribution of random phases, being general for steep slopes as well.By comparing the two formulations, it can be explained why the classical calculation underpredicted the set-down over steep slopes, being recovered by the new slope-dependent model for slopes milder than 1/100.However, the radiation stress can still not lead to a continuous formulation between set-down and set-up, whereas the stochastic approach succeeds in this matter.Despite the usage of a simple model for the distribution of random phases, this work demonstrates that wave statistics are useful for computing primary variables of the water wave solutions, paving the way for the realization that physically visible effects can be handled by the intangible randomness measures of a sea state.
Declaration of Interests.The author reports no conflict of interest.
The second term is at least one order of magnitude smaller than the first, and I can neglect it.I compute the exponential of the potential energy variation, obtaining an expansion for the assumed sea conditions (400ε 2 ∼ 1): Expanding the exponential 20 ⩽ µ 4 , 0 /π 2 ε 2 ⩽ 40 (see figure 9), I conclude the proof for the second-order small amplitude wave theory in steepness (ε ⩽ 1/20): It can be seen in figure 9 that for low steepness one should probably use the coefficient lower bound of 40 to compensate for the gap between solid and dashed curves, whereas for the bulk of the range in steepness the upper bound is the best approximation.

Eq. ( 28 ) -ε = 1 FIG. 2 :
FIG.2:(a) Set-down/up normalized by significant wave height as a function of water depth as computed from the non-ergodic formula eq.(28) and the classical approach of eq.(8) adapted to irregular waves in eq.(29) for a fixed steepness ε = 1/20.(b) Excess set-down in percentage of H s by the non-ergodic computation as compared to the classical counterpart with varying pre-shoal steepness ε 0 and relative water depth.