# Irrotationality of Water Waves and Topography

Alan Riquier, Emmanuel Dormy

► To cite this version:

Alan Riquier, Emmanuel Dormy. Irrotationality of Water Waves and Topography. 2024. hal-04782997v1

HAL Id: hal-04782997

https://hal.science/hal-04782997v1

Preprint submitted on 14 Nov 2024 (v1), last revised 20 Nov 2025 (v2)

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 Attribution 4.0 International License
# Irrotationality of Water Waves and Topography

Alan Riquier & Emmanuel Dormy

## Abstract

Most theoretical and numerical results regarding water waves rely on the irrotationality of the flow in order to simplify significantly the mathematical formulation. This assumption is rarely discussed in details. In this work we investigate the well-foundedness of this important hypothesis using numerical simulations of the free-surface Navier-Stokes equation, using a scheme introduced in Riquier and Dormy (2024). We show that, in the presence of an irregular bottom, a gravity wave of non-negligible height can effectively destabilise the bottom boundary layer, stripping off vortices into the main flow. As a vortex approaches the surface the solution of the Navier-Stokes flow in the limit of vanishing viscosity is shown to differ from the inviscid Euler solution.

## 1 Introduction

The mathematical description of water waves remains one of the oldest yet most challenging topics in fluid mechanics, even to this day. The highly nonlinear nature of the phenomenon, coupled with the very high values of the Reynolds number of the flow, strongly restrain the use of the full equations in both analytical and numerical works. Hence, through the ages, many simplified models have been proposed in order to reduce these difficulties. Most of these models correspond to asymptotic regimes of the Water Waves equations (Lannes 2013), which are themselves derived from Euler equations assuming the flow is irrotational (Craig and Sulem 1993; Zakharov 1968). This assumption, while being widely accepted, is rarely justified. We investigate here its validity.

In the last twenty-five years, a wide range of numerical studies were able to solve for the full non-linear incompressible free-surface Navier-Stokes equations and observed, in light of the analytical result of Lundgren and Koumoutsakos (1999), the appearance of a boundary layer below the free-surface (Raval et al. (2009); it is also visible in Deike et al. (2015)). Recently (Riquier and Dormy 2024), we showed that the effects of this vorticity sheet are only noticeable on long time-scale or on small length-scale (i.e. when the curvature becomes comparable to the Reynolds number). We will show here that the situation significantly differs in the case of a bottom topography with no-slip boundary conditions.

In the absence of spilling or breaking, the free surface contributes little to the vorticity generation. The present work is thus concerned with the influence of the bottom boundary of the fluid domain on the irrotationality. We show, through a numerical study of the incompressible free-surface Navier-Stokes equations, that a wave passing over an obstacle can strip vorticity from the boundary layer, leading to the formation of vortices in the fluid. The effect of these vortices on the free surface do not decrease with the viscous dissipation. It is shown that this mechanism does not only arise around sharp obstacles but also next to smooth boundaries with a sufficient curvature. The presence of these vortices has a non-negligible effect on the surface evolution, which seriously questions the use of irrotational models for waves above sharp topographies as well as smooth, i.e. infinitely differentiable, topographies with regions of large curvature where the boundary layer can separate.

## 2 Mathematical description and numerical framework

We solve the full non-linear Navier-Stokes equation with a free surface (see Riquier and Dormy 2024), and consider a non-flat bottom topography on which we apply the no-slip boundary condition. The fluid domain $\Omega(t)$ is a two-dimensional domain, periodic in $x$ with period $L$, bounded by a rigid bottom $\Gamma_b$ and a moving interface $\Gamma_i(t)$. We introduce $\gamma_b(s)$ a parametrisation of $\Gamma_b$ and $\gamma_i(t, s)$ a time-dependant parametrisation of $\Gamma_i(t)$. In order to carry out smoothing of the topography in section 5, we further assume that $\gamma_b(s)$ corresponds to the arc-length parametrisation. A schematic representation of this configuration is presented in figure 1.
<figure><img src="image_2.png" /><figcaption>Figure 1 – Geometry of the domain $\Omega(t)$ in dimensional units. The bottom topography can be chosen arbitrarily.</figcaption></figure>

The domain $\Omega(t)$ is filled with a fluid (usually water) whose kinematic viscosity is denoted by $\nu$. We neglect the effects of an upper fluid (usually air) on the dynamics of the flow (in effect setting its density to zero). We also disregard the effects of surface tension, even though they could easily be included in this formulation. The characteristic length of the flow is chosen as the distance $h_0$ between the average height of the interface and the minimum height of the bottom (see figure 1). The chosen typical velocity $U$ corresponds to the phase velocity of an inviscid gravity wave in shallow water, i.e. $U = \sqrt{gh_0}$. The non-dimensional governing equations then take the form

$$ \frac{\partial \mathbf{u}}{\partial t} + (\mathbf{u} \cdot \nabla)\mathbf{u} = -\nabla p + \frac{1}{Re} \Delta \mathbf{u} - \hat{\mathbf{y}} \quad (1a) $$

$$ \nabla \cdot \mathbf{u} = 0, \quad (1b) $$

where the Reynolds number defined as

$$ Re = \frac{h_0 \sqrt{gh_0}}{\nu}, \quad (2) $$

and where $\hat{\mathbf{y}}$ is the unit vector along the y-axis (figure 1). It appears through the non-dimensionalisation of the gravity term in the dimensional Navier-Stokes equation.

The no-slip condition $\mathbf{u} = 0$ is imposed on the bottom boundary $\Gamma_b$. This is in contrast with Riquier and Dormy (2024) where impenetrability in the normal direction and vanishing stress in the tangential direction were chosen instead. As we shall see, the formation of a boundary layer near $\Gamma_b$ can have a significant effect on the overall flow and on the interface.

On the interface $\Gamma_i(t)$ the stress-free condition,

$$ p\hat{\mathbf{n}} - \frac{1}{Re} \left[ \nabla \mathbf{u} + (\nabla \mathbf{u})^t \right] = 0, \quad (3) $$

is enforced. Each fluid particle on the interface is advected in a lagrangian manner. This is achieved through the following kinematic condition

$$ \frac{\partial \gamma}{\partial t}(t, s) = \mathbf{u}(t, \gamma(t, s)). \quad (4) $$

We are interested in comparing our results with those obtained solving the irrotational Euler equation. To that end, we rely on the code of Dormy and Lacave (2024). Hence an irrotational initial condition is chosen. The initial interface $\gamma(0)$ is chosen as a simple cosine wave of finite amplitude $A$ and the initial velocity $\mathbf{u}_0$ corresponds to the first order solution of the water waves equations (see e.g. Johnson (1997), chap. 2),

$$ \gamma(0, x) = \begin{bmatrix} x \\ h_0 + A \cos(kx) \end{bmatrix} \quad \text{and} \quad \mathbf{u}_0 = A \sqrt{gk \tanh(kh_0)} \cdot \begin{bmatrix} \left(\tanh(kh_0)\right)^{-1} \cos(kx) \\ \sin(kx) \end{bmatrix}, \quad (5) $$

where $h_0 = 1$. The solution to this problem is computed numerically using the FreeFem finite-element library on unstructured meshes (Hecht 2012). A crank-Nicholson time-stepping scheme is used for all terms but the non-linearity, which is treated in a mixed explicit-implicit manner in order to prevent the use of a nonlinear solver. In order to advect the interface points with the velocity of the fluid, we employ
the arbitrary lagrangian-eulerian method of Hirt et al. (1974). The entire mesh is advected at each time step with a velocity $w$, obtained after numerically solving the following problem

$$ \left\{ \begin{array}{ll} \Delta w = 0 & \text{in the fluid } \Omega(t) \\ w = u & \text{on the interface } \Gamma_i(t) \\ w = 0 & \text{on the bottom } \Gamma_b \end{array} \right. \qquad (6) $$

This is a lagrangian advection scheme for the interface points while preventing spurious triangle elongation near the bottom boundary layer. The code is made multi-threaded through the use of an additive Schwarz solver combined with a geometric multigrid preconditioner in order to maintain the convergence speed as the number of MPI cores is increased.

## 3 Rectangular step and non-convergence of the viscous solution

<figure><img src="image_3.png" /><figcaption>Figure 2 – Evolution of the interfaces for different values of Re and the inviscid (Euler) solution. A video animation is available as a supplementary material.</figcaption></figure>
Simulations have been carried out with values of the Reynolds number ranging from $Re = 10^2$ to $Re = 10^5$. The length of the domain is set to $L = 2\pi$ and the amplitude of the initial cosine wave is $A = 0.1$. A larger value of $A$ would lead to a breaking wave in short time. A rectangular step of length $\ell_b = 2\pi/3$ and height $h_b = 0.5$ characterises the bottom $\Gamma_b$, otherwise flat, see figure 1.

We perform comparisons with the inviscid and irrotational solution of the Euler equation, computed using the dipole method of Dormy and Lacave (2024). The interfaces for different values of $Re$ and the irrotational solution are shown in figure 2. In order to emphasize the logarithmic nature of the convergence, the Reynolds number has been increased in multiplicative steps of size $\sqrt{10}$.

The solution pictured in figure 2 exhibits an interesting behaviour. A difference between the viscous and inviscid solution is observed, which does not disappear as the viscosity is decreased. This difference stems from the stripping of vortices from the bottom boundary layer, as shall be emphasised in section 4.

In order to understand whether the difference between the Navier-Stokes solution and the Euler solution might disappear should we increase the Reynolds number any further, we need to investigate the convergence in a quantitative manner. To do so, since the interface remains a graph, i.e. no overturning occurs, we have computed the $L^2$ distance between the two interfaces for different values of the Reynolds. Let $\eta_E(x, t)$ be the Euler interface represented as a graph, and $\eta_{Re}(x, t)$ the viscous solution with Reynolds $Re$. Then the $L^2$ distance between both curves is defined as

$$ \delta_{Re}(t) = \left[ \int_{0}^{L} (\eta_{Re}(x, t) - \eta_{E}(x, t))^2 dx \right]^{1/2} . \qquad (7) $$

<figure><img src="image_4.png" /><figcaption>Figure 3 – Evolution of the $L^2$ distance between viscous interfaces and either the inviscid (Euler) interface (left) or the highest Reynolds number $Re = 10^5$ interface (right).</figcaption></figure>

As can be seen in figure 3 (left), the distance between the viscous and inviscid free surfaces converges to a non-vanishing value which slowly increases with time. Figure 3 (right) presents the distance (7) between the $Re = 10^5$ interface and those corresponding to lower values of $Re$. The asymptotically vanishing distance between the viscous solutions as the Reynolds number is increased is a key feature of convergent sequences. Hence both plots lead us to conclude that the viscous solution converges to something that does not correspond to the irrotational inviscid solution, even though the initial condition was irrotational. This is due to the appearance of vortices in the fluid domain, in contrast with the flat bottom case we studied previously (Riquier and Dormy 2024). These do not vanish with the viscosity.

## 4 Vortices emission

The origin of the non-vanishing difference highlighted in the previous section can be understood by computing the two-dimensional vorticity $\omega = \partial_x u_y - \partial_y u_x$ of the flow. This is shown in figure 4 for $Re = 10^5$. A comparison between the solutions for various values of $Re$ at time $t = 10$ is shown in figure 5. Note that since the value of $\omega$ in the vortices increases linearly with the Reynolds number (as the vortices get more and more localised), the color scheme has been truncated in order to highlight regions of strong vorticity.
<figure><img src="image_5.png" /><figcaption>**Figure 4** – Evolution of the vorticity at $Re = 10^5$. The color scheme has been centered on the range [-5, 5] to highlight the vortices. A video animation is available as a supplementary material.</figcaption></figure>

The appearance of vortices near the corners of the rectangular step is observed. The velocity associated with the passage of the wave induces an instability in the bottom boundary layer which, ultimately, leads to a pair of counter-rotating vortices being stripped off. The globally averaged vorticity remains zero throughout the simulation but the presence of vortices far from the boundary, whose effects do not vanish with the viscosity, prevent the solution to converge to the irrotational Euler solution as $Re$ tends to infinity. Instead the solution appears to converge toward a flow characterised by localised vortices.

At $Re = 10^5$, vortices are transported in the vicinity of the free surface before having decayed (fig. 4 at $t = 20$). Should we increase the Reynolds number any further, a vortex pair could eventually reach the interface with a high vorticity. This has however not been observed so far.

Figure 4 highlights another surprising behavior of the boundary layer: the first pair of vortices emitted on the left corner ($t = 4$) translates to the bottom boundary layer ($t = 14$) and strips away two new vortices ($t = 17.5$) and even more afterward ($t = 20$).

In order to understand the vortex-generating mechanism, it is useful to observe the streamlines near a corner as the first pair of vortices is generated. It is presented in figure 6 for $Re = 10^4$. As usual when considering water waves, the streamlines globally wrap around the crests and the troughs of the wave. Locally, in the vicinity of the obstacle the underlying irrotational flow gets through three different stages. First, near $t = 2$ it corresponds to an up-going flow near a corner, whose velocity potential is well-known (see e.g. Batchelor 1967, sec. 6.5). At a later time (about $t = 4$), as the inflection point of the interface passes above the salient edge, the direction of the stream changes. The streamlines no longer wrap around the topography. Instead the horizontal left-coming flow meets the vertical up-going one at the corner in a fashion similar to that of the rear tip of an aircraft's wing. These streamlines effectively guide the pair of vortices out of the boundary layer. Finally, at subsequent times (around $t = 6$ in fig. 6),
<figure><img src="image_6.png" /><figcaption>Figure 5 – The vorticity at $t = 10$ for Reynolds numbers $Re$ ranging from $Re = 10^{3.5}$ to $Re = 10^5$. The color scheme has been centered on the range $[-5, 5]$ to highlight the vortices.</figcaption></figure>

the underlying irrotational flow wraps once again around the corner, but now in the opposite direction.

Throughout the passage of a crest, the direction of the flow outside the boundary layer corresponds to the wave propagation direction and the vorticity generated in the boundary layer has a negative sign (fig. 6 right). As a trough passes above the edge, the stream’s direction is opposite to the wave propagation so the vorticity has a positive sign (fig. 6 left). This yields, as the vortices are separated from the boundary layer, the correspondence crest ↔ negative vortex and trough ↔ positive vortex. However the complex vortex dynamics happening once the vortices have been separated prevents the Navier-Stokes solution to converge to the Euler flow.

# 5 Smooth topography and influence of the curvature

We observed the emission of vortices near sharp corners following from vorticity stripping from the boundary layer. We now want to investigate whether this behavior persists around edges of finite curvature. To that end, we consider smoothed out versions of the rectangular topography considered in the previous section. We use the method of mollifiers, introduced originally by Friedrichs (1944) and widely used in the mathematical study of partial differential equations.

Let $\gamma_b(s)$ the parametrisation of the rectangular step. This function is non-differentiable at the edges. We thus introduce the following bump function,

$$ \phi(x) = \begin{cases} C \exp\left(-\frac{1}{1-x^2}\right) & \text{for } -1 < x < 1, \\ 0 & \text{otherwise.} \end{cases} \qquad (8) $$
<figure><img src="image_7.png" /><figcaption>Figure 6 – The streamlines (white lines with arrows) and vorticity $\omega$ near the top right corner of the obstacle at $t = 2, 4$ and $t = 6$ for $Re = 10^4$. The color scheme has been centered on the range [-5,5] to highlight the vortices.</figcaption></figure>

where $C$ is a normalisation factor ensuring that the integral of $\phi$ from $x = -1$ to $x = 1$ is equal to 1. This function is infinitely differentiable. Let us introduce a parameter $\epsilon > 0$ which will correspond to the curvature of the obstacle. The mollifier $\phi_{\epsilon}$ is then defined as

$$ \phi_{\epsilon}(x) = \frac{1}{\epsilon} \cdot \phi\left(\frac{x}{\epsilon}\right), \qquad (9) $$

and can be used to smooth out the function $\gamma_b(s)$ using a convolution

$$ \gamma_{\epsilon}(s) = \int_{-1}^{1} \phi_{\epsilon}(t) \gamma_{b}(s-t) dt. \qquad (10) $$

The resulting parametrisation $\gamma_{\epsilon}(s)$ is smooth and differentiable. It converges to $\gamma_b(s)$ as $\epsilon \to 0$.

We now use the parametrisation $\gamma_{\epsilon}(s)$ constructed above to define the bottom topography $\Gamma_b$ in three new simulations with $Re = 10^5$ and $\epsilon = 1, 0.5$ and $0.1$. The initial condition corresponds to the one described in section 2.

The associated vorticity is presented in figure 7. At $\epsilon = 1$ the boundary layer is stable. At $\epsilon = 0.5$ a vortex is detached on the right corner only. Finally at $\epsilon = 0.1$ another one is ripped off near the left corner.

The presence of these vortices is reminiscent of the von Kármán (1911) votrex street, in which an adversary pressure gradient in the boundary layer leads to the emission of vortices past a cylinder when the Reynolds number is increased. We are, in fact, considering a rather similar situation except that vorticity stripping is not here associated with a steady stream, but instead results from an oscillation of the main stream in time.

# 6 Conclusion

In the present work we use numerical simulations to investigate the separation of boundary layers and vortex emission associated with water waves. We identify cases for which the Navier-Stokes solution does not converge to the irrotational Euler solution in the limit of vanishing viscosity. We observe that vortices are emitted near regions of negative curvature.

When using a water wave model derived from Euler's equation with irrotationality assumed, as is often the case, one must keep in mind that such models bypasses a major physical phenomenon, which does not disappear as the Reynolds number is increased and whose effects on the shape of the interface
<figure><img src="image_8.png" /><figcaption>**Figure 7** – Emission of vortices at $Re = 10^5$ and $t = 10$ before and after an obstacle as its smoothness $\varepsilon$ varies $\varepsilon = 1$ (up), $\varepsilon = 0.5$ (middle), $\varepsilon = 0.1$ (down). The color scheme has been centered on the range [-2, 2] to highlight the vortices. A video animation is available as a supplementary material.</figcaption></figure>

is not negligible. This is particularly true in the near-shore analysis of the ocean when these vortices are emitted near *e.g.* immersed sharp boulders.

## Acknowledgments

This work has been presented during the 1st European Fluid Dynamics Conference (EFDC1) in Aachen in
september 2024. The authors wish to thank Christophe Lacave, Geoffrey Beck, Erwan Faou and Vincent
Duchêne for useful discussions, as well as Pierre Jolivet, Bertrand Maury and Georges Sadaka on the use
of FreeFem. The simulations have been carried out with help from the MeSU supercomputing center at
Sorbonne Université, Paris.

## References

* Batchelor, G. K. (1967). *An Introduction to Fluid Dynamics*. Cambridge University Press.

* Craig, W. and C. Sulem (1993). Numerical simulation of gravity waves. In: *Journal of Computational Physics* **108**, pp. 73–83.

* Deike, L., S. Popinet, and W. K. Melville (2015). Capillary effects on wave breaking. In: *J. Fluid Mech.* **769**, pp. 541–569.

* Dormy, E. and C. Lacave (2024). Inviscid Water-Waves and interface modeling. In: *Quarterly of Applied Mathematics*, in press, available online. DOI: https://doi.org/10.1090/qam/1685 and arXiv:2306.02363.

* Friedrichs, K. O. (1944). The identity of weak and strong extensions of differential operators. In: *Trans. Amer. Math. Soc.* **55**, pp. 132–151.

* Hecht, F. (2012). New development in FreeFem++. In: *J. Numer. Math.* **20.3-4**, pp. 251–265. URL: https://freefem.org/.

* Hirt, C.W, A.A Amsden, and J.L Cook (1974). An arbitrary Lagrangian-Eulerian computing method for all flow speeds. In: *Journal of Computational Physics* **14.3**, pp. 227–253.

* Johnson, R. S. (1997). *A Modern Introduction to the Mathematical Theory of Water Waves*. Cambridge University Press.
<figure><img src="image_9.png" /><figcaption>**Figure 8** – The streamlines (white lines with arrows) and vorticity $\omega$ near the top right corner of the smoothed obstacle ($\epsilon = 0.1$) at $t = 2, 4$ and $t = 6$ for $Re = 10^5$.</figcaption></figure>

Lannes, D. (2013). *The Water Waves Problem: Mathematical Analysis and Asymptotics*. Vol. 188. *Mathematical Surveys and Monographs*. American Mathematical Society.

Lundgren, T. and P. Koumoutsakos (1999). On the generation of vorticity at a free surface. In: *J. Fluid Mech.* **382**, pp. 351–366.

Raval, A., X. Wen, and M. H. Smith (2009). Numerical simulation of viscous, nonlinear and progressive water waves. In: *Journal of Fluid Mechanics* **637**, pp. 443–473.

Riquier, A. and E. Dormy (2024). Numerical study of a viscous breaking water wave and the limit of vanishing viscosity. In: *Journal of Fluid Mechanics* **984**, R5. DOI: 10.1017/jfm.2024.208.

von Kármán, T. (1911). Ueber den Mechanismus des Widerstandes, den ein bewegter Körper in einer Flüssigkeit erfährt. In: *Nachrichten von der Gesellschaft der Wissenschaften zu Göttingen, Mathematisch-Physikalische Klasse*, pp. 509–517.

Zakharov, V. E. (1968). Stability of periodic waves of finite amplitude on the surface of a deep fluid. In: *Journal of Applied Mechanics and Technical Physics* **9**, pp. 190–194.

**Alan Riquier**

DÉPARTEMENT DE MATHÉMATIQUES ET APPLICATIONS (DMA), UMR 8553
ÉCOLE NORMALE SUPÉRIEURE - PSL
45 rue d'Ulm, 75005 Paris, France
Email adress: alan.riquier@ens.fr

**Emmanuel Dormy**

DÉPARTEMENT DE MATHÉMATIQUES ET APPLICATIONS (DMA), UMR 8553 AND CNRS
ÉCOLE NORMALE SUPÉRIEURE - PSL
45 rue d'Ulm, 75005 Paris, France
Email adress: emmanuel.dormy@ens.fr