# Pressure analysis in nonlinear waves by revisiting the breaking wave onset

Yves-Marie Scolan, Stéphane Etienne

► To cite this version:

Yves-Marie Scolan, Stéphane Etienne. Pressure analysis in nonlinear waves by revisiting the breaking wave onset. European Journal of Mechanics - B/Fluids, 2023, 101, pp.246-256. 10.1016/j.euromechflu.2023.06.002. hal-04193152

HAL Id: hal-04193152

https://ensta.hal.science/hal-04193152v1

Submitted on 1 Sep 2023

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons CC BY 4.0 - Attribution - International License
Pressure analysis in nonlinear waves by
revisiting the breaking wave onset.

Yves-Marie SCOLAN(1), Stéphane ETIENNE(2)

(1) ENSTA Bretagne IRDL UMR 6027, Brest, France

(2) Ecole Polytechnique de Montréal, Canada¹

**Abstract**

Breaking waves are numerically simulated by solving the two dimensional nonlinear free surface boundary conditions in potential theory. Onset of breaking wave is determined for two configurations : a focused wave in intermediate depth and a dambreaking flow leading to a soliton. Two criteria are identified: first, a kinematic criterion follows from the ratio of the horizontal velocity at the crest to the phase velocity; second, a dynamic criterion is formulated in terms of the vertical pressure gradient. Thresholds of these criteria are proposed in the literature. For the two configurations studied here, it is shown that these thresholds are not universal but they are reasonable estimates of the onset of breaking wave. To facilitate this study, the Hessian matrix of the pressure is analyzed. It is shown that there exists a continuous curve that links the crest of the wave to the sea bottom. This line follows from the colinearity of the pressure gradient and one of the two eigenvectors of the Hessian matrix. The analysis of the pressure, its gradient and its second derivatives along this line gives more insights into the definition of the thresholds associated to the onset of breaking wave.

*Keywords:* Potential flow, breaking wave onset, pressure analysis.

This paper is dedicated to Lili (Olivier) Kimmoun who was a great friend and a great researcher as well. She participated to this paper by providing us in May 2022 with the focused wave data.

**1. Introduction**

It is well known that there exists a kinematic criterion that determines the onset of wave breaking. The parameter is usually expressed as the ratio of the horizontal velocity at the crest to the phase velocity. Whatever the type of breaker and the water depth, its lowest value is found to be 0.85 but it can be higher. This value depends on the point where the fluid velocity is calculated. It is not the purpose here to discuss the contributions of all researchers in that field. The most recent papers already yield pertinent comments on these results; among them Barthelemy *et al.* (2018), Derakhti *et al.* (2020) or Varing *et al.* (2021).

Much earlier, Stokes (1847) followed by Havelock (1918), Davies (1951) determine the limit of the inner angle in the crest of a regular wave before breaking. This angle is 120°. As

Email address: yves-marie.scolan@ensta-bretagne.fr, stephane.etienne@polymtl.ca
(Yves-Marie SCOLAN(1), Stéphane ETIENNE(2))

(1) ENSTA Bretagne IRDL UMR 6027, Brest, France

(2) Ecole Polytechnique de Montréal, Canada)
a consequence Longuet-Higgins (1963) shows that in a progressive wave "*the acceleration has a magnitude* $\frac{1}{2}g$ *and is directed everywhere outwards from the vertex.*" This also means that in a close vicinity of the crest, the nondimensional vertical gradient of the pressure $p_{,y}/(\rho g)$ is negative and close to $-\frac{1}{2}$ (see Nadaoka et al., 1997) with $\vec{g} = -g\vec{y}$. It is less known that the analysis of the spatial and temporal variations of the pressure $p$ also offers some insights into the conditions that precede the breaking of a wave. The present paper focuses on this analysis.

As the pressure is nil at the free surface, the pressure gradient has only a normal component which is negative all along the free surface if the normal vector points outside the fluid. Inside the fluid, the pressure is positive and then vanishes more or less abruptly at the free surface. In a recent paper (see Scolan, 2023), it is shown that the Hessian matrix<fnref n="1" /> of the pressure has characteristics that enable a clear analysis about the crest of the wave where breaking is expected. In particular, a continuous curve can be determined inside the fluid that links the wave crest to the sea bottom. This line is the collection of the points where one of the two eigenvectors of the Hessian matrix is parallel to the pressure gradient. The retained eigenvector (denoted here $\vec{v}_1$) gives also the direction along which the curvature radius of the pressure surface (denoted $r_1$) may change sign. The other eigenvector is associated with a curvature radius that never changes sign and remains negative. The nondimensional vertical pressure gradient $|p_{,y}|/(\rho g)$ decreases almost monotonically from unity at the sea bottom to the crest. However in the close vicinity of the crest, it has a minimum close to $\frac{1}{2}$ if the crest breaks. It is also shown that the curvature radius $r_1$ is mainly positive along this line but changes sign in the vicinity of the crest, meaning that the determinant of the pressure Hessian matrix (known as the Gaussian curvature) is positive and consequently the pressure vanishes very smoothly at the crest. The present paper aims to describe the spatial and temporal variations of the pressure in the crest when reaching its limit of stability.

Two configurations are studied in the next two sections. Section §2 is concerned with the last stages of a focused wave that may (or not) break in intermediate depth. The second configuration, described in section §3, is a dambreaking flow leading to a soliton that may (or not) break in shallow depth. Conclusion and perspectives of future works are drawn in section §4. The numerical model and the theoretical developments regarding the pressure field are described in Appendix A.

## 2. Focused wave in intermediate depth

The wave train that leads to a focused wave is yielded by a Boussinesq model (see Kimmoun and Scolan, 2009). Figure 1 shows the free surface elevation and the distribution of the velocity potential along the horizontal position from the left wall to the right wall of the tank. The tank is rectangular with length $L = 15m$ and the mean water depth is set to $h = 0.7m$. In order to cover parametrically the breaking and the non breaking configurations, the initial free surface elevation measured from the mean water level is multiplied with a factor; this factor is denoted $\chi$. The free surface profile shown in figure 1 corresponds to the factor $\chi = 1.9205$. This value is the highest scale factor before breaking. Only the free surface elevation is affected; that means that the initial potential energy is varying whereas the initial kinetic energy hardly varies with the factor $\chi$. Figure 2 shows the successive free surface profiles starting from the initial one up to $t = 3.5s$ of time simulation. The corresponding scale factor is $\chi = 1.9205$. It is worth noting in figure 2 the change of asymmetry of the

<footnote n="1">¹The governing equations and all numerical developments are detailed in Appendix A.</footnote>
crest over the spatial interval $x \in [8 : 13]$ in meters. This is consistent with the observations made in Banner *et al.* (2014) or Fedele *et al.* (2020).

The crest is always the location where the velocity is maximum. The breaking is simply detected when the front of the wave becomes vertical. Figure 3 shows the free surface profiles for a breaking case with the factor set to $\chi = 1.95$. A closer view is made in the area drawn by a rectangle. A vertical front is clearly detected before breaking.

The parametrical study is achieved in terms of the factor $\chi$. For each simulation, the time variation of the velocity $u$ in the fluid at the maximum free surface elevation is computed. The phase velocity $c$ is determined by tracking the velocity of this maximum. The velocity ratio $u/c$ is obtained and compared to the critical value $u/c = 0.85$ given by Barthelemy *et al.* (2018). Figure 4 shows the time variation of $u, c$ and $u/c$ for $\chi = 1.92$. For this non breaking wave, the $u/c$ ratio becomes greater than the threshold 0.85 over a significant duration. The maximum $u/c$ ratio occurs while a minimum phase velocity $c$ is being reached. It is also worth noting in the time sequence shown in figure 4 that the crest slows down from $c = 2\,m/s$ to the minimum $c = 1.5\,m/s$. Just after the crest speeds up over a short duration to reach $c = 3\,m/s$. The characteristics of the time variation of the phase velocity (see figure 4) are similar to the observations made in Banner *et al.* (2014) or Fedele *et al.* (2020).

Simulations are performed by varying $\chi$ in the interval $[1 : 1.95]$. Above $\chi = 1.9205$, breaking occurs. We record the maximum of $u/c$ in the time interval when the horizontal velocity reaches its maximum for each factor $\chi$. We denote $(u/c)_m$ this maximum value. We also record the values of $u/c$ when the horizontal velocity $u$ is maximum. The corresponding ratio is denoted $(u/c)_{u_{max}}$ while $u_{max}$ is the maximum horizontal velocity. Figure 5 shows the variation of $(u/c)_m$ and $(u/c)_{u_{max}}$ in terms of $u_{max}$ for non breaking waves only. It appears that a threshold $u/c = 0.92$ is better suited to the present configuration than the standard $u/c = 0.85$. However, we better agree with the standard value if $u/c$ is determined at the instant when the horizontal velocity is maximum.

As mentioned in section 1, other criteria can be proposed to identify the separation between breaking and non-breaking configurations. For that we analyze the pressure, its gradient and its Hessian matrix. There exists a line which yields all these quantities of interest. It can be shown that this line connects the tip of a crest to the sea bottom, even if the crest is overturning like in a plunging breaker, as shown in Scolan (2023). This line is described by the colinearity of the pressure gradient and the eigenvector $\vec{v}_1$ (see eq. A.10), therefore this is the line where $v$ in equation (A.11) vanishes. This continuous curve is described over time $t$ by means of an arclength $\zeta$ starting from the sea bottom and ending almost at the tip of the crest; it is denoted $v(\zeta, t) = 0$. Figure 6 shows all the computed lines at each instant of the simulations for two values of $\chi \in \{1.92, 1.95\}$; those correspond to non-breaking and breaking configurations respectively. It is worth noting that the lines $v = 0$ computed at each time step ($\Delta t = 0.01\,s$ here) are almost superposed whatever the factor $\chi$. Discrepancies become noticeable when breaking occurs. In figure 6 a small difference of crest elevation occurs slightly above the horizontal position $x \approx 11\,m$.

The locations of the main crest are also tracked over time with marks in figure 6. They are located at the top end of these lines $v = 0$. That clearly shows that the lines $v = 0$ indeed arrive at the crest. They are not always vertical lines as illustrated in figure 6. In the present computations, the highest wave is reached at $x \approx 10.5\,m$ along the tank whatever the factor $\chi$.

We now consider the temporal-spatial variations of the wave by analyzing its dynamic properties. The pressure is analyzed along the line $v = 0$ after the maximum horizontal velocity is reached; here we select the instant $t = 2.38\,s$ for a near breaking case associated to the factor $\chi = 1.92$. The location of the corresponding line $v = 0$ is emphasized in figure
6 by a black line; this line arrives almost vertically at the free surface. The pressure follows from Bernoulli equation (A.3). The inertia component ($\phi_{,t}$) and kinetic energy component $\frac{1}{2}(\nabla\phi)^2$ are opposite in sign. In the present case the latter component contributes to the total pressure much less than the former component. But the pressure is mainly hydrostatic so that the pressure monotonically decreases from the sea bottom towards the free surface.

Figure 7 shows the spatial variation of the pressure gradient (modulus) and the La-grangian acceleration (modulus). These quantities are made nondimensional with $\rho g$ and $g$ respectively. It is worth reminding that the vectorial sum of these last two contributions is exactly the gravity. The limiting value 1/2 is almost reached at the crest. This is consistent with the theoretical results by Longuet-Higgins (1963) even for the present near breaking case. However the spatial variations of these quantities are not monotonical when arriving at the free surface. Indeed there is a local extremum that is close to 0.5. For all the studied cases here these spatial variations are identical to the one shown in figure 7. That's why, in the sequel, we record the value of this extremum at some distance of the free surface to be compared to 0.5. Figures 8 show the contours of $p_{,x}/(\rho g)$ and $p_{,y}/(\rho g)$ at the tip of the crest. The horizontal pressure gradient changes sign in a close vicinity of the line $v = 0$. It is confirmed that the extremum close to 0.5 is indeed local and localized very close to the tip. It is also observed that, along the line $v = 0$, the horizontal component $p_{,x}$ is negligible compared to the vertical one $p_{,y}$; the ratio $|p_{,y}/p_{,x}|$ is about 10 in the present simulation. The critical state of the crest can be assessed by examining the Gaussian curvature of the pressure and alternatively the eigenvalues of the pressure Hessian matrix. The spatial variation of these quantities are plotted in figures 9. The variations of $\lambda_1$ and $\lambda_2$ along the line $v$ are similar for two scale factors $\chi = 1.92$ at time $t = 2.38s$ and $\chi = 1.95$ at time $t = 2.29s$. The inverse of the two eigenvalues ($\lambda_1, \lambda_2$) are also the curvature radii along the principal directions given by the eigenvectors, respectively $(\vec{v}_1, \vec{v}_2)$. One of them, $\lambda_2$ is always negative and the other one, $\lambda_1$ can change sign. Along $\vec{v}_2$ the pressure is hence concave. This means that the line $v = 0$ can be viewed as a mountain crest line.

Along most of the line $v = 0$ from the sea bottom, the pressure being hydrostatic, the eigenvalues are almost zero. This also means that the dynamic of the fluid is very small, in other words the fluid is almost at rest. The increase of the eigenvalue $\lambda_1$ up to a positive maximum means that the pressure decreases towards the crest with a positive curvature; in other words the function $p_{v=0}(\zeta)$ is convex. The greater the positive maximum, the flatter the pressure before the inflection point and then before vanishing at the free surface. The change of sign of the eigenvalue $\lambda_1$ in the vicinity of the crest means that the pressure has a local inflection point and the pressure vanishes at the free surface with a non zero normal gradient. The spatial variations of the pressure and its derivatives (first and second) described above are quite standard in a near breaking crest. In particular the threshold $|p_{,y}|/(\rho g) \approx 0.5$ can be retained as a dynamic criterion to define the onset of wave breaking. This is illustrated in figure 10 which shows the time variation of the ratio $u/c$ and the minimum of $|p_{,y}|/(\rho g)$ along the line $v = 0$ for all factors $\chi$.

It should be noted that the duration over which the minimum of $|p_{,y}|/(\rho g)$ remains below 0.5, is similar to the duration over which the ratio $u/c$ is above the threshold 0.85. In addition, a large positive eigenvalue $\lambda_1$ before its change of sign can also be considered as a criterion. This last statement is maybe less intuitive. Figure 11 shows the time variation of the maximum positive $\lambda_1$ along the line $v = 0$ for all factors $\chi$. As we approach the breaking limit, this maximum increases significantly. This suggests to analyze more deeply the time variation of the line $v = 0$ for an overturning crest at a larger scale, for example a plunging breaker in shallow water.
### 3. Dambreaking flow in shallow water

We consider here a simple dam breaking case yielding a nonlinear wave that possibly breaks. The initial free surface position is a part of a Gaussian centered at the right vertical wall

$$y = h + ae^{-2(x-L)^2}, \quad \text{at initial time} \qquad (1)$$

This equation is defined in a coordinate system centered at the left bottom corner of a rectangular tank. The velocity potential on this initial free surface is set to zero. The resulting wave thus propagates in the decreasing $x$ axis direction. The water level $h$ is of $0.5m$ and the length $L$ is chosen depending on the needed time to yield a breaking wave. $h$ is not exactly the mean water level at rest of the total fluid mass but it is the water depth over which the soliton propagates. A parametrical analysis is achieved in terms of the height $a$ setting the initial potential energy. The selected interval of variation of $a$ provides an easy way to cover a range of breaking/non-breaking waves.

By increasing the initial amplitude $a$, we arrive at the breaking threshold. The following example shows how the wave breaks when we choose an initial amplitude $a$ large enough. Here $a = 0.72m$ is the lower limit of the breaking wave configuration; slightly below, say $a = 0.71971m$, the wave does not break. A convergence analysis is achieved in terms of number of markers and time steps. Basically the number of markers $N = 400$ (in tank length $L = 10m$) and a time step $\Delta t = 0.01s$ are chosen. It is shown that the relative errors on mass and energy conservations are bounded by $10^{-4.6}$ and $10^{-3.5}$ respectively. Increasing the number of markers to $N = 1000$ decreases the errors with a factor 10. It is also checked that increasing the discretization $N = 1000$ and $\Delta t = 0.0001s$ does not modify the transition between breaking and non breaking limit. The free surface profiles for a breaking case are plotted for $a = 0.72m$ in figure 12. A global view and a closer view of the overturning wave show that the overturning wave occurs over a tiny area. On the basis of this analysis, $N = 400$ and $\Delta t = 0.01s$ have been selected as it leads to a sufficiently well captured transition.

The relative importance of $u$ and $c$ are shown in figure 13 for the near breaking case with $a = 0.71971m$. This value of $a$ corresponds to the highest initial free surface deformation that leads to a non breaking case. The variations of $u$, $c$ and $u/c$ are similar to those described in figure 4 and which concerns the focused wave. For the present parametrical analysis, the maximum horizontal velocity is reached after the ratio $u/c$ becomes higher than 0.85. It should be noted that the rapid variations of $u$ and $c$ after $t \approx 3s$ are due to the beginning of the reflection on the right wall. This has no interest in the present analysis. Increasing the length of the tank would avoid these disturbances.

Figure 14 shows the variation of the maximum ratio $u/c$ for $a \in [0.7m : 0.71971m]$. In this range of $a$, breaking does not occur. The ratio $u/c$ is plotted in terms of the maximum horizontal velocity. The variation of $(u/c)_m$ and $(u/c)_{u_{max}}$ with $u_{max}$ are identical to the those shown in figure 5 of section §2. Here again it is shown that the threshold 0.85 is not universal.

The analysis of the dynamical properties of the flow along the line $v = 0$ shows exactly the same general variations of the pressure and its derivatives as described in figures 7 to 9. Figure 15 shows the time variation of the minimum of $|p_y|/(\rho g)$ along the line $v = 0$ for all the amplitudes $a$ whatever the breaking or non breaking cases. This quantity drops below 0.5. This result hence proves than the threshold $p_y/(\rho g) = -0.5$ is not universal. However, the fact that the eigenvalue $\lambda_1$ increases substantially as we reach a near breaking configuration remains a clear feature of the limit between breaking and non breaking cases as illustrated in figure 15.
## 4. Conclusion

On the one hand, onset of wave breaking can be defined in terms of kinematic and dynamic criteria. The kinematic criterion $u/c$ is the ratio of velocity at the crest ($u$) to the phase speed ($c$). A threshold is proposed in the literature $u/c = 0.85$. The present analysis shows that this threshold is not universal. In particular, for a focused wave in intermediate depth and for a dambreaking flow leading to a soliton, the ratio $u/c$ would range from 0.9 to 1.

On the other hand, the pioneering works by Stokes (1847) lead to a dynamic criterion formulated in terms of the vertical gradient of the pressure. It is shown that $p_{,y}/(\rho g) = -1/2$ is the critical value for a regular wave. Here we revisit this criterion by examining the properties of the flow along a line that links the seabottom to the crest of a propagating wave packet. This continuous curve is the location where the pressure gradient is colinear to one of the eigenvectors of the Hessian matrix of the pressure. Indeed the critical value $p_{,y}/(\rho g) = -1/2$ can be considered as a threshold, but not universal threshold. Future works should confirm this finding for other configurations.

The two application cases studied here are not very severe because velocities and accelerations remain of a reasonable magnitude. This is not the case in Scolan (2023) where the same analysis is performed by identifying the line $v = 0$. Much higher accelerations (than here) are reached in critical jets suddenly appearing at the free surface. Indeed this analysis reveals premonitory signs of a localized, highly nonlinear free surface flow associated with large velocities or accelerations. That draws the perspectives of future works. One of them is the analysis of the time variation of the line $v = 0$ when a crest overturns at a large scale like in a plunging breaker. It is expected to better describe the transition from the threshold $p_{,y}/(\rho g) = -1/2$ to the free fall of the crest where the fluid acceleration is the acceleration of gravity and hence the pressure gradient must vanish.

## References

* BANNER, M. L., BARTHELEMY, X., FEDELE, F., ALLIS, M., BENETAZZO, A., DIAS, F. AND PEIRSON, W. L. 2014, Linking reduced breaking crest speeds to unsteady nonlinear water wave group behavior, *PhysRevLett.*, .112.114502.

* BARTHELEMY, X., BANNER, M., PEIRSON, W., FEDELE, F., ALLIS, M., & DIAS, F., 2018, On a unified breaking onset threshold for gravity waves in deep and intermediate depth water. *J. Fluid Mech.*, **841**, 463-488.

* CAO, Y., SCHULTZ, W.W., BECK, R.F., 1991, A three-dimensional desingularized boundary integral method for potential problems, *International Journal Numerical Method in Fluids*, 11, 785-803.

* DAVIES, T. V. 1951, The theory of symmetrical gravity waves of finite amplitude. I. *Proc. Roy. Soc A*, 208, 475-86.

* DERAKHTI, M., KIRBY, J.T., BANNER, M.L., GRILLI, S.T., THOMSON, J., 2020, A unified breaking onset criterion for surface gravity water waves in arbitrary depth. *Journal of Geophysical Research: Oceans* e2019JC015886.

* FEDELE, F., BANNER, M., & BARTHELEMY, X., 2020, Crest speeds of unsteady surface water waves. *J. Fluid Mech.*, **899**, A5. doi:10.1017/jfm.2020.424
HAVELOCK, T. H. 1918, Periodic irrotational waves of finite height. *Proc. Roy. Soc. A*, **95**, 38-51.

KIMMOUN O., AND SCOLAN Y.-M., 2009, Generation of focalized wave packet *Proc. 24th International Workshop on Water Waves and Floating Bodies*, Zelenogorsk, Russia.

LONGUET-HIGGINS, M., 1963. The generation of capillary waves by steep gravity waves. *Journal of Fluid Mechanics*, **16**(1), 138-159.

LONGUET-HIGGINS, M.S. 1980, A technique for time-Dependent free-surface flows. *Proceedings of the Royal Society of London. Series A, Mathematical and Physical Sciences*, **371**, 441-451.

NADAOKA, K., ONO, O. AND KURIHARA, H., 1997, Analysis of near-crest pressure gradient of irregular water waves as a dynamic criterion of breaking, *Proc. The Seventh International Offshore and Polar Engineering Conference*, Honolulu, Hawaii, USA, May 1997.

SCOLAN, Y.M., 2010, Some aspects of the flip-through phenomenon: A numerical study based on the desingularized technique *J. Fluid Struc.*, **26**, Issue 6, 918-953.

SCOLAN, Y.M., 2023, Some aspects of the pressure field preceding the onset of critical jets in a breaking wave. *J. Eng. Math.*, **138**, https://doi.org/10.1007/s10665-022-10246-3.

STOKES G. G., 1847, On the theory of oscillatory waves. *Trans. Camb. Phil. Soc.*, **8**, 441-55.

TUCK E.O., 1998, Solution of nonlinear free-surface problems by boundary and desingularised integral equation techniques, Invited Lecture, *Computational Techniques and Applications: CTAC'97 Eds J. Noye, M. Teubner, A. Gill, World Scientific, Singapore*, 11-26.

VARING A., FILIPOT J.F., GRILLI S., DUARTE R. AND ROEBER V., 2021 A new definition of the kinematic breaking onset criterion validated with solitary and quasi-regular waves in shallow water. *Coastal Engineering*, **164**, 103755.
<figure><img src="image_2.png" /><figcaption>Figure 1: Free surface elevation (green and blue) and distribution of the velocity potential on the free surface (red) along the horizontal position from left wall to right wall of the tank. Length of the rectangular tank $L = 15m$. Mean water depth $h = 0.7m$. Blue line : free surface elevation of reference corresponding to $\chi = 1$. Green line : free surface elevation corresponding to $\chi = 1.9205$.</figcaption></figure>
<figure><img src="image_3.png" /><figcaption>Figure 2: Successive free surface profiles up to $t = 3.5s$ for the scale factor $\chi = 1.9205$. Time lag between each profile : 0.1s. The wave train propagates along the increasing $x$ axis. Initial profile in black.</figcaption></figure>
<figure><img src="image_4.png" /><figcaption>Figure 3: Successive free surface profiles for a breaking case with scale factor $\chi = 1.95$. The inserted figure is a closer view in the area drawn by a black rectangle. Only the free surface profiles computed in the time interval $t \in [1.5s, 2.29s]$ are plotted. In the main part of the figure, just 1 in every 5 of the times computed is drawn as a profile. In the inserted graph all profiles are plotted with a time step : 0.01s.</figcaption></figure>
<figure><img src="image_5.png" /><figcaption>Figure 4: Time variation of $u$, $c$ and $u/c$ for scale factor $\chi = 1.92$ corresponding to a near breaking wave. Threshold $u/c = 0.85$: black line.</figcaption></figure>

<figure><img src="image_6.png" /><figcaption>Figure 5: Ratio $(u/c)_{u_{max}}$ calculated at the maximum horizontal velocity $u_{max}$. Maximum ratio $(u/c)_m$ reached over the simulation. Variation of $(u/c)_m$ and $(u/c)_{u_{max}}$ in terms of $u_{max}$ for non breaking waves $\chi < 1.9205$. Threshold $u/c = 0.85$: black line.</figcaption></figure>
<figure><img src="image_7.png" /><figcaption>Figure 6: Plot of the lines $v(\zeta, t) = 0$ for two scale factors : $\chi = 1.92$ (nonbreaking case in red) and $\chi = 1.95$ (breaking case in green). The black line corresponds to $v(\zeta, t) = 0$ at $t = 2.38s$ for $\chi = 1.92$. Superposed marks indicate the location of the crest. The wave propagates along the increasing $x$ axis. Time lag between each line : 0.01s. The coordinate system is centered at the left bottom corner of a rectangular tank.</figcaption></figure>
<figure><img src="image_8.png" /><figcaption>Figure 7: Spatial variation along the curve $v = 0$ of the pressure gradient (in green modulus non-dimensionalized with $\rho g$) and the Lagrangian acceleration (in blue modulus non-dimensionalized with $g$). The origin of the arclength is at the sea bottom. Scale factor $\chi = 1.92$ at time $t = 2.38s$. The corresponding line $v = 0$ is drawn with a black line in figure (6).</figcaption></figure>
(a)

<figure><img src="image_9.png" /></figure>

(b)

<figure><img src="image_10.png" /></figure>

Figure 8: Contours of the horizontal (a) and vertical (b) components of the pressure gradient in the tip of the crest for the scale factor $\chi = 1.92$ at time $t = 2.38s$. The pressure gradient components are non-dimensionalized with $\rho g$.
<figure><img src="image_11.png" /><figcaption>Figure 9: Spatial variation of the two eigenvalues $\lambda_1$ and $\lambda_2$ along the line $v = 0$. The inner graph is a closer view of the region where the largest variations occur. Thick line : scale factor $\chi = 1.92$ at time $t = 2.38s$. Light line : scale factor $\chi = 1.95$ at time $t = 2.29s$ corresponding to the last computed profile plotted in figure (3). Non-dimensionalized with $\rho g/h$.</figcaption></figure>
<figure><img src="image_12.png" /><figcaption>Figure 10: Time variation of the ratio $u/c$ (a) and the minimum $|p_y|/(\rho g)$ (b) along the line $v = 0$ for all scale factors $\chi$ ranging from 1 to 1.95. Blue curve: greatest factor $\chi$ without breaking $\chi = 1.9205$. Arrow: increasing factor $\chi$.</figcaption></figure>
<figure><img src="image_13.png" /><figcaption>Figure 11: Time variation of the maximum positive $\lambda_1 h/(\rho g)$ along the line $v = 0$ for all factors $\chi$. Arrow : increasing factor $\chi$.</figcaption></figure>
<figure><img src="image_14.png" /><figcaption>Figure 12: Successive free surface profiles for a breaking case $a = 0.72m$. The inserted figure is a closer view in the area drawn by a rectangle. Initial free surface deformation given by equation (1) with $L = 10m$, $h = 0.5m$. Number of markers $N = 1000$ and time step : $\Delta t = 0.01s$ for $t \le 2s$, $\Delta t = 0.001s$ for $t \le 2.1s$, $\Delta t = 0.0001s$ for $t > 2.1s$. In the main graph the profiles are not all plotted. In the inserted graph 1 profile over 20 is plotted.</figcaption></figure>
<figure><img src="image_15.png" /><figcaption>Figure 13: Time variation of $u$, $c$ and $u/c$. Initial free surface deformation given by equation (1) with $a = 0.71971\text{m}$, $L = 10\text{m}$, $h = 0.5\text{m}$. Threshold $u/c = 0.85$: black line.</figcaption></figure>

<figure><img src="image_16.png" /><figcaption>Figure 14: Ratio $(u/c)_{u_{max}}$ calculated at maximum horizontal velocity $u_{max}$. Maximum ratio $(u/c)_m$ reached over the simulation. Variation of $(u/c)_m$ and $(u/c)_{u_{max}}$ in terms of $u_{max}$ for non breaking wave. Threshold $u/c = 0.85$: black line.</figcaption></figure>
<figure><img src="image_17.png" /><figcaption>Figure 15: Time variation of the minimum of $|p_y|/(\rho g)$ (a) maximum positive $\lambda_1 h/(\rho g)$ (b) along the line $v=0$ for all amplitude $a$. Initial free surface deformation given by equation (1) with $a \in [0.7m : 0.73m]$, $L = 10m$, $h = 0.5m$.</figcaption></figure>
## Appendix A. Numerical modelling

The fully nonlinear free surface boundary conditions are solved in potential theory. The velocity potential is formulated in the frame of a desingularized technique (see Tuck, 1998). The fluid is contained in a fixed tank filled over a mean depth $h$. Here the tank is rectangular with length $L$ and the equations are posed in a twodimensional space. The variables are calculated in a coordinate system centered at the bottom left corner of the tank. A conformal mapping turns the inner tank (infinite half vertical strip) described in the physical plane $z$ into a half space in the transformed plane $w$. As a consequence the impermeability conditions are easily accounted for. The complex potential $F$ is expanded at time $t$ as a finite sum of Rankine source singularities

$$F(z,t) = \sum_{j=1}^{N} q_j(t)G(w, W_j(t)) \qquad (A.1)$$

The sources are located at the complex coordinates $Z_j$ in the physical plane and $W_j$ in the transformed plane; they are mutual images by using the conformal mapping function $\ell$ such that $W_j = \ell(Z_j)$. The intensity of the source $\# j$ is $q_j$. $G$ is the complex potential of this source plus its mirror image with respect to the axis $\Im(w) = 0$,

$$G(w, W_j) = \log(w - W_j) + \log(w - \overline{W_j}) \qquad (A.2)$$

The overline denotes the conjugate of the complex variable. The singularities are located at a short distance from the actual free surface as proposed in Cao et al. (1991). A Lagrangian marker belonging to the free surface is attached to each singularity. The velocity potential $\phi = \Re(F)$ is calculated at these Lagrangian markers. The time differential system is formulated for the velocity potential $\phi$ and the Cartesian coordinates of the Lagrangian markers. A fourth order Runge-Kutta algorithm is used to update these quantities over time (see Scolan, 2010).

For the applications considered in section 2 the accuracy of the computations is checked over time in terms of conservation of mass and energy. This is illustrated in figure (A.16) which shows the time variations of the relative errors on the mass and energy conservations corresponding to the simulations described in section §2. Here the mass conservation is

<figure><img src="image_18.png" /><figcaption>Figure A.16: Time variation of the relative errors on mass (left) and energy (right) conservations for all scale factors. Number of markers $N = 1000$. Time step $\Delta t = 0.01s$.</figcaption></figure>

ensured with a maximum relative error of $10^{-5.5}$ and the energy conservation is met with a maximum relative error of $10^{-3.8}$.
Pressure is calculated a posteriori from Bernoulli equation

$$p = -\rho \left( \phi_{,t} + \frac{1}{2} \left( \vec{\nabla} \phi \right)^2 + g(y-h) \right) \quad (A.3)$$

Pressure is always positive in the fluid and cannot be expected to drop below the ambient pressure. Without loss of generality this ambient pressure is set arbitrarily to zero.

Pressure is zero at the free surface. The pressure gradient $\vec{\nabla}p$ is always oriented in the direction normal to the free surface, pointing inside the fluid. The Euler equation relates the pressure gradient $\vec{\nabla}p$ to the velocity $\vec{u}$. The Euler equation is formulated here in its complex form

$$F_{,zt} + F_{,zz} \overline{F_{,z}} + \frac{p_{,z}}{\rho} - ig = 0 \quad (A.4)$$

where $p_{,z} = p_{,x} - ip_{,y}$ and the sum of the first two terms is the complex Lagrangian acceleration. The spatial and temporal derivatives of $F$ are obtained in closed form starting from equation (A.1).

In order to better characterize the spatial variation of the pressure, the Hessian matrix of the pressure is introduced. It is denoted $\mathbf{H}(x, y, t)$ and its expression is

$$\mathbf{H} = \begin{pmatrix} p_{,xx} & p_{,xy} \\ p_{,xy} & p_{,yy} \end{pmatrix} \quad (A.5)$$

This is a symmetric matrix and its elements are (see Longuet-Higgins, 1980a)

$$p_{,xx} = -\rho (|F_{,zz}|^2 + \Re(F_{,zzz} \overline{F_{,z}} + F_{,zzt})) \quad (A.6)$$

$$p_{,yy} = -\rho (|F_{,zz}|^2 - \Re(F_{,zzz} \overline{F_{,z}} + F_{,zzt})) \quad (A.7)$$

$$p_{,xy} = \rho \Im(F_{,zzz} \overline{F_{,z}} + F_{,zzt}) \quad (A.8)$$

By using the notations

$$f = \rho |F_{,zz}|^2, \quad A = \rho (F_{,zzz} \overline{F_{,z}} + F_{,zzt}), \quad (A.9)$$

the eigenvalues $\lambda$ and the eigenvectors $\vec{v}$ of the matrix $\mathbf{H}(x, y, t)$ are expressed as follows

$$\lambda = -f \pm |A|, \quad \vec{v} = \begin{pmatrix} \Im(A) \\ \Re(A) \pm |A| \end{pmatrix}. \quad (A.10)$$

The eigenvectors are normalized by dividing their components with their norm. The Gaussian curvature $\Omega$ is the determinant of the matrix $\mathbf{H}$. This is the product of the eigenvalues, $\Omega = \lambda_1\lambda_2$. The inverse of each eigenvalue is also the curvature radius of the surface $p(x, y)$ along their corresponding principal directions defined by the eigenvectors. The eigenvalue (denoted $\lambda_2$) is always negative meaning that the surface $p(x, y)$ is concave along the direction of $\vec{v}_2$. The other eigenvalue (denoted $\lambda_1$) can change sign when $f = |A|$. This also corresponds to the change of sign of the Gaussian curvature $\Omega = |A|^2 - f^2$.

In the fluid region we seek the points where $\vec{\nabla}p$ and $\vec{v}_1$ are colinear. This occurs when the following quantity changes sign

$$v = \Im(A)p_{,y} - (\Re(A) + |A|)p_{,x} \quad (A.11)$$