2.1. Flume configuration

Two test conditions were specified for the test to meet the requirements set out above.

1. (i) Constant water depth (C01 in table 1) of 0.6 m, which corresponds to depth (non-dimensional) $k_mh\approx 3$ , where $k_m$ is the mean wavenumber of the dispersive wave group. The condition is thus in the deep to intermediate water depth range.

2. (ii) Variable water depth (C02). A 4 m long trapezoidal-shaped shoal is placed on the seabed below the measurement probes, as illustrated in figure 1. The shoal covers the entire width of the wave flume and rises 40.8 cm above the seabed, giving 19.2 cm of water over the top of the shoal in still water. This inhomogeneous condition has a depth range $0.9\lt k_mh\lt 3$ , which spans from deep water to shallow. The slope of the shoal facing the wave maker has an angle of $\approx 30$ °, while the back slope is $\approx -50$ °.

Figure 1.

Wave flume test set-up.

2.2. Wave input

Dispersive wave groups of limited duration were generated by frequency focusing (see Rapp & Melville Reference Rapp and Melville1990). As shown in Lindgren (Reference Lindgren1970) and discussed by Tromans, Anaturk & Hagemeijer (Reference Tromans, Anaturk and Hagemeijer1991), the expected shape of the steep crest tends to the autocorrelation function of the process (surface elevation). Using focused wave groups also has several other practical benefits. The focus point $(x_p)$ can easily be set anywhere in the basin. Away from the focus point, the wave group will be out of focus, avoiding steep waves and unwanted breaking upstream where the wave is generated. This makes the waves easy to generate (experimentally and numerically) and very repeatable. Secondly, in steep conditions, sideband instabilities can cause amplitude modulation due to energy leakage to the peak frequency sidebands (see Benjamin & Feir Reference Benjamin and Feir1967). These instabilities can be difficult to control in an experimental setting, and therefore a challenge to repeatability. The transient nature of focused waves, however, confines the nonlinear interactions in time and space near the focus point, improving repeatability, whether the waves be strongly nonlinear or breaking. Finally, focused wave groups may have a short duration, eliminating the possibility of contamination by reflecting waves. The linear surface elevation of a focused wave train, given time $t$ and position $x$ may be described by

(2.1) \begin{equation} \eta (x,t)=\sum _{i=1}^N a_i \cos\left (k_i(x-x_{f})-\omega _i(t-t_{f})\right ), \end{equation}

where $a_i$ , $k_i$ and $\omega _i$ are the amplitude, wavenumber and frequency of each frequency component of the irregular wave train, and the two constants $x_{f}$ and $t_{f}$ are used to control the position and time of linear focus.

Any spectral shape may be used to generate the focused wave train; however, the spectral bandwidth and decay rate will influence the degree of nonlinear interaction observed (see Johannessen & Swan Reference Johannessen and Swan2001). In these tests, we chose to use a discrete amplitude spectrum based on a simple signal decay formulation, as previously applied by Johannessen (Reference Johannessen1997) and Johannessen & Swan (Reference Johannessen and Swan2001) due to its simplicity:

(2.2) \begin{equation} a_i=\omega _i^{-2}. \end{equation}

The amplitude spectrum is normalised and scaled to achieve the desired target amplitude $A_{targ}$ at the focus point:

(2.3) \begin{equation} a_i = A_{targ}\frac {\omega _i^{-2}}{\sum _j \omega _j^{-2}} .\end{equation}

The discrete frequency range is given by

(2.4) \begin{equation} \omega _i=\frac {2\pi }{64}i\sqrt \alpha ,\end{equation}

where $i=46{-}106$ . The spectral shape and frequency band, which corresponds to spectrum B in the work of Johannessen (Reference Johannessen1997), were chosen because it has previously been shown to exhibit a substantial nonlinear behaviour as wave steepness increases. Thus, it may serve as a challenging test case for a numerical wave model. The scaling factor $\alpha$ was added to the original definition, to scale the frequency range to match the effective range of the wave tank. This was achieved with a factor $\alpha =3/4$ , which gave a frequency range $f \approx 0.9{-}1.9\,\text{Hz}$ .

The resulting amplitude spectrum has a mean wavenumber of $k_m\approx 6.4 \,\text{m}^{-1}$ , obtained using spectral weighted averaging (see Drazen et al. Reference Drazen, Melville and Lenain2008). This gives a mean wave period $T_m\approx 0.79$ s, which is in good agreement with the measured trough-to-trough period. Four test conditions were selected to cover the steepness span. The same four wave packets were run for both bathymetry conditions described above. The full range and combination of the set-up and conditions tested are presented in table 1.

Table 1.

Test matrix of focused wave groups, run in the wave flume. Asterisk denotes conditions where PIV measurements were taken. Consequently, the ID tag for each condition is used in the results section and should be read: C01A25F01 = condition 1 (flat bed), linear target amplitude 25 mm, focus point 1.

The linear target amplitude $A_{targ}$ serves as a gain to scale the steepness of the wave group, $A_{targ} k_m$ . Note that the measured peak steepness of the wave groups will generally be higher than the reported linear steepness in table 1, due to nonlinear wave interaction (see e.g. Longuet-Higgins Reference Longuet-Higgins1962; Sharma & Dean Reference Sharma and Dean1981). Secondly, wave interaction and modulation instabilities will result in a shift of the focus point downwave (Rapp & Melville Reference Rapp and Melville1990; Johannessen & Swan Reference Johannessen and Swan2001), compared with the linear description. The generation of the amplitudes determined above is done through the wave-maker transfer function, which creates the demand signal for the wave paddle. The linear focus point was calibrated to $x=8.1\,\text{m}$ , measured from the position of the wave paddle at rest (which is consistently used as the reference point $(x=0)$ in all tests). However, for the steepest conditions, the focus point was set at $x=7.1\,\text{m}$ , to ensure that all wave breaking had occurred before reaching the rear wave gauge (gauge 4), so the energy dissipation of wave breaking can be studied.

The travel distance from the wave maker to the point of linear focusing corresponds to approximately $8.3\lambda _m$ , where $\lambda _m$ is the mean wave length of the wave group. The long propagation distance and the nonlinear behaviour of the steeper wave groups make the set-up ideal for the validation of nonlinear dispersive effects. In addition to the well-defined discrete definition of the input wave, the horizontal motions of the piston wave maker were recorded, such that its position and velocity may be used as an alternative way of prescribing boundary input to wave models.

2.3. Wave probes

Four highly accurate non-intrusive acoustic wave probes from General Acoustics (Ultralab USS 02/HFP) were used in the experiments, which can sample at 250 Hz at 0.1 mm resolution/accuracy. The probes do not disturb the wave field as traditional wire gauges. The system was configured to operate at a sampling frequency of 100 Hz, which exceeds the minimum requirement necessary to accurately resolve surface elevations. The noise levels (standard deviation) of each wave gauge were quantified, measuring between 0.026 and 0.04 mm in still-water conditions. Wave gauge 1 was placed near the wave maker, while the remaining probes were carefully placed around the point of linear wave focus, as illustrated in figure 1. The locations of the probes are specified relative to the position of the wave maker at rest.

It is important to note that in some instances of steep wave slopes, minor signal dropouts may occur when utilising acoustic probes. However, these dropouts can typically be identified and rectified using standard interpolation techniques. This approach was also applied during the steepest conditions test in this experiment.

2.4. Particle image velocimetry

Particle image velocimetry (PIV) measurements were made for the constant-depth condition under wave gauge 3 ( $x=8.15\,\text{m}$ ). Two cameras were used, vertically stacked and mounted 60 cm from the glass sidewall of the flume. The top camera was mounted at the height of the still-water line to obtain the sharpest possible surface interface. The bottom camera was mounted to cover the lower part of the wave profile with an overlapping field of view to the top camera. The two cameras covered a range from $-0.4\,{\rm m}$ to the surface. The following equipment was used:

1. (i) Camera: 2× Photron WX100 with resolution of 2048 × 2048 pixels at 125 frames per second.

2. (ii) Lens: Zeiss.

3. (iii) Light source: red LED light focused through a lens, mounted underneath the wave flume.

4. (iv) Particles: $50\,{\unicode{x03BC}}\text{m}$ PSP polyamid seeding particles.

Post-processing was performed using MATLAB, with the HydrolabPIV toolbox developed and maintained by UiO (Kolaas Reference Kolaas2016).

2.5. Assessment of model test dataset

2.5.1. Non-breaking energy dissipation

In a physical wave flume, there are several sources of wave energy loss in addition to wave breaking. Viscous dissipation, bottom friction, wall friction and surface contact line damping along the walls may all be contributing factors which drain energy from a wave group as it propagates down the flume. Quantifying the loss of energy in the experiments is critical if we wish to reproduce the results numerically.

Following the methodologies outlined by Rapp & Melville (Reference Rapp and Melville1990), Drazen et al. (Reference Drazen, Melville and Lenain2008), Tian et al. (Reference Tian, Perlin and Choi2010) and Craciunescu & Christou (Reference Craciunescu and Christou2020), our objective is to quantify energy dissipation by performing a temporal integration of the measured surface elevation time series. At each gauge location, the time-integrated potential energy per unit length can be calculated:

(2.5) \begin{equation} E_p=\frac {1}{2}\rho _w g \int _{t_1}^{t_2}\eta ^2 \,\text{d}t ,\end{equation}

where $\eta$ is the elevation of the surface, g is the acceleration of gravity, $\rho _w$ is the density of water and $[t_1,t_2]$ is the integration length, which encompasses the entire wave group. Following the approach of Drazen et al. (Reference Drazen, Melville and Lenain2008) and Tian et al. (Reference Tian, Perlin and Choi2010), under the linear assumption of equipartition of energy, the total energy $E$ is given as twice the potential wave energy $E_p$ . As such, the energy loss $\Delta E$ may therefore be estimated from surface elevation measurements only as

(2.6) \begin{equation} \Delta E = \frac {E}{E_{0}}, \end{equation}

where $E_{0}$ is taken as the energy estimated in wave gauge 1. Note, however, that the equipartition assumption between the kinetic and potential energy is strictly a linear wave assumption, and is only reasonable to use for non-breaking waves.

The resulting energy loss as a function of position in the tank can be seen in figure 2, where non-breaking focused waves of variable steepness were run. Note that in these particular tests, wave gauges 2–4 were evenly distributed in the upper half of the tank. For our particular set-up of focused wave groups with a constant depth of water of 0.6 m, the energy loss was found to be $\approx 1.5\,\,\%$ per wavelength ( $\lambda _m$ ). The dissipation rate is higher than expected; however, similar values have been reported in similar laboratory experiments (see Tian et al. Reference Tian, Perlin and Choi2010). The water depth used in the experiments is deep relative to the wave frequencies used, and the longest wave components are only marginally affected by the finite depth. This eliminates bottom friction as a potential source of energy loss. The most plausible cause is wall friction and surface contact line damping. During the experiments, it was indeed observed that small-amplitude, high-frequency waves (ripples) were generated along the glass walls as the focused wave group passed. A more detailed assessment of the wall-friction contribution and its regime of validity is provided in Appendix A. Note that this test was only carried out for the constant-water-depth condition, since the shoal will introduce reflecting waves due to wave diffraction.

Figure 2.

Wave energy dissipation ( $E/E_0$ ) plotted for various positions in the tank (x axis), for non-breaking focused wave groups of different steepness, at constant water depth. The dissipation rate is observed to be linear and independent of amplitude.

2.5.2. Repeatability

Laboratory experiments are inherently subject to uncertainties arising from factors such as measurement noise, signal dropout, equipment precision limits, mechanical constraints, wave reflections and more. For the dataset generated in this experiment to be used effectively for validation purposes, it is crucial to conduct repeated trials and document the variability in the results. Consequently, all test datasets were replicated a minimum of three times. Figure 3 illustrates three repetitions of event C01A65F02, which corresponds to a breaking-wave event where we anticipate significant variability in the measurements due to the effects of wave breaking. Both the measured wave elevation (figure 3 a) and the velocity profile extracted from the PIV measurements beneath the peak crest (figure 3 b) demonstrate excellent repeatability.

Figure 3.

Three repetitions of the same condition (C01A65F02), showing measured surface elevation (a) and the measured velocity profile (b) underneath the passing crest at wave gauge 3 (at $t=0.544 s$ ).

2.6. Applicability and limitations

The experimental dataset comprises eight selected wave events, deliberately chosen to span a broad range of wave steepness and water depths, extending from deep to shallow conditions. The wave maker was carefully calibrated to ensure accurate reproduction of linear wave components prior to generating steeper, nonlinear wave states. The dataset includes both the prescribed linear frequency components and the corresponding paddle motion signals, thereby enabling numerical reconstruction of the wave fields using different modelling approaches. Consequently, the dataset is suitable for validation and performance assessment across a wide class of wave models.

The experimental configuration is deliberately narrow and idealised, with compactness and numerical reproducibility prioritised over physical generality. Although the focused wave groups exhibit dispersive characteristics and include wave breaking, key features of natural ocean wave fields, such as directional spreading and stochastic phase variability, are excluded.

Given that the test cases involve strongly dispersive wave conditions, non-dispersive or weakly dispersive formulations (e.g. shallow-water equations or certain Boussinesq-type models) may face inherent limitations in accurately reproducing these events. Furthermore, models constrained to constant water depth cannot represent the spatially varying bathymetry considered in the experiments. In particular, the steepest inhomogeneous wave conditions investigated in these experiments are expected to be reliably represented only by high-fidelity, phase-resolving numerical models capable of resolving wave-breaking processes. Accordingly, nonlinear wave modelling approaches, such as potential flow models (e.g. boundary element methods and higher-order spectral methods) or layered grid-based formulations (e.g. multi-layer non-hydrostatic models), would be anticipated to reproduce a substantial portion of the dataset, although they may encounter limitations in accurately predicting the onset of breaking. The dataset is therefore well suited for systematically evaluating and screening the capabilities of different phase-resolving modelling frameworks.

It should be emphasised that the experiments were performed in a wave flume and are therefore restricted to long-crested wave conditions. The finite lateral width of the facility may suppress fully three-dimensional dynamics, including lateral energy spreading, oblique breaking and cross-wave instabilities.

Nevertheless, the flume width was considered sufficiently large to prevent substantial confinement effects on the breaking process itself and on the associated turbulence generation, since the characteristic turbulent length scales in the post-breaking region are significantly smaller than the basin width. However, a wider basin could potentially reduce the effect of non-breaking energy dissipation on the measurements, and loosen the three-dimensional dynamics restrictions of the test. Consequently, conducting equivalent experiments in a wider basin and extending the dataset to include short-crested wave conditions constitute logical and valuable directions for future work. Similarly, incorporating additional and more diverse bathymetric configurations such as shorelines and bar formations into the experimental set-up would further broaden the applicability and relevance of the dataset.

2.7. Availability

The experimental data presented are made available open access in the University of Oslo’s data repository DataverseNO (Lande Reference Lande2026), and can be accessed at https://doi.org/10.18710/FZN3AL.

3.1. Numerical domain set-up

A two-dimensional numerical set-up was selected to accurately model the entire wave flume used in the experiments at high resolution. The domain dimensions are 25 m in length and 1.56 m in height, with a water depth of 0.6 m. At the far end of the domain, a wave-absorbing damping zone is implemented, corresponding to the location of the dissipating beach used in the experiments.

The total simulation time was set to 50 s, which is sufficient to ensure complete generation of the prescribed wave group at the inflow boundary and its subsequent propagation through the central measurement region of the basin.

Although the selected CFD framework supports fully three-dimensional simulations, a two-dimensional configuration was adopted to achieve high spatial grid resolution within feasible computational limits. This approach necessarily neglects inherently three-dimensional processes, but was considered justified given the prohibitive cost of equivalently resolved three-dimensional computations.

The AMR was updated at each time step based on local velocity gradients within the computational domain. In addition, the highest level of refinement was imposed along the entire free surface to ensure a consistent and accurate extraction of surface elevation data. This strategy results in a substantial increase in the total number of computational cells and, consequently, the computational cost. However, such an increase remains computationally manageable in two-dimensional simulations. Relaxing this refinement constraint would make three-dimensional simulations of the same conditions feasible, although they would remain computationally demanding.

3.2. Boundary conditions

An open inflow boundary condition (Dirichlet-type velocity condition) is prescribed to facilitate the generation of the focused wave groups. The prescribed flow velocity is derived from the time derivative of the paddle’s position. The linear focus point of the wave group is located 8.1 m downwave, resulting in the wave group being out of focus at the wave paddle. The stroke length required to generate these waves is therefore relatively small compared with the wavelength, ensuring that no wave breaking occurs in the immediate vicinity of the wave paddle or inflow boundary. This supports the adoption of an open stationary flow boundary to simulate the motion of the wave paddle’s moving wall.

As an alternative to prescribing inflow kinematics based on the wave paddle signal, linear or second-order wave theory may also be used, employing the spectral definitions provided in § 2.2. Both approaches have been shown to yield results that align well with the experimental data.

To maintain consistency between the pressure gradient and the fluid acceleration at the boundary, the pressure gradient is computed using the following equation:

(3.6) \begin{equation} \boldsymbol{\nabla }p = -\rho \left (\frac {\partial u_{n}}{\partial t}\right ), \end{equation}

where ${\partial u_n}/{\partial t}$ represents the normal velocity gradient at the boundary face. An open boundary condition is applied at the top of the flume to ensure a balanced inflow and outflow, achieved by imposing a Dirichlet-type pressure boundary condition on the top boundary along with a Neumann condition for velocity.

The shoal is modelled using an embedded boundary with a no-slip condition applied parallel to the slope.

3.3. Wall friction

To account for energy loss due to wall friction in the experiments, we introduced a simple kinetic energy dissipation model, defined as follows:

(3.7) \begin{equation} \boldsymbol{u}^*=\frac {\boldsymbol{u}}{1+c_{f}\Delta t}, \end{equation}

where the linear wall friction coefficient $c_f \approx 0.0119$ was determined to yield the same amplitude-independent dissipation rate as observed in the experimental measurements for non-breaking waves (see figure 2). Here $\Delta t$ is the time-step size that varies during simulation and $\boldsymbol{u}$ and $\boldsymbol{u}^*$ represent the velocity fields before and after the application of the wall friction model. Consequently, the wall friction coefficients found for non-breaking waves were used for all simulations.

3.4. Scales, energetics and grid resolution

The numerical modelling of wave evolution, particularly in regimes involving wave breaking, constitutes a challenging multi-scale problem. Accurate representation of the relevant physical processes requires sufficient resolution across a range of spatial and temporal scales. Consequently, the choice of grid resolution can substantially influence the quality of the simulation results. To evaluate the adequacy of the numerical discretisation and to characterise the dominant physical scales, a systematic analysis of characteristic scales and energetics has been conducted.

In addition, simulation of wave propagation in two-phase flows with large density contrasts between air and water introduces further numerical complexity. The high density ratio across the interface, coupled with shear flow along the free surface, imposes stringent requirements on interface capturing, stability and resolution. These factors make the problem particularly demanding from a numerical perspective, as discussed by Desjardins & Moureau (Reference Desjardins and Moureau2010) and Remmerswaal (Reference Remmerswaal2023). A discussion of scales and energetics related to grid resolution is presented in Appendix B.

4.1. Linear dispersion and energy conservation

We first assess the ability of the numerical model to retain the dispersion properties of the wave groups. Even though the chosen model is fully nonlinear and should incorporate all physical processes of steep and breaking waves, the accuracy and performance of a CFD wave model is influenced by the choices of numerical schemes and the various parameters that control the numerical simulation. We start by comparing surface elevation measurements for the lowest-amplitude case, C01A07F01, in constant deep water. Linear dispersive wave behaviour is expected. As such, most wave models for irregular waves will be able to reproduce this event. The results for the four wave gauges are presented in figure 5(a). All wave elevation comparisons have been normalised using the target amplitude $A_{targ}$ when presented. Consequently, the red lines represent the CFD wave model, the blue lines represent the linear wave theory and the black dots represent the experimental measurements. There is agreement between the two wave models, albeit the low steepness of this condition being close to the measurement accuracy of the acoustic wave probes. Given the propagation distance of $\approx 8.5 \lambda _m$ (measuring from wave paddle to wave probe 4), both the amplitude and the phase of the waves are very well reproduced by the CFD wave model, as well as the linear wave theory. Although the results are as expected, they serve as confirmation of measurement accuracy as well as conservation of energy in our CFD model.

Figure 5.

Comparison of normalised wave elevation measurement for linear (a) and weakly nonlinear (b) test case, run for deep-water conditions with constant depth.

4.2. Nonlinear dispersive effects

As steepness increases, nonlinear wave interaction will occur, which changes the shape and modulates the phases of the dispersive wave group. This is indeed a well-known effect on gravity waves. Reference is made to the pioneering work of Hasselmann (Reference Hasselmann1962, Reference Hasselmann1963), Longuet-Higgins & Stewart (Reference Longuet-Higgins and Stewart1960), Phillips (Reference Phillips1960), Longuet-Higgins (Reference Longuet-Higgins1962) and Longuet-Higgins & Phillips (Reference Longuet-Higgins and Phillips1962). Figures 5(b) and 6(a) compare the surface elevation for our next two conditions, which are nonlinear but non-breaking waves. For these cases, we expect that linear and weakly nonlinear wave models will struggle. Although changes in amplitudes and wave shape to a large extent are well modelled by weakly nonlinear wave theory, such as second-order irregular wave theory, the nonlinear modulation that occurs is not. Linear wave theory, which we in this case compare against, does not do particularly well, as expected. Our CFD wave model, on the other hand, is a fully nonlinear wave model, and as such we expect the same modulation and amplification as the experiments. This is indeed confirmed when we compare the surface elevations for the non-breaking moderate test cases (figure 5 b) and steep test cases (figure 6 a). The nonlinear amplification and modulation which occur as the waves steepen are well captured by the CFD model, whereas the deviation between the experiments and linear wave theory increases with steepness.

Figure 6.

Comparison of normalised wave elevation measurement for steep test case at the limit of breaking (a) and the severely breaking case (b) for deep-water conditions with constant depth.

4.3. Breaking event

We then look at the steepest wave condition tested, where the wave energy is sufficient to cause successive breaking. Snapshots from simulation of the test case C01A65F02 are shown in figure 7. Visually we observe a substantial increase of vorticity in both phases and air–water mixing, leading to energy dissipation, which is studied later. The corresponding wave elevation measurements for this event are compared in figure 6(b). All probes show excellent agreement between experiments and the CFD model, before and after wave breaking. As expected, this cannot be said for the linear wave model, which does not account for nonlinear interaction, nor wave breaking.

Figure 7.

Snapshots from simulation of test C01A65F02 during the onset of wave breaking. The simulation was run with grid refinement level 14. In the water, the colouring indicates the velocity magnitude; in the air, the vorticity is shown (in greyscale).

4.4. Inhomogeneous seabed

The same four wave groups are repeated, introducing the shoal in the centre of the flume, as shown in figure 1. With a slope of $\approx 0.3$ , the incoming waves will transcend from deep water $kh\approx 3$ to shallow water $kh\approx 0.9$ in a short distance of $\approx 1.5 \lambda _m$ . The combination of dispersive waves and a reasonably quick reduction in water depth makes these test cases challenging for wave models, which are limited to mild slopes, weakly dispersive or weakly nonlinear wave models. Our CFD wave model should fully incorporate these properties and is thus expected to perform well. Figures 8(a) and 8(b) compare surface elevation measurements from the two steepest conditions, C02A45F01 and C02A65C02, which are both breaking-wave conditions as the wave reaches the shoal. Note that the linear wave model has been omitted for this comparison, since the test condition is beyond what can realistically be modelled using linear wave theory. As a reference, we have instead plotted the time series from the constant-water-depth condition (C01) as grey dots to show the substantial influence of the shoal on the wave train. Again we see that our chosen numerical wave model is capable of reproducing the measured surface elevations from the experiments quite well, considering that the waves are breaking as they pass gauges 2 and 3.

Figure 8.

Comparison of wave elevation measurement for case C02A07F01 (a) and C02A65F02 (b). Grey dots show the time series from the corresponding experiment without the shoal, to illustrate the large influence of the shoal on the wave elevations.

4.4.1. Velocity profiles and momentum flux

From processing of the recorded PIV measurements, velocity profiles under wave gauge 3 are obtained for three of the eight test conditions. A comparison of the velocity profile underneath the largest crest is presented for two of the events in figures 9(b) and 10(b). For the larger-amplitude case, we observe a substantial deviation of velocity. For the moderate-amplitude condition the agreement is quite good for the entire profile. However, for the large-amplitude condition, the measured velocities are lower than the numerical model close to the free surface. This is caused by difficulty in obtaining accurate measurements from the PIV recording due to overexposure near the surface. The upper camera was mounted with a lens centre height level to still water. As a result, a larger projected free-surface area is exposed for large surface elevations, and particle velocities in the upper $\approx 10\,\%$ – $20\,\%$ of the larger wave crest either could not be obtained from the post-processing or were tainted by noise levels too high to filter out or reconstruct reliably. The affected vertical extent therefore scales with crest amplitude: the moderate-amplitude condition C01A25F01 is only marginally affected, whereas the breaking case C01A65F02 shows the effect clearly in the upper part of the velocity profile. Away from the surface zone, the kinematics obtained from the model and experiments agree well.

In figures 9(a) and 10(a), the momentum flux is compared, obtained by integrating the velocity profile from a depth of $z=-0.4\,\text{m}$ , which is the lower limit of the PIV image frame, to the free surface over a range of time steps. The agreement is good, as we expect from the velocity profile comparison, although the peaks are underestimated in the experiments because the missing or noisy upper-crest velocities discussed above are excluded from the integration. The effect is small for C01A25F01 but clearly visible at the peak of C01A65F02, consistent with the amplitude-dependent extent of the affected near-surface band. The PIV measurements should therefore be regarded as reliable below this band, while the upper-crest kinematics under steep and breaking waves are best inferred from the validated CFD solution. The PIV set-up could likely have been improved by adjusting the camera configuration to each individual wave condition, by using multiple cameras with overlapping fields of view or by reducing the shutter speed; however, this was not attempted during the present measurement campaign.

Figure 9.

Comparison of momentum flux (a) for test case C01A25F01. (b) A comparison of the instantaneous velocity profile underneath wave gauge 3, at the time of maximum momentum.

Figure 10.

Comparison of momentum flux (a) for C01A65F02. (b) A comparison of the instantaneous velocity profile underneath wave gauge 3, at the time of maximum momentum.

4.5. Non-breaking energy dissipation

Energy conservation poses a challenge for iterative numerical fluid flow models. It is well established that the choice of numerical schemes and simulation parameters can significantly influence energy conservation (see Buist et al. Reference Buist, Sanderse, Dubinkina, Henkes and Oosterlee2022). The Navier–Stokes equations, which our solver is based upon, describe the conservation of mass and momentum; however, as noted by Mostert & Deike (Reference Mostert and Deike2020), the methods utilised in Basilisk are not energy-conserving by design. Consequently, accurately recovering the full energy budget can be difficult, despite energy conservation following from the mass and momentum conservation equations as a secondary conserved quantity (Buist et al. 2022).

In this context, our objective is to propagate waves over considerable distances, which makes it particularly important to document energy conservation. Figure 11(a) illustrates the loss of potential energy along the length of the wave flume for different levels of grid refinement. The energy ratio ( $E/E_0$ ) has been calculated based on wave elevation measurements, as discussed in § 2.5.1. The simulations were conducted for condition C01A25F01.

For the solid-line cases, the wall friction model (described in § 3.3) has been disabled to illustrate that numerical dissipation within the simulation model is minimal, resulting in negligible energy loss. It is important to emphasise that achieving such accurate energy conservation required setting the Courant–Friedrichs–Lewy (CFL) number, which limits the time-step size $\delta t$ based on cell size and fluid velocity, to 0.1. The dashed lines depict the comparison between the experimental data and the numerical results when the wall friction model is enabled. Dissipation rates are reported for various grid refinement levels, demonstrating that non-breaking energy conservation is largely insensitive to mesh refinement.

Figure 11(b) compares the potential energy loss across different amplitudes for condition C01. The three non-breaking cases $A_{targ}=7$ , $25$ and $45\,\text{mm}$ closely align with the dissipation rate obtained from the model test. However, for the largest wave group $A_{targ}=65\, \text{mm}$ , breaking occurs, resulting in a substantial energy drop, which is discussed further in the next section.

Notably, a fluctuating component is observed, attributed to kinetic–potential energy fluctuations that directly stem from the use of focused wave groups. As the waves focus, the assumption of equipartition of energy is violated, leading to the observed fluctuations in potential energy. These fluctuations are amplitude-dependent and diminish as one moves away from the focus point. A similar inverse fluctuating component can also be detected by analysing the kinetic energy as a function of position.

This observation highlights that the assumption of equipartition of energy may be a poor approximation for assessing total energy when using focused wave groups, even in the case of non-breaking waves.

Figure 11.

Potential energy loss $E/E_0$ , plotted as a function of position X in the wave flume. (a) Simulation results for condition C01A25F01 for different grid refinement levels. Solid lines show the CFD model run without a wall friction model, while the dashed lines are with the friction model included and are comparable to the experiments. (b) Comparison of dissipation for all four waves of condition C01, run with the wave model (solid lines), with the energy dissipation of the non-breaking experiment waves (dashed line).

4.6. Breaking-wave energy dissipation

Three of the eight test conditions resulted in significant wave breaking due to energy focusing and depth-induced breaking (specifically, for condition C02). For all three breaking events, the majority of the onset of breaking occurs before reaching wave gauge 4. Although the linear assumption of equipartition of energy is not applicable to steep breaking waves or shallow-water conditions, for the sake of comparison and considering our focus on relative values, the variation of potential energy along the length of the tank is deemed relevant and is anticipated to exhibit consistency between the model and the experimental results.

Figure 12.

(a) Potential energy dissipation for breaking wave events C02A45F01 and C02A65F02 for refinement level 14. Parameter $E_0$ , used to normalise the energy loss, is obtained from time integration of wave gauge 1. (b) Estimated energy loss due to wave breaking.

The resulting energy loss $\Delta E$ as a function of position $x$ is presented in figure 12(a) for the breaking-wave events C01A65F02, C02A45F01 and C02A65F02, using grid refinement level 14. The star markers indicate the energy loss calculated from the experimental data. A notable decrease in wave energy is observed for both the model and the experiments within the transitional breaking zone, which spans approximately from $x=6$ to $x=8$ m.

From the model results, which feature a higher density of positional measurements, we observe fluctuations consistent with the expectation that we are only examining potential energy. Therefore, we should anticipate a comparable degree of variability in the experimental results; however, the limited number of wave gauges is insufficient to fully capture this variability.

In addition, the energy loss attributed to wall friction is illustrated as a black dotted line. This information, shown in figure 2, was previously determined to be independent of amplitude. To estimate the energy dissipation due to breaking, we consider the difference in energy between non-breaking dissipation (represented by the dashed black line) and the energy associated with breaking wave events, evaluated at wave gauge 4. The resulting estimates are presented in figure 12(b) for our three breaking conditions across different levels of grid refinement.

The agreement between the experimental measurements and the numerical wave model is generally quite good, even at the coarser refinement levels of 12 and 13. However, at level 11, the grid becomes too coarse to accurately model the breaking waves.

In general, these results are encouraging, especially considering the limitations in turbulent scales resolved by the simulation model, as discussed in Appendix B. This is particularly relevant given the substantial computational costs often associated with increasing mesh resolutions in CFD wave models. While grid refinement level 14 is feasible for use in a two-dimensional domain, the associated costs would likely be prohibitive when extending to three-dimensional domains and short-crested seas. Consequently, the satisfactory performance observed at grid refinement level 12 is indeed a positive finding.

4.7. Benchmarking the wave models

To evaluate the performance of various wave models under differing wave conditions, or to assess the effectiveness of a mesh-based iterative wave model across different refinement levels, a benchmark is essential. Using the normalised surface elevation ( $\hat {\eta } (x,t)$ ), the relative error between the experimental data and the wave model can be determined using a weighted least-squares approach:

(4.1) \begin{equation} \epsilon =\sum _{j=1}^M \sum _{i=1}^N w_{ij}\big(\hat {\eta }^{exp}_{ij}-\hat {\eta }^{model}_{ij}\big)^2, \end{equation}

where the weighting function ( $w_{ij}=|\hat {\eta }^{exp}|$ ) is applied to penalise errors associated with larger wave excursions. In this equation, the indices i and j represent the discrete time and position, respectively. Consequently, the computed least-squares errors are treated as an equally weighted sum of the errors from the four wave gauges.

The notation of the hat in $\hat {\eta }$ indicates that the surface elevation time series have been normalised against the target amplitude ( $A_{targ}$ ) for each test condition. Thus, the calculated errors are relative to this target amplitude. The resultant single number error ( $\epsilon$ ) serves as a performance indicator for the wave model under the specified test condition, with lower values indicating better performance. Although we do not expect the error to approach zero, we anticipate it to stabilise at a low value. Exact matches cannot be realised due to inherent uncertainties in both the model tests and the numerical wave models.

Benchmarking is conducted on four levels of mesh refinement, as described in table 7. As additional reference, we also report the performance of three further wave models: linear irregular wave theory, second-order irregular wave theory and a depth-averaged shallow-water equation (SWE) solver. Linear and second-order theory are applied only to the constant-depth test cases (C01), as their underlying assumptions do not accommodate the shoal geometry used in the C02 series. The SWE model is applied to both series. The corresponding surface-elevation time series are shown in Appendix D, and the benchmark results are summarised in table 2.

Table 2.

Benchmark error $\epsilon$ for different grid refinement levels (L11–L14) and reference models (linear, second-order, SWE). Smaller values indicate better agreement with the experiment.

Generally, the differences observed among the three highest levels of refinement are minimal, confirming convergence regarding surface elevation. However, at refinement level 11, a noticeable decline in performance is observed. The best results are obtained for the 25 mm linear amplitude test conditions, which represent weakly nonlinear waves. Although it might be anticipated that the linear test cases (7 mm amplitude) would perform equally well, the very low amplitude renders these conditions more challenging due to a low signal-to-noise ratio in the model test measurements. Given that benchmark errors are reported relative to the amplitude, these results are not surprising.

For the steepest conditions tested, the benchmark results are slightly inferior to those for the weakly nonlinear case, yet remain robust across the three highest levels of refinement, particularly considering that these conditions are strongly nonlinear and involve wave breaking. Comparisons of wave elevation measurements for different refinement levels under test conditions C02A25F01 and C02A65F02 are illustrated in figure 13, visually corroborating the findings indicated by the benchmark results.

Figure 13.

Comparison of wave elevation $\hat {\eta }$ , computed using different grid refinement levels, for almost linear (a) and severely breaking (b) conditions.

In contrast, the performance of linear wave theory demonstrates a marked deterioration as the amplitude increases, as previously observed in the wave elevation comparisons shown in figures 5 and 6.

The benchmark errors obtained with second-order irregular wave theory are very similar to those of the linear theory across the C01 series, and likewise deteriorate sharply with increasing amplitude. Although second-order theory produces a more realistic crest–trough asymmetry and bound-wave structure than the linear model, it does not capture the nonlinear amplitude and frequency modulation that develops as a dispersive wave group steepens and approaches focus. In the benchmark metric, where the elevation is sampled at fixed gauge positions and weighted by the experimental amplitude, the residual error is therefore dominated by a phase offset between the numerical and experimental wave groups, and is only weakly reduced by the improved wave shape. This is visible in the surface-elevation comparisons in figure 15.

The SWE model, by contrast, produces benchmark errors that are two to three orders of magnitude larger than those of the Navier–Stokes model and do not improve with increasing amplitude. Because the shallow-water equations assume a hydrostatic pressure distribution, they are non-dispersive, and all frequency components propagate at the long-wave celerity $c = \sqrt {gh}$ . For the deep- and intermediate-depth conditions considered here this assumption is fundamentally violated: the individual spectral components no longer travel at their respective linear phase speeds, and the focused wave group cannot be reproduced. The effect is already evident at wave gauge 1, only $\approx 1\,\textrm{m}$ from the paddle, where the SWE solution has already lost its prescribed amplitude and phase structure, as shown in figures 15 and 16. The SWE entries in table 2 are included as a baseline illustrating the penalty incurred when a non-dispersive model is applied outside its intended regime of validity.

Finally, we evaluated the contribution of the various wave gauges to the reference error $\epsilon$ , as shown in table 3. The percentage contributions are averaged across all eight test conditions and were found to be relatively consistent across the different cases. As anticipated, the contribution from wave gauge 1 is negligible due to its proximity to the wave maker. The most significant contributions arise from wave gauges 2 and 3, which can be attributed to the occurrence of wave breaking in that area, making it difficult to accurately measure surface elevation in the experiments using acoustic probes and in the CFD model due to air–water mixing. However, once the wave groups reach wave gauge 4, measurements become easier, resulting in a lower error.

Table 3.

The contribution to the benchmark error $\epsilon$ from different wave probes. The result is found to be quite consistent for all levels of refinement.

4.8. Computational performance

The computational cost of CFD wave models depends on the level of simulation refinement, as well as the resulting velocities in the fluid flow. Our wave model is no different. The time step of the iterative solver is controlled by the CFL number, which may be seen as an expression of the maximum distance that a particle is allowed to move relative to the size of the cell through which it moves in a time step. The CFL number for our simulations was set to 0.1 for all cases to ensure energy convergence. The simulations were run on FRAM, one of Norwegian Research Infrastructure Services (NRIS) high-performance computational clusters, with Intel E5-2683v4 2.1 GHz nodes, each of which has 32 threads available. Since the number of cells in the domain varies over time due to AMR, the number of threads utilised varies between 2 and 32, during runtime, to optimise performance. The total simulation time accumulated (wall clock) is calculated as

(4.2) \begin{equation} t_{awc} = t_{wallclock} \bar {N}_{thread}, \end{equation}

where $t_{wallclock}$ is the simulation time (wall clock) and $\bar {N}_{thread}$ is the average, time-weighted, number of threads used during the simulation. The accumulated simulation times (in hours) for all the simulated cases are presented in table 4.

Table 4.

Accumulated simulation time $t_{awc}$ , computed for all simulated cases. Values are given in number of hours.

It is evident that the computational cost of simulating each event varies significantly, both across different refinement levels and among individual cases, depending on the steepness of the conditions. Although the three highest grid refinement levels generally yield very similar results, it is noteworthy that the cost of running the steepest conditions at level 12 is less than 5 % of that at level 14.