The Efficient Application of an Impulse Source Wavemaker to CFD Simulations

: Computational Fluid Dynamics (CFD) simulations, based on Reynolds-Averaged Navier–Stokes (RANS) models, are a useful tool for a wide range of coastal and offshore applications, providing a high ﬁdelity representation of the underlying hydrodynamic processes. Generating input waves in the CFD simulation is performed by a Numerical Wavemaker (NWM), with a variety of different NWM methods existing for this task. While NWMs, based on impulse source methods, have been widely applied for wave generation in depth averaged, shallow water models, they have not seen the same level of adoption in the more general RANS-based CFD simulations, due to difﬁculties in relating the required impulse source function to the resulting free surface elevation for non-shallow water cases. This paper presents an implementation of an impulse source wavemaker, which is able to self-calibrate the impulse source function to produce a desired wave series in deep or shallow water at a speciﬁc point in time and space. Example applications are presented, for a Numerical Wave Tank (NWT), based on the open-source CFD software OpenFOAM, for wave packets in deep and shallow water, highlighting the correct calibration of phase and amplitude. Furthermore, the suitability for cases requiring very low reﬂection from NWT boundaries is demonstrated. Possible issues in the use of the method are discussed, and guidance for accurate application is given. C.W. and J.D.; validation, P.S., C.W. and J.D.; formal analysis, P.S., C.W. and J.D.; investigation, P.S., C.W. and J.D.; resources, P.S., C.W.; data curation, P.S., C.W.; writing—original draft preparation, P.S., C.W. and J.D.; writing—review and editing, P.S., C.W., J.D., J.V.R. and T.W.; visualization, C.W.; project administration, J.V.R. and T.W.; funding acquisition, J.V.R. and T.W.


Introduction
While CFD enables the simulation of complex flow phenomena, such as two-phase flows and gravity waves, setting up a simulation of coastal or marine engineering processes, with correct wave creation and minimum reflection, can be a challenge. A wide range of Numerical Wavemaker (NWM) methods exist in the literature for the generation and absorption of waves in a CFD simulation. Wave generation methods can be broadly grouped into two main categories: 1.
Replication of physical wavemakers, such as oscillating flaps, paddles or pistons 2.
Implementation of numerical/mathematical techniques, introducing source terms or similar into the governing equations The first category of wave generation methods directly simulates the wavemaker by a moving wall, requiring mesh motion/deformation in the CFD domain [1], which can be a complex and computationally-expensive task. Additionally, because these methods replicate the real processes in a physical laboratory tank, they also incur the well-known drawbacks of physical wavemakers encountered in an experimental facility, such as evanescent waves, imperfect wave absorption and the required implementation of a wavemaker control system. Another drawback of conventional wave pistons is that they are typically limited to low order waves (first or second order).
In contrast to directly mimicking physical wavemakers, the second category of wave generation methods is generally more computationally efficient, using numerical algorithms to create the desired flow conditions by manipulating the field variables. In principle, there are four types of methods in the second category of wave generators: • relaxation method: relaxes the results of the simulation inside the domain with the results given by wave theory [2][3][4][5] (linear or higher order) or other numerical models like Boundary Element Methods (BEM) for example [6,7]. • boundary method: creates and absorbs waves on boundary patches through Dirichlet boundary conditions [8][9][10] • mass source method: the wave is created by adding a source term to the continuity equation [11][12][13][14] • impulse source method: the wave is created by adding a source term to the impulse equation [15][16][17][18] Correspondingly, wave absorption can also be achieved using a number of methods. The trivial solution, of using a very long tank and limiting the simulation duration to a small number of wave cycles, is certainly not efficient or elegant and will therefore not be considered herein. The relaxation and boundary wave generation methods allow for wave absorption to be achieved actively by the NWM, whereas the mass and impulse source function methods can be used in combination with numerical beaches to dissipate waves and eliminate reflections.
While a comparison of different wavemakers is not within the scope of this paper, and has been addressed by Windt et al. [19], we want to highlight some common issues in NWTs. Wavemakers relying on open boundaries or algorithms that set the species variable, such as relaxation methods, are not necessarily mass conservative. Compared to an experimental facility, where waves are typically created in an enclosed volume of water, wave run-up along a sloping bottom, for example, might not be recreated correctly. Numerical beaches, while in many cases expected to be more computationally expensive than boundary conditions, are able to absorb waves to any desired level of reflection. The user thus has the possibility to balance computational burden against permissible reflections. Mass source methods are very similar to impulse source methods, but problems can arise when generating waves of significant height in shallow water. Once the surface level falls below the source region, wave generation fails. It should be noted that, at the wavemaker location, the surface elevation is twice the wave amplitude of the target wave train, since two wave trains travelling in opposite directions are created. The present paper thus focuses on the impulse source wavemaker.

Impulse Source Wavemaker
The difficulty in utilising an impulse source as an NWM lies in calculating the required source function to obtain a desired target wave series. Since the free surface is not a variable in the RANS equations, there is no direct expression relating the impulse source function to the resulting generated wave series. However, for depth-integrated equations, the free surface does appear as a variable, which Wei et al. [20] utilised to derive a transfer function between the source amplitude and the surface wave characteristics, generating mono-and poly-chromatic waves in Boussinesq-type wave models. Following this approach, generating waves by manipulating the impulse term in the Navier-Stokes equation has successfully been applied to shallow water waves in CFD simulations by Choi and Yoon [18] and Ha et al. [15]. Ha et al. [15] describe several issues, especially when creating deep water waves, originating from the Boussinesq simplifications used in the source term derivation. They also discussed limitations inducing numerical errors when applying a momentum source function to random wave generation. Despite those promising first steps, internal wavemakers have not seen widespread application.
In theory, limitations due to the Boussinesq simplification might be overcome by employing stream functions or other high order wave theories. In practice, evaluating higher order wave theories at multiple spatial positions (across all faces of the patch or all cells in the wavemaker region) can be computationally expensive. Furthermore, in many practical CFD applications, the accurate wave trace description is required at some distance away from the wave maker, typically in the middle of the domain. Schmitt [21] presented some progress in the use of neural networks for calibrating extreme waves in shallow water conditions, and another interesting application of two internal wavemakers in a single NWT in deep and intermediate water can be found in McKee et al. [22].
The present paper follows an alternative, more generalised, method to determine the required source function. The proposed method follows from standard calibration procedures utilised for physical wavemakers in real wave tanks, which iteratively tune the input signal applied to the wavemaker, in order to minimise the error between the measured and targeted wave series. Unlike previous methods to determine the source function, which are restricted to shallow water waves, the method proposed herein is applicable from deep to shallow water conditions, able to generate any realistic wave series at any desired location in an NWT, without explicitly employing any wave theory. Additionally, while dissipation can reduce a theoretically-correct wave height from the source region as it travels to the target location, the proposed method overcomes this problem since it is designed to obtain the desired wave signal at the target location.
Iterative calibration methods will of course increase the computational burden when compared directly to an analytically-derived accurate wavemaker function. However, by performing the calibration runs on a two-dimensional slice, the computational overhead is almost negligible compared to fully-three-dimensional simulations of realistic test cases. For many scientific or industrial applications, wavemakers with sophisticated theoretical solutions are used, although transient time traces require calibration, practically making the details of the wavemaker transfer function irrelevant and even creating an unnecessary overhead [23]. Applying a calibrated two-dimensional result to a corresponding three-dimensional domain is efficient and is also the recommended approach by Perić and Abdel-Maksoud [24].
The paper will demonstrate that, by utilising the proposed method, the source function can consist of a purely horizontal impulse component that varies in time. The time evolution of the magnitude and direction of this simple impulse source term can be calibrated using a standard spectral analysis method, which is commonly used in physical wave tanks. The remaining parameters of the impulse source wavemaker that need to be chosen are then the geometrical size and shape of the source function inside the NWT domain. An investigation of the effect of these geometric parameters will be presented to offer guidance on the selection of those values.

Outline of the Paper
The remainder of the paper is organised in the following four sections. In Section 2, the implementation of the impulse source wavemaker and the theoretical background is discussed in detail. The following two sections then present the application of the impulse source wavemaker to a wave packet in deep and shallow water conditions, with Section 3 detailing the particulars of the case study and Section 4 presenting the results and a discussion. Finally, in Section 5, a number of conclusions are drawn.

Implementation
This section details the implementation of the impulse source wavemaker, within a CFD solver, based on RANS models. First, the modification of the impulse equation to produce the wavemaker is detailed. Next, the definition of the parameter values, associated with the modified impulse equation, is discussed. Finally, the calibration procedure used to generate a target wave at a specified location is presented.

The Impulse Equation
A Volume Of Fluid (VOF) RANS model for incompressible flow includes the following impulse equation: where t is time, U the fluid velocity, p the fluid pressure, ρ the fluid density, T the stress tensor and F b the external forces, such as gravity. The stress tensor can include terms from a turbulence model; however, for this paper, all cases presented use a laminar flow model, since only non-breaking waves are considered. The proposed impulse source wavemaker is implemented by adding two terms to Equation (1): • rρa wm : This is the source term used for wave generation, where r is a scalar variable that defines the wavemaker region and a wm is the acceleration input to the wavemaker at each cell centre within r. r is dimensionless, while a wm is given in the units of m s −2 . • sandρU: This describes a dissipation term used to implement a numerical beach, where the variable field sand controls the strength of the dissipation, equalling zero in the central regions of the domain where the working wavefield is required and then gradually increasing towards the boundary over the length of the numerical beach [25]. sand is given in units of s −1 .
Introducing these two terms to Equation (1) yields the adapted impulse equation:

Setting Parameter Values
The values for the parameter fields r and sand are set during preprocessing, and their definition is crucial for the correct functioning of the method.

•
r is set to one in the region where the wavemaker is acting and zero everywhere else in the NWT domain. Therefore, the size of the wavemaker and its position within the NWT must also be selected. To offer guidance on the selection of the wavemaker size and position, the case study presented in Sections 3 and 4 investigates the effect of these parameters on the wavemaker performance. • sand is initialised using an analytical expression relating the value of sand and the geometric coordinates of the NWT. The simplest expression would be a step function, where the value of sand is constant inside the beach and is zero everywhere else in the NWT. However, such a sharp increase in the dissipation will cause numerical reflections. Instead, the value of sand should be increased gradually from the start to the end of the numerical beach. Equation (3) is used in the current implementation, which has been shown to result in sufficient absorption [25]: where l beach is the length of the numerical beach, x is the position within the numerical beach, equalling zero at the start and increasing to l beach at the NWT wall, and sand Max stands for the maximum value of sand. Guidance on the selection of the parameters l beach and sand Max is also given in the case study (see Section 3.4.4).
Perić and Abdel-Maksoud [26] recently derived an analytical solution describing the ideal setting for sand and validated the method in a numerical experiment, removing the need for parameter studies in the future.

Calibration Procedure
The source acceleration time trace, a wm (t), is calibrated, using a standard spectral analysis method, based on work presented by Masterton and Swan [27], to produce a target wave series at a desired position within the NWT. Figure 1 shows a schematic of the calibration procedure, which comprises the following steps: 1.
Define target wave series at desired NWT location, η T (t), with a signal length L and N samples 2.
Perform a Fast Fourier Transform (FFT) on η T (t), to obtain the amplitudes, A T ( f j ), and phase components, Generate an initial time series for the wavemaker source term, a wm,1 (t) (can be chosen randomly or informed by η T (t)) 4.
Perform an FFT on a wm,1 (t), to obtain the amplitudes, A a,1 ( f j ), and phase components, φ a,1 ( f j ), for each frequency component of the input source term 5.
Run the simulation, for iteration i, using the wavemaker source term a wm,i (t), and measure the resulting free surface elevation at the chosen NWT location, η R,i (t) 6.
Perform an FFT on η R,i (t), to obtain the amplitudes, A R,i ( f j ), and phase components, φ R,i ( f j ), for each frequency component of the generated wave series 7.
Calculate the new amplitudes for each frequency component of the input source term, A a,i+1 ( f j ), by scaling the previous amplitudes, A a,i ( f j ), with the ratio of target surface elevation amplitude, A T ( f j ), and the generated surface elevation amplitude from the previous run, A R,i ( f j ): Calculate the new phase components, φ a,i+1 ( f j ), by summing φ a,i ( f j ) with the difference between the target elevation phase, φ T ( f j ), and the measured surface elevation phase from the previous run, φ R,i ( f j ): Generate the new time series for the wavemaker source term, a wm,i+1 (t), by performing an inverse Fourier transform on A a,i+1 ( f j ) and φ a,i+1 ( f j ) 10.
Repeat Steps 5-9 until either a maximum number of iterations or a threshold for the Mean-Squared Error (MSE) between the target and resulting surface elevation is reached.
It should be noted that, although the method is based on the work described by Masterton and Swan [27], we have further simplified it. Instead of evaluating the phase-shift caused by the distance between the wavemaker and the point of interest, which would require the use of some wave theory to find the wavenumber k, here, the phase between the wavemaker and the target is directly adjusted (Step 8).

Calibration Procedure for Regular Waves
The calibration method detailed above might fail to resolve the single peak in the frequency domain for short time traces of regular waves. At the same time, in regular waves, the phase of the wave is generally not of interest, allowing for a simplification of the calibration. Schmitt and Elsaesser [28], Schmitt et al. [29], Vorhoelter et al. [30] and Windt et al. [19] used the same solver formulation and simply tuned the amplitude of an oscillating source to create monochromatic waves of a desired height. The following calibration method is used for regular waves: The source term is set to oscillate in the horizontal direction with the desired wave frequency.

2.
The amplitude A a,i is initialised with an arbitrary value. 3.
After the initial, and each subsequent run, the surface elevation is analysed in the time domain. The mean is removed from the surface elevation. The mean wave height, H R,i , is then obtained as the difference between the mean of the positive and mean of the negative peaks.

4.
A new wave maker amplitude, A a,i+1 , is obtained by linearly scaling the previous value with the ratio of target H T and resulting wave height, H R , as follows Steps 3 and 4 are repeated until the desired wave height is achieved.

Case Study
In this section, a case study is presented with two main objectives: firstly, to demonstrate the capabilities of the impulse wavemaker and the self-calibration procedure in producing realistic deep, intermediate or shallow water wave series, at a specified location in a wave tank. The second objective is to provide guidance on the selection of the wavemaker source region, by investigating the effect that the size and position of the source region has on the resulting waves. It might be expected that a very long source area will decrease accuracy, whereas a very short source region will require very large velocity components to create a target wave. For the shallow water impulse wavemaker presented by Choi and Yoon [18], a source length of about a quarter to half a wavelength is recommended.
In shallow water, the entire water column, from the sea floor up to the water surface, performs an oscillating motion, whereas in deep water, only part of the water column is affected, and the wave motion does not extend to the sea floor. It is thus an interesting question how to choose the length and depth of the source region to achieve optimal results.

Target Waves
The case study considered two types of target waves, multi-frequency and regular waves, to demonstrate the calibration procedures outlined in Sections 2.3 and 2.4, respectively.

Multi-Frequency Wave Packet
To demonstrate the ability of the impulse wavemaker, the case study created a unidirectional multi-frequency wave packet at a specified location in the NWT. The wave packet considered was a realisation of the NewWave formulation as presented by Ning et al. [31]. A NewWave wave packet comprises a summation of all the frequency components of a given spectrum, such that the largest amplitude wave crest is created in the temporal centre of the packet. The governing equations for the surface elevation, η T (t), and amplitude of each wave component, a k , with k = 1, ..., K , are given in Equations (4) and (5), respectively: In Equation (4), x 0 represents the spatial focal location and t 0 the temporal focal instant. In Equation (5), S( f ) is the spectral density and ∆ f the frequency step. A 0 represents the amplitude of the largest wave crest.
For this case study, three different wave packets, for deep, intermediate and shallow water conditions, were considered. The wave characteristics for the different waves are listed in Table 1. Although the parameters used in this case study were, to a certain extent, arbitrarily chosen, wave packets are generally well suited for testing the calibration of amplitudes and phases of different frequency components [32] and are also of practical relevance in industrial applications [33]. While for the deep and intermediate water case, the wavelength was kept constant, the water depth was equal for the deep and shallow water case. Choosing the parameters for the shallow water case requires care, in order to avoid breaking waves. Time traces for the surface elevation, η(t), and plots of the spectral density, SDF( f ), are shown in Figure 2, for all three wave packets.

Regular Waves
Two cases demonstrating the creation of regular waves in intermediate and deep water are presented. The target wave height was set to 0.037 m and the period to 1.1 s for the shallow water case and 0.977 s for the deep water case, replicating the peak frequency and maximum wave height of the corresponding wave packets described in Section 3.1.1. If not mentioned specifically, the settings, found to be suitable for the corresponding wave packets, have been used for the setup of these regular wave cases.

Source Region
A rectangular-shaped source region was used for all simulations. Length L, height H and positions were varied to investigate the effect of the source region shape and position on the wave maker performance.

Simulation Platform
The NWT implementation for this case study was based on the interFoam solver from the OpenFOAM toolbox. The interFoam solver uses a VOF approach for modelling the two different fluid phases, air and water, in order to capture and track and the free surface. More details on this solver were presented by Rusche [34]. Although the case study employed OpenFOAM as the CFD solver, the method can easily be applied to any CFD software that allows user coding. For example, Choi and Yoon [18] implemented the shallow water impulse source wavemaker in ANSYS Fluent, using its user-defined functions capability.
The calibration function was implemented in the free scientific programming software, GNU Octave [35]. This results in all of the software utilised in this case study being open-source or free. The source code for the case study setup and implementation has been shared by the authors on the CCP-WSIrepository [36], allowing anyone to access and utilise the developed impulse source wavemaker easily.

Numerical Wave Tank Setup
Since the case study considered a unidirectional input wave, a Two-Dimensional (2D) NWT was implemented to simplify the setup and reduce computational overheads. The NWT setup is described below, detailing the geometry, boundary conditions, mesh and calibration of the absorption beach.

Geometry
The NWT geometry is depicted in Figure 3. The depth of the NWT was set to 0.74 m for the shallow water case, 0.25 m for the intermediate water case and 0.74 m for the deep water case. The target location for the input wave packet was located 2λ downwave from the centre of the wavemaker source region. Two absorption beaches were then located at 1λ upwave and 3.7λ downwave from the source centre.

Boundary Conditions
The front and back faces of the NWT, lying in the x-z plane, were set to empty, indicating a 2D simulation. The bottom, left and right boundaries were set to a wall. The top boundary was set to atmospheric inlet/outlet condition.

Mesh
The mesh is depicted in Figure 3. The mesh was one cell thick in the y direction, for the implementation of the 2D simulation. Mesh refinement has been employed in the interface region leading to a cell size of eight cells per A 0 and 100 cells per wavelength. The spatial discretisation has been chosen based on prior user experience. Furthermore, a grid convergence study has been performed to determine these sizes; based on the case of the wave packet in intermediate water depth. The resulting grid size was consistent with commonly-used discretisations found in the literature [37]. Adjustable time stepping, based on a maximum allowable Courant number of 0.9, was used in all simulations.

Calibration of the Numerical Beach
Perić and Abdel-Maksoud [26] recently published a method based on analytical theory to find the ideal magnitude of the damping parameter sand for monochromatic waves. In the future, it will thus be possible to set the ideal damping parameters prior to each simulation. Herein, a parameter study is performed to investigate the absorption performance for the wave packet. Three different pairs of l beach and sand Max were tested for the intermediate water case, i.e.: l beach = λ p and sand Max = 5; l beach = λ p and sand Max = 7; l beach = 2λ p and sand Max = 6. The first test, using sand Max = 5 s −1 and l beach = λ p , yielded a reflection coefficient of 3.5%, with reflection evaluated using the method described by Mansard and Funke [38]. Increasing sand Max to 7 s −1 resulted in a reflection coefficient of 2.5%; using sand Max = 6 s −1 and l beach = 2λ p reduced the reflection coefficient significantly to 0.2%. The latter beach configuration was then used in all subsequent cases. While even the first configuration, with a resulting reflection coefficient of 3.5%, was better than for example the experimental facilities described in [27] and many numerical methods, it highlights the flexibility of the approach. Larger beaches will invariably come with a greater computational cost, but for cases where very low reflection over a larger frequency range is required, they seem the only viable method. Figure 3 includes the trajectory of the increasing damping parameter, as used in the simulations.

Source Region Shape and Position
To investigate the influence of the size and location of the source region, 120 multi-frequency wave packet experiments were run, for 90 intermediate depth cases and 30 deep water cases, with varying source region layouts. For the intermediate water cases, simulations were run with the source region progressively centred at a third, a half, and two thirds of the water depth. For deep water cases, the centre of the source region was located at a depth of λ p 4 , or half the water depth. For each different source centre location, 30 experiments were run with varying the source length region, L, of 0.125λ p , 0.25λ p , 0.5λ p , 1λ p , 1.25λ p and varying the height of the source region, H, of 0.125d, 0.25d, 0.5d, 0.75d, 1d and 1.25d. For brevity, the parametric study on the source region size and location has only been performed for the deep and intermediate water depth case. Since the results were delivered in a parametrised form, as functions of the wave characteristics, the optimal setting found for the deep and intermediate water depth case can directly be applied to the shallow water case.
The performance of the various source region layouts was evaluated as the MSE between the target surface elevation, η T , and the achieved result, η R : where K is the total number of time steps in the simulation.

Results and Discussion
In this section, the case study results are presented and discussed. First, an example result from a single experiment (deep water case with a source height of d /4 and a source length of λ p/4) is presented in Section 4.1. The complete set of the results from all the deep, intermediate and shallow water experiments is summarised in Sections 4.2-4.4, respectively. Application of the simplified calibration method for regular waves is discussed in Section 4.5. Figure 4 displays an example of the typical decrease in MSE for increasing calibration iterations during an experiment. Note that the evolution is non-monotonic, which is explained by Figure 5, showing the corresponding surface elevation time traces from a selection of these calibration iterations. The first iteration was initialised with random input of small amplitude, and the resulting surface elevation was almost zero, while the MSE of 3 × 10 −5 was large. Obviously, convergence could be improved by choosing better initial values, but even for random noise as start values, this behaviour has been observed as typical for all cases. The second iteration yielded several waves with similar amplitude and frequency to the target wave packet, but out of phase, and the resulting MSE almost doubled. The fourth and fifth iterations decreased the MSE by orders of magnitude to 1 × 10 −5 and 1 × 10 −6 . The main peak of the wave packet was now already correctly resolved, and most of the error stemmed from some small high frequency waves after 11 s. Further iterations decreased the MSE to about 2.8 × 10 −7 , and agreement was now further improved for the small ripples after 11 s. For all the results presented Sections 4.2-4.4, nine calibration iterations were used.

Wave Packet: Deep Water Waves
Deep water waves only affect the water column up to a depth of half a wavelength; it thus seems more appropriate to base the vertical position of the source region on that parameter. Test cases were run for varying source heights and lengths, with the source centre position fixed at λ p /4. Figure 6 shows the results from the deep water experiments; data are also provided in Table A1. The minimum error was 2.8 × 10 −7 for a source length of 0.25λ p and a height of 0.25-times the water depth. A wide range of source heights yielded correct results; only for small values of source heights, where the source region did not reach the surface, did the error increase significantly.

Wave Packet: Intermediate Water Waves
The results from the intermediate water experiments are summarised in Figures 7-9 and also in tabular form in Tables A2-A4. Simulations were run with the centre of the source region at one third, one half and two thirds of the water depth, respectively. The results showed that decreasing the length of the source region, from 1λ p to 0.3λ p or less, reduced the MSE by over two orders of magnitude. In contrast, the height of the source region was seen to have a much smaller influence on the wavemaker performance, with the best results obtained when the source region spanned over the entire water depth. Overall, the smallest MSE occurred for the experiment with the source region centred at half the water depth, the source length of 0.25λ p and a source height of 1.25-times the water depth.  Figure 10 shows the target wave time series and the results of the experiments that produced the best results for the three different source centre depths. Variations between simulations were minimal. The centre of the wave packet was very well reproduced around its peak and the surrounding crests and troughs, with the main discrepancies occurring at the beginning and the end of the wave packet, with the CFD simulations presenting small amplitude ripples when the target wave was already reaching still water conditions.

Wave Packet: Shallow Water Waves
For the shallow water wave experiment, the source region stretched over the whole water column and only different lengths were tested. This seems appropriate considering that depth-independent velocity profiles are a defining feature of shallow water waves. The length of the source was set to λ p/4, λ p/2 and 0.75λ p . Figure 11 shows the resulting MSE values for the three source region lengths, data are also presented in Table A5. The best results were obtained with a length of λ p/2 with an MSE  Figure 12 shows the target and resulting surface elevation time traces after nine calibration iterations for the overall best setting, that is a source length of λ p/2.  Figure 13 shows the surface elevation during the first four calibration steps for deep and shallow water waves. In both cases, the results converged rapidly; after the second iteration, the results were almost identical for subsequent runs. Figure 14 shows the absolute difference between target wave height and the current value. For the second iteration errors were already within the millimetre range and decreased almost another order of magnitude in the third iteration, which can be deemed sufficient for any practical application. 12

Conclusions
A NWM, based on the impulse source method, was implemented in the OpenFOAM framework. In combination with a calibration procedure, complex irregular wave patterns can be recreated in conditions ranging from deep to shallow water using a simple impulse source term acting in the horizontal direction. While the simple formulation of the source term facilitates implementation in the flow solver and reduces computational overheads, the calibration procedure, which is often required in such cases independent of the wavemaker used, ensures that the target wave is created at the desired position and time in the NWT. Furthermore, reflection analysis demonstrates the transparency of the wavemaker region and the ability of the numerical beach to achieve arbitrarily low reflection. In the future, the work of [26] will allow setting the ideal beach parameters without the need to run parameter or calibration studies.
Parameter sensitivity studies, investigating the effect of the shape and position of the wavemaker region, show that correct results can be achieved over a wide range of parameters. A wavemaker region length of less than a quarter of a wavelength is suitable, which is less than recommended in previous work [18], and might enable the use of smaller computational domains.
The vertical position of the centre of the wavemaker has overall little effect in shallow water conditions, and good results were achieved if it was placed at half the water depth. For deep water cases, good results were obtained with the centre of the wavemaker at a quarter of a wavelength below the surface. Even when starting from a poor initial parameter specification, the method was shown to converge to an accurate solution within a few iterations.
For monochromatic or regular waves, a simplified calibration method can be used. Tests showed that the amplitude can typically be found within four iterations.
Although the calibration method is standard in experimental test facilities, it must be expected to become inaccurate in more nonlinear wave conditions. A possible way forward might be the use of neural networks as described by [21].

Conflicts of Interest:
The authors declare no conflict of interest. The founding sponsors had no role in the designof the study; in the collection, analyses, or interpretation of data; in the writing of the manuscript, and in thedecision to publish the results.