Solitary wave propagation and wave breaking: two benchmark problems for low Mach number solvers

"Mathematical and Numerical Aspects of Low Mach Number flows", Porquerolles, France, June 21-25 2003 

Wave computation and animation by F. Golay and P. Helluy.


We propose a workshop dedicated to the modelling and simulation of  the breaking due to gravity of a solitary wave over a step-like reef. It is classical in that simulation to suppose that the whole flow (air plus water) is incompressible. This workshop aims at evaluating the behavior of low Mach number solvers applied to wave breaking. Because this simulation is difficult, the benchmark is splitted into two parts: propagation of a solitary wave, and then, the breaking.

Figure: Geometry of the step-like reef.

The benchmark problem we have selected is described in the paper [1] . That paper provides physical experiments, including PIV velocity measurements  and also numerical results, using an incompressible model, before the breaking.


 The computation has to be performed in both air and water. The two fluids are supposed to be inviscid, but compressibility has to be modelled. The air will be assimilated to a perfect gas with a pressure law (r is the density and e is the internal energy)


The liquid is almost incompressible. For example, it can be supposed to satisfy a stiffened-gas equation of state ([3] ).


Because the flow is incompressible, the participants are allowed to modify the pressure laws in order to obtain shorter computations. The modelling can be adapted to the numerical method used by the participants. A possible modelling is the following.

The flow is subject to the Euler equations (u and v are the two components of the velocity, Y is the mass fraction of air) with a source term due to gravity


with the pressure law


For the breaking simulation, the computation domain can be started at the position of the gauge P1. It has to enclosed a sufficient quantity of air (let us say around one meter over the free surface).


At time t=0, the density of the water is 1000 kg/m3, the pressure and the velocity in water are known and are those of a solitary wave, computed by the Tanaka method ([2] ). The density of air is 1 kg/m3, the pressure is of 105 Pa and the velocity can be neglected.

Here is a FORTRAN program to compute the Tanaka solution soliton.f

The input of this program is the free surface profile computed in tanaka.dat, and a list of points given by the participants where they wish to compute the initial condition pcal.dat. The program return a list of points with the density, the two components of the velocity and the pressure in plotuv.dat.

The values given by “soliton.f” have a dimension: meters for the positions x and y, kg/m3 for the density r, m/s for the components of the velocity u and v and Pascal (kg/m/s2) for the pressure p. The crest of the solitary wave is supposed to be at the position x=0. The participants are invited to translate the initial condition according to their needs (for comparisons for example). The velocity of the solitary wave is .

For information: here is the program (written by  S. T. Grilli) that computes the surface profile. Remarks and suggestions are welcome ( ).


On the solid boundaries: top, bottom, left and right of the computation domain, the participants will use solid wall boundary conditions.


For this part, we suggest to use a domain –13h1<x<40h1 and –h1<y<2h1. A the initial time, the crest of the solitary wave is at x=0 and propagates to the right at a velocity of . The participants are allowed to modify the size of the computation domain if the CPU time is too long. The objective is to evaluate the attenuation of the numerical solitary wave during the propagation.

The “exact” translated wave should be compared with the numerical one when it reaches x=20h1. The total kinetic energy in the domain should also be plotted during time (if the flow is supposed to be perfectly incompressible, it should be constant).

For information: here (balot.f) is the program that I intend to use to perform this part of the workshop.


The participants should plot the elevation of the free surface during time at the gauges P2, P3 and P4. The free surface profile at the breaking point will also be compared with experimental data.

The experimental time evolution of the water elevation is stored in the following files: p1_exp.txt , p2_exp.txt , p3_exp.txt , p4_exp.txt . In that files, the height is expressed in meters and the time in seconds.


The experimental surface profiles (a) and (b) are also given in the following files: prof1_exp.txt , prof2_exp.txt . The distances are in meters.

The positions of the artificial boundaries of the computation domain can be chosen by the participants, but a mesh example can be found in balot.f.


April 30th 2004: intention to participate.

May 28th 2004: send to the organizer a short paper (about 5 pages) describing the required results.

June 23st 2004: presentation of the results and discussion in the afternoon.


Preliminary results of the participants to the workshop are given here: soliton/solitonpaper/solitonpaper.htm


[1]   Yasuda, T., Mutsuda, H. and Mizutani, N., Kinematic of overtuning solitary waves and their relations to breaker types, Coastal Engineering, 29 (1997), 317–346.

[2]   Tanaka, M., The stability of solitary waves. Phys. Fluids, 29 (1986), 650–655.

[3]   Cocchi, J. P. and Saurel, R., A Riemann problem based method for the resolution of compressible multimaterial flows. Journal of Computational Physics, 137 (1997) n2, 265–298.



Philippe HELLUY