# Numerical strategies to speed up CFD computations with free surface

Alban Leroyer, Jeroen Wackers, Patrick Queutey, Emmanuel Guilmineau

► To cite this version:

Alban Leroyer, Jeroen Wackers, Patrick Queutey, Emmanuel Guilmineau. Numerical strategies to speed up CFD computations with free surface: Application to the dynamic equilibrium of hulls. Ocean Engineering, 2011, 10.1016/j.oceaneng.2011.09.006. hal-01145125

HAL Id: hal-01145125

https://hal.science/hal-01145125

Submitted on 26 Apr 2021

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
# Numerical strategies to speed up CFD computations with free surface—Application to the dynamic equilibrium of hulls

Alban Leroyer *, Jeroen Wackers, Patrick Queutey, Emmanuel Guilmineau

Laboratoire de Mécanique des Fluides, Ecole Centrale de Nantes, CNRS UMR 6598, BP 92101, 44321 Nantes Cedex 3, France

This article presents two numerical procedures to speed up computations when dealing with a Reynolds Averaged Navier Stokes (RANS) solver based on the Volume of Fluid (VoF) or multifluid method to treat the free surface. The first one is a time-splitting procedure for the volume fraction equation, enabling the use of larger time steps for the resolution of the flow, without penalizing accuracy. However, these large time steps destabilize the coupling with the ship motion simulation when computing a dynamic equilibrium position in marine applications. The second procedure is therefore a quasi-static approach to solve the coupled problem of dynamic equilibrium. A comparison of these procedures with classical simulations shows that numerical solutions of realistic problems can be obtained up to four times faster.

## 1. Introduction

Free-surface capturing methods, such as the Volume of Fluid (VoF) or level-set formulation, have become more and more popular among the CFD developers involved in viscous marine hydrodynamics. The reason for this increasing interest is that these approaches are more robust than those based on a free-surface fitting methodology since no regridding is necessary. Moreover, the merging or breakup of the interface is also handled in a natural way. To achieve computations with VoF-like approach, specific compressive discretisation schemes are used to solve the volume fraction transport equation and keep the sharpness of the interface (Ubbink, 1997; Queutey and Visonneau, 2007). Even if the capacities and the flexibility of such an approach are unquestionable, two drawbacks can be highlighted:

* The formulation is intrinsically unsteady since the volume fraction is convected by the flow. Up to now, no steady formulation has ever been successful for non-academic test cases, to the knowledge of the authors.

* The compressive property imposes a numerically severe Courant number limitation.

It seems inappropriate and wasteful to use such a complete unsteady approach when dealing with physically steady cases. It is all the more a pity since for an implicit solver, the volume

fraction equation is the only equation to have such a Courant number limitation! This issue was first underlined by Ubbink (1997) in the conclusions of his PhD thesis: "the Courant limitation is not insurmountable because it should be possible to apply a technique of sub-cycling where the time step of the main loop is divided into smaller steps in order to advance the volume fractions...".

In this paper, to reduce this Courant number limitation, such an original time-splitting (also called time subcycling) procedure for the volume fraction equation is developed and validated for steady state cases. For instance, it enables to increase the global time step while keeping the Courant number constant. The procedure creates no problem for fixed bodies. However, when coupling the flow solver with a resolution of Newton's laws for the ship body, in order to reach the dynamic equilibrium position of a hull, the large time steps lead to a divergent flow-motion coupling, due to the added mass effects (Söding, 2001). To avoid this problem, a quasi-static approach has been implemented to govern the motion of the body in time. This technique has been successfully combined with the time-splitting procedure for the volume fraction, removing the flow/motion instabilities.

Compared to a classical unsteady approach using the resolution of Newton's law, this new numerical procedure to deal with steady cases for hydrodynamics applications enables to reduce significantly the CPU time.

After a brief description of the RANS solver in which this work has been implemented (Section 2), this article outlines the time-splitting procedure (Section 3) as well as the quasi-static approach (Section 4). Then, in Section 5 a test-case is shown demonstrating the capability and the efficiency of these techniques.

<footnote>* Corresponding author. Fax: +33 2 40 37 25 23.
E-mail address: alban.leroyer@ec-nantes.fr (A. Leroyer).</footnote>
## 2. The ISIS-CFD flow solver

ISIS-CFD is developed by the EMN group (Equipe Modélisation Numérique) of the Fluid Mechanics Laboratory of the Ecole Centrale Nantes. The solver is based on the Finite-Volume method to build a second-order accurate discretisation of the RANS equations. The flow equations are constructed face by face, which means that cells having an arbitrary number of arbitrarily shaped faces can be accepted. This enables the simulation of flows around complex geometries. The velocity field is obtained from the momentum conservation equations and the pressure field is extracted from the mass conservation constraint, or continuity equation, transformed into a pressure equation. In the case of turbulent flows, additional transport equations for modelled turbulent parameters are solved in a form similar to the momentum equations: they are discretised and solved using the same principles.

Free-surface flows are computed through an interface capturing method: the flow phases are modelled with a transport equation for the volume fraction of water $c_w$ in each cell: $c_w=1$ means that the cell is completely filled with water, $c_w=0$ means that only air is present in the considered cell. The interface between air and water is represented by the numerical discontinuity in the solution for $c_w$. The effective flow physical properties (dynamic viscosity $\mu$ and density $\rho$) are obtained from the properties for each phase, $(\mu_w, \mu_a)$ and $(\rho_w, \rho_a)$, respectively, for water and air, with the following constitutive relations: $\mu = c_w\mu_w + (1-c_w)\mu_a$ and $\rho = c_w\rho_w + (1-c_w)\rho_a$.

Special attention has to be paid to preserve the sharpness of the interface when solving the transport equation of $c_w$. Therefore, this equation is discretised with a scheme using specific anti-diffusive properties, which are fully activated when the amount of fluid drained out of a cell during one time step does not exceed the cell volume. This implies a constraint in the local Courant number (the Courant number is an adimensional parameter roughly defined by: $\Delta tV/\Delta x$, where V is the velocity through the considered cell, $\Delta x$ is the size of the cell, and $\Delta t$ is the global time step of the temporal discretisation). The latter should not exceed 1 and the usual target value retained is 0.3 for 3D flows. Here, the Blended Reconstructed Interface Capturing Scheme (BRICS) scheme is used (Wackers et al., 2011).

An ALE (Arbitrary Lagrangian Eulerian) approach is used to deal with moving bodies (Hirt et al., 1974; Donea et al., 2004; Leroyer, 2004). All configurations of motion (up to six solved or imposed DOF) can be applied. Analytical weighted deformation techniques have been developed to preserve a mesh fitted to the body during its motion (Leroyer, 2004; Leroyer and Visonneau, 2005).

## 3. The time-splitting procedure for the volume fraction

### 3.1. Description

With an interface capturing approach, we have to solve the ALE convection equation for the volume fraction ($c_w$ is simply denoted by $c$ in the following):

$$ \frac{\delta}{\delta t} \int_V c \, dV + \oint_S c(\vec{U} - \vec{U}_d) \cdot \vec{n} \, dS = 0, \quad (1) $$

where $V$ is the domain of interest, or control volume, bounded by the closed surface $S$ moving at the velocity $\vec{U}_d$ with a unit normal vector $\vec{n}$ directed outward. The time derivative following the moving grid is written $\delta/\delta t$.

Finite-Volume discretisation leads to the following discretised form:

$$ \frac{c(t_c)V(t_c)-c(t_p)V(t_p)}{\Delta t} + \sum_{\text{faces } S_f} c_f(\mathcal{F}\vec{U} - \mathcal{F}\vec{U}_d) = 0, \quad (2) $$

where $t_p$ and $t_c$ mean, respectively, the previous and the current instant, and $\Delta t = t_c - t_p$ is the current time step. $c_f$ denotes the reconstruction of the volume fraction at the centre of the face, whereas $\mathcal{F}\vec{U}$ and $\mathcal{F}\vec{U}_d$ represent, respectively, the velocity flux and the grid displacement velocity flux through the considered face $S_f$. Since steady configurations are investigated here, accurate resolution in time is not required, so the basic first order Euler implicit scheme is applied for the time derivative.

The spirit of the time-splitting approach is to enhance the fulfilment of the Courant-Friedrichs-Lewy (CFL) condition related to Eq. (2), by using a specific time step for the volume fraction, which is a fraction of the time step associated with the global simulation. In other terms, the global time step $\Delta t$ is split into a sequence of smaller ones leading naturally to smaller Courant numbers (Fig. 1). As a consequence, the volume fraction equation is solved sequentially several times during a single global time step, for only one update of the other flow equations; each resolution is called a "subcycle". If we note $N$ the number of "subcycles", the split time step $\Delta t_i$ is equal to $\Delta t/N$. The intermediate volumes $V^i$, for which no mesh is reconstructed, are linearly interpolated between $t_p$ and $t_c$. Thus, the convection equation for $c$ becomes

$$ \frac{c^i V^i - c^{i-1} V^{i-1}}{\Delta t_i} + \sum_{\text{faces } S_f} c_f^i (\mathcal{F}\vec{U} - \mathcal{F}\vec{U}_d) = 0. \quad (3) $$

Instead of solving Eq. (2), we solve $N$ times Eq. (3), going forward in time progressively. From one equation with a typical Courant number $C_0$ for the classical approach, the time-splitting approach leads to $N$ equations to solve but with typical Courant numbers around $C_0/N$. Finally, by summing Eq. (3) from $i$ equal 1 to $N$, we obtain the following Eq. (4), similar to Eq. (2):

$$ \frac{c(t_c)V(t_c)-c(t_p)V(t_p)}{\Delta t} + \sum_{\text{faces } S_f} \tilde{c}_f(\mathcal{F}\vec{U} - \mathcal{F}\vec{U}_d) = 0 \quad \text{with} \quad \tilde{c}_f = \sum_{i=1}^{N} \frac{c_f^i}{N}. \quad (4) $$

The free-surface moves little by little during these small time steps. As the CPU time related to the resolution of the volume fraction equation is not large compared with other parts of the solver (especially the computation of the pressure), the global CPU time of the simulation is strongly reduced. When the steady state of the flow is reached, the solution is obviously the same as the one obtained with a classical approach, since the temporal derivative term vanishes. Then, solving either once or several times a convection equation with a null temporal derivative does not change anything. With the time splitting approach, the way to reach a steady state remains physically correct (even if we try to remove the transitory state quickly without solving it accurately). This is probably the reason why this procedure is very robust.

<figure><img src="image_2.png" /><figcaption>Fig. 1. Temporal diagram.</figcaption></figure>
## 3.2. Example

To illustrate the time splitting procedure, a 2D configuration of a submerged NACA-0012 hydrofoil is considered. This well-known test-case experimentally investigated by Duncan (1983) is often used for the validation of numerical methods developed to compute free-surface flows. In the experimental setup, the hydrofoil, whose chord is $c=20.3$ cm, is towed in a tank with speed $U=0.8$ m/s, with an angle of attack $\alpha=5^{\circ}$ (see Fig. 2). The relevant non-dimensional parameters based on the chord length and the free-stream velocity are $Fr=0.5672$ and $Re=1.423e5$. The distance between the profile and the bottom of the basin is kept fixed ($H=17.5$ cm from the mid-chord of the profile), whereas the depth of submergence $s$ is varying. Here, we are only interested in the case $s=23.6$ cm.

The simulations are performed using a mesh of 38 000 cells. The two-equation $k-\omega$ SST closure (Menter, 1993) is used to take into account the turbulence phenomena. Boundary conditions are imposed as follows: on the top and on the outlet, a hydrostatic pressure is imposed. On the inlet, the velocity is imposed at its far field value. A slip condition is applied on the bottom of the tank, whereas wall-function boundaries are applied on the whole surface of the foil.

Different simulations are performed with various numbers of subcycles $N$. For all the simulations, the time step is such that the Courant number $Co$ for the split volume fraction equation remains around 0.3. Table 1 summarizes the numerical setup and results in terms of CPU time. As expected, the CPU time is greatly reduced using the time-splitting approach, and the final solutions (forces on the body and free-surface deformation) are nearly identical: Fig. 3 shows that the converged value of the lift is the same for all computations, but the path to convergence is slightly different. The classical computation reaches the mean value quickly but with oscillations, which take a long time to damp out. When the number of subcycles is increased, the physical time to obtain the converged value increases too, but the gain in CPU time remains significant (the CPU time for 10 subcycles is four times shorter than for classical approach). Above $N=10$, the efficiency of the

approach saturates but robustness is kept even if a large break is visible at $t \approx 9$ s for $N=20$. Fig. 4 confirms that the converged free surface is not influenced by the time-splitting approach.

## 4. The quasi-static approach

When dealing with potential flows, Delhommeau (1987) underlined that the loads on a hull can be divided into four parts: external forces (due to mechanical bindings or mooring lines), gravity force, hydrostatic loads, and hydrodynamic loads. A first-order evaluation of the hydrostatic loads around an equilibrium position can be obtained as a function of the sink $\delta T_z$, the roll $\delta R_x$ and the trim $\delta R_y$, as well as geometric and inertia characteristics

<figure><img src="image_3.png" /><figcaption>Fig. 3. Comparison of lift force in time.</figcaption></figure>

<figure><img src="image_4.png" /><figcaption>Fig. 2. Experimental setup.</figcaption></figure>

<table><caption>**Table 1**

Numerical setup and CPU time (*): with some residual oscillations.</caption><tr><td>Time step (s)</td><td>0.0025</td><td>0.0125</td><td>0.025</td><td>0.05</td></tr><tr><td>Number of subcycles, N</td><td>0</td><td>5</td><td>10</td><td>20</td></tr><tr><td>Duration of the velocity ramp (s)</td><td>1</td><td>2</td><td>3</td><td>4</td></tr><tr><td>Duration of the simulation (s)</td><td>20</td><td>30</td><td>40</td><td>50</td></tr><tr><td>Total CPU time (min)</td><td>258</td><td>112</td><td>103</td><td>102</td></tr><tr><td>CPU time per second of physical time (min)</td><td>13</td><td>3.7</td><td>2.5</td><td>2.06</td></tr><tr><td>Physical time to reach the converged lift (min)</td><td>15 (*)</td><td>17</td><td>18</td><td>32</td></tr><tr><td>CPU time to reach the converged lift (min)</td><td>195 (*)</td><td>63</td><td>45</td><td>66</td></tr></table>
(Delhommeau, 1987). The relation is given with the following equation:

$$ \begin{bmatrix} F_x \\ M_x \\ M_z \end{bmatrix} = \underbrace{\begin{bmatrix} S_{33} & S_{34} & S_{35} \\ S_{43} & S_{44} & S_{45} \\ S_{53} & S_{54} & S_{55} \end{bmatrix}}_{\mathcal{S}} \begin{bmatrix} \delta T_z \\ \delta R_x \\ \delta R_y \end{bmatrix}, \qquad (5) $$

$$ S_{33} = \rho g \int_{S_w} \underbrace{(\vec{n} \cdot \vec{z})}_{dS_z} dS, $$

$$ S_{34} = S_{43} = \rho g \int_{S_w} (Y - Y_G) dS_z, \quad S_{35} = S_{53} = \rho g \int_{S_w} (X - X_G) dS_z, $$

$$ S_{44} = \rho g \int_{S_w} (Y - Y_G)^2 dS_z + V_i(Z_c - Z_G), $$

$$ S_{55} = \rho g \int_{S_w} (X - X_G)^2 dS_z + V_i(Z_c - Z_G), $$

$$ S_{45} = S_{54} = -\rho g \int_{S_w} (X - X_G)(Y - Y_G) dS_z. $$

The coefficients of the matrix $\mathcal{S}$ depend on the density of water ($\rho$), the magnitude of gravity ($g$), the geometry of the wetted surface of the hull ($S_w$), the immersed volume ($V_i$), the centre of gravity (coordinates ($X_G, Y_G, Z_G$)) and the vertical position of the centre of buoyancy ($Z_c$).

By computing $\mathcal{S}^{-1}$, Eq. (6) can be used to provide the hydrostatic position of any hull through an iterative process, assuming

<figure><img src="image_5.png" /><figcaption>Fig. 4. Comparison of free-surface for $N=0$ and $N=10$.</figcaption></figure>

that the mass and the centre of gravity position are known:

$$ \begin{bmatrix} \delta T_z \\ \delta R_x \\ \delta R_y \end{bmatrix} = \mathcal{S}^{-1} \begin{bmatrix} F_x \\ M_x \\ M_z \end{bmatrix}. \qquad (6) $$

Coupled to a flow solver, it enables to predict a dynamic equilibrium position. In fact, let us consider $\mathcal{L} = (F_z, M_x, M_y)$ the vertical force and torque (evaluated at the centre of gravity $G$) along the $X$ and $Y$-axis, respectively, acting on the hull (including gravity, fluid force and possible external forces), in a given imposed position. By computing the coefficients of the matrix $\mathcal{S}$ using this current position of the ship, a prediction of the new position can be obtained by solving the linear system equation (6), giving $\delta \mathcal{P} = (\delta T_z, \delta R_x, \delta R_y)$ as a result.

This amounts to extrapolating an equilibrium position with a first order method while keeping the hydrodynamic load and the external forces constant. This inferred new position is relaxed and then applied progressively during the kinematic transient period $\Delta T_i$. After some time steps, a new prediction of the equilibrium position using the current fluid forces applied to the hull can be obtained. Then the procedure can go on up to convergence. The procedure is summed up as follows ($t_e$ means the time when an evaluation is carried out, $\Delta T_h$ is the temporal interval between two evaluations, with $\Delta T_h \ge \Delta T_i$):

① $t_e=t_c$,

② computation of the current forces $\mathcal{L}$ acting on the hull,

③ prediction of $\delta \mathcal{P}$ to reach the new equilibrium position using Eq. (6),

④ underrelaxation of the new position with a coefficient $\alpha \in [0, 1]$: $\mathcal{P}^{\text{new}} = \mathcal{P}^{\text{old}} + \alpha \delta \mathcal{P}$,

⑤ gradual variation of the position to reach $\mathcal{P}^{\text{new}}$ within $\Delta T_i$,

⑥ when $t_c \ge t_e + \Delta T_h$, go to ①.

## 5. An application test-case: DTMB 5415 with free trim and sinkage

Both the techniques described previously are applied here to the bare hull DTMB5415 test case (see Fig. 5) in the following conditions: $Fr=0.28$ and $Re=1.26e7$, trim and sinkage free, $L_{pp}=3.048$ m. This test case was used for the CFD Workshop Tokyo 2005 (CFDWS2005, case 1.3) and is the subject of experimental and numerical studies (see Longo and Stern, 2005; Carrica et al., 2007 for example). A mesh of 1.6 M cells around the half body is used to compute the flow (see Figs. 6 and 7). Boundary conditions are imposed as follows: on the top and the bottom of the fluid domain, a hydrostatic pressure is imposed, with a Neumann condition for velocity. On the inflow, outflow and side boundary, the velocity is imposed at its far field value, as is the volume fraction (in this case, a Neumann condition is prescribed

<figure><img src="image_6.png" /><figcaption>Fig. 5. The DTMB5415 bare hull.</figcaption></figure>
<figure><img src="image_7.png" /><figcaption>Fig. 6. Fluid domain.</figcaption></figure>

<figure><img src="image_8.png" /><figcaption>Fig. 7. Zoomed view of the mesh at the bow.</figcaption></figure>

for the pressure). Mirror symmetry is used for the plane $Y=0$. On the body, a wall function is used everywhere, except on the deck where a slip condition is applied.

A classical computation Ref (no subcycling, time-accurate resolution of Newton's laws) and a computation TS-QS combining the time-splitting and the quasi-static approaches are performed. For the latter, 10 subcycles are applied, using a time step equal to 0.05 s, i.e. 10 times those of the classical computation. The ship is accelerated from rest to its nominal velocity ($=1.531 \text{ m/s}$) in 2 s and 4 s, respectively, for the classical and quasi-static computation. Acceleration ramp follows a quarter of sine law to reach the nominal velocity smoothly. The convergence time is not very sensitive to the duration of this ramp, within a certain range. However, to avoid convergence troubles, it is recommended not to go far beyond an acceleration of $1 \text{ m/s}^2$ for model scale hulls. For the quasi-static computation, we prefer to double this transient stage since the time step is far larger.

The parameters for the quasi-static approach are as follows: a new prediction of the equilibrium position is computed each $\Delta T_h = 0.3$ s and the relaxation value $\alpha$ of the displacement is set to 0.3. While the sinkage is updated from the beginning, the trim is released only after reaching the nominal velocity, i.e. after 4 s. Linear laws for trim and sinkage are imposed to the ship to reach the new predicted position within a temporal interval of $\Delta T_i$. For this case, $\Delta T_i$ is chosen equal to $\Delta T_h$, i.e. as soon as the ship reaches its new position, a new prediction occurs, without any latency period in which the boat remains fixed. This latency period (i.e. $\Delta T_i < \Delta T_h$) is needed for stiffer cases like planing hulls, to let the hydrodynamic force stabilize before a new prediction of the equilibrium.

Figs. 8 and 9 show the evolution of trim and sinkage for both configurations. It takes more time to reach the equilibrium position for the simulation TS-QS, but since the time step is multiplied by 10, this is not really a problem. However, contrary
<figure><img src="image_9.png" /><figcaption>**Fig. 8.** Trim as a function of time.</figcaption></figure>

<figure><img src="image_10.png" /><figcaption>**Fig. 9.** Sinkage as a function of time.</figcaption></figure>

to the classical approach in which the oscillations of the motion are difficult to damp out quickly (even using an artificial damping term), the quasi-static approach gives a very stable equilibrium position. The values are in quite good agreement with the experimental data: differences are comparable to those obtained in Carrica et al. (2007). Here again, the free-surface is the same for both results (see Fig. 10). The gain in CPU time is similar to that exposed in Section 3.2.

## 6. Conclusions

This article describes first an algorithm related to a time-splitting technique, in which the advance in time is not the same for the volume fraction and for the other variables: the equation of the volume fraction is solved using a smaller time step which reduces significantly the Courant number. A test-case is presented showing that the optimal number of splits is around 10. Then, a quasi-static approach is described to be used instead of the

<figure><img src="image_11.png" /><figcaption>**Fig. 10.** Elevation of free surface. —: positive value, - - -: negative value, $ΔZ/L_{pp} = 0.001$ between two successive lines.</figcaption></figure>

resolution of the Newton's laws to reach dynamic equilibrium. The convergence is carried out through a succession of predicted positions reached with an imposed motion. Compared to the dynamic approach, the transients are easier to remove, the dynamic equilibrium is far less oscillating. As a consequence, a criterion can be added to test automatically several velocities the ones after the others. Moreover, the large time steps associated with the time-splitting technique can then be used without any flow/motion instabilities. As a consequence, both techniques together enable to speed up computations for reaching a dynamic equilibrium position, without penalizing accuracy. A test-case is shown demonstrating that the results are similar to a classical approach, but with a gain of CPU time around 4.

These techniques have been applied on various kinds of hulls. For some of them (like planing boats), for which the resolution of the equilibrium position is more stiff, a latency period is needed between two predictions to let the fluid forces to be stabilized.

This approach can also be used together with an adaptive time step law to the Courant number (Hay et al., 2006). In this case, a target Courant number can be specified for both the split fraction volume equation and the global time step (the ratio of these both values give approximatively the number of subcycles). This avoids thinking about which time step should be used to have enough accuracy for the interface capturing.

The time-splitting method can also be used with a standard small time step to reduce the Courant number and improve the accuracy of the interface capturing. As a matter of fact, in real cases, the criterion of having the Courant number lower than 0.3 is often difficult to reach everywhere.

The next step is to investigate this sub-cycling method for unsteady physical problems, even if the accuracy of such an approach needs to be addressed in this case.

## Acknowledgements

This work was granted access to the HPC resources of IDRIS under the allocation 2009-000129 made by GENCI (Grand Equi-pement National de Calcul Intensif).

## References

* Carrica, P.M., Wilson, R.V., Noack, R.W., Stern, F., 2007. Ship motions using single-phase level set with dynamic overset grids. Comput. Fluids 36, 1415-1433.

* Delhommeau, G., 1987. Les problèmes de diffraction-radiation et de résistance de vagues: étude théorique et résolution numérique par la méthode des singularités. Ph.D. Thesis. Ecole Nationale Supérieure de Mécanique de Nantes, France.
Donea, J., Huerta, A., Ponthot, J., Rodríguez-Ferran, A., 2004. Arbitrary Lagrangian-Eulerian methods. Encyclopedia of Computational Mechanics.

Duncan, J., 1983. The breaking and non-breaking wave resistance of a two-dimensional hydrofoil. J. Fluid Mech. 126, 507-520.

Hay, A., Leroyer, A., Visonneau, M., 2006. H-adaptive Navier-Stokes simulations of free-surface flows around moving bodies. J. Mar. Sci. Technol. 11, 1-18.

Hirt, C., Amsden, A., Cook, J., 1974. An arbitrary Lagrangian-Eulerian computing method for all flow speeds. J. Computat. Phys. 14, 227-253.

Leroyer, A., 2004. Etude du couplage écoulement/mouvement pour des corps solides ou à déformation imposée par résolution des équations de Navier-Stokes. Contribution à la modélisation numérique de la cavitation. Ph.D. Thesis. Ecole Centrale de Nantes, France.

Leroyer, A., Visonneau, M., 2005. Numerical methods for RANSE simulations of a self-propelled fish-like body. J. Fluids Struct. 20, 975-991.

Longo, J., Stern, F., 2005. Uncertainty assessment for towing tank tests with example for surface combatant DTMB model 5415. J. Ship Res. 49, 55-68.

Menter, F., 1993. Zonal two-equation $k-\omega$ turbulence models for aerodynamic flows. In: AIAA 24th Fluid Dynamic Conference, Orlando, FL.

Queutey, P., Visonneau, M., 2007. An interface capturing method for free-surface hydrodynamic flows. Comput. Fluids 36, 1481-1510.

Söding, H., 2001. How to integrate free motions of solids in fluids. In: 4th Numerical Towing Tank Symposium, Hamburg, Germany.

Ubbink, O., 1997. Numerical prediction of two-fluid systems with sharp interfaces. Ph.D. Thesis. Imperial College, University of London.

Wackers, J., Koren, B., Raven, H., Ploeg, A., Starke, A., Deng, G., Queutey, P., Visonneau, M., Hino, T., Ohashi, K., 2011. Free-surface viscous flow solution methods for ship hydrodynamics. Arch. Comput. Methods Eng. 18, 1-41.