Published in Encyclopedia of Computational Mechanics, Encyclopedia of Computational Mechanics, E. Stein, R. de Borst and T.J.R. Hughes (Eds.), John Wiley & Sons Ltd, Vol. 3, Chapter 18, pp. 579  607, 2004
DOI: 10.1002/9781119176817
This chapter presents an overview of some computational methods for analysis of ship hydrodynamics problems. Attention is focused on the description of a stabilized finite element formulation derived via a finite calculus procedure. Both arbitrary LagrangianEulerian (ALE) and fully Lagrangian forms are presented. Details of the treatment of the freesurface waves and the interaction between the ship structure and the sea water are also given. Examples of application to a variety of ship hydrodynamics problems are shown.
Keywords Ship hydrodynamics, finite element method, finite calculus, fluidstructure interaction
Accurate prediction of the hydrodynamic forces on a ship in motion is of paramount importance in ship design. The water resistance at a certain speed determines the required engine power and thereby the fuel consumption. Minimization of the hydrodynamic forces is therefore an important issue in ship hull design. Further, excitation of a wave pattern by ship motion not only induces wave resistance but may also limit the speed in the vicinity of the shore for environmental reasons, which must also be taken into account in ship design.
The usual simplification in ship hydrodynamics design is to separately consider the performance of the ship in still water and its behaviour in open sea. Hydrodynamic optimisation of a hull primarily requires the calculation of the resistance in a calm sea and the open sea effects are generally taken into account as a wave added resistance.
The resistance of a ship in still water can be considered as the sum of several contributions: a viscous resistance associated with the generation of boundary layers, the wave resistance, the air resistance on the superstructure and the induced resistance related to the generation of lift forces.
Wave resistance in practical cases amounts to to of the total resistance of a ship in still water (Raven, 1996 [84]). It increases very rapidly at high speeds dominating the viscous component for high speed ships. Furthermore, wave resistance is very sensitive to the hull form design and easily affected by small shape modifications. For all these reasons, the possibility to predict and reduce the wave resistance is an important target.
The prediction of the wave pattern and the wave resistance of a ship has challenged mathematicians and hydrodynamicists for over a century. The Boundary Element Method (BEM) is the basis of many computational algorithms developed in past years. Here the flow problem is solved using a simple potential model. BEM methods, termed by hydrodynamicists as Panel Methods may be classified into two categories. The first one uses the Kelvin wave source as the elementary singularity. The main advantage of such scheme is the automatic satisfaction of the radiation condition. The theoretical background of this method was reviewed by Wehausen (1970) [103], while computational aspects can be found in Soding (1996) [87] and Jenson and Soding (1989) [58]. The second class of BEM schemes uses the Rankine source as the elementary singularity. This procedure, first presented by Dawson (1977) [29], has been widely applied in practice and many improvements have been addressed to account for the nonlinear wave effects. Among these, a succesful example is the Rankine Panel Method (Xia, 1986 [106]; Jenson and Soding, 1989 [58]; Nakos and Sclavounos, 1990 [109]).
In addition to the important developments in potential flow panel methods for practical ship hydrodynamics analysis during the period 19601980, much research in the second half of the twentieth century was oriented towards the introduction of viscosity in the CFD analysis. In the 1960's the viscous flow research was mainly focused in 2D boundary layer theory and by the end of the decade several methods for arbitrary pressure gradients were available. This research continued to solve the 3D case during the following decade and an evaluation of the capability of the new methods to predict ship wave resistance was carried out at different workshops (Bai and McCarthy, 1979 [7]; Larsson, 1981 [61]; Noblesse and McCarthy, 1983 [69]). Here application to some well specified test cases were reported and numerical and experimental results compared acceptable well for most part of the boundary layer along the hull, while wrong results were obtained near the stern. This prompted additional research and by the end of the 1980's a number of numerical procedures for solving the full viscous flow equation accounting for simple turbulence modes based on Reynolds averaged NavierStokes (RANS) equations were available. Considerable improvements for predicting the stern flow were reported in subsequent workshops organized in the 1990's (Kim and Lucas, 1990 [59]; Reed et al. , 1990[85]; Beck et al., 1993 [8]; Raven, 1996 [84]; Soding, 1996 [87]; Janson and Larsson, 1996 [57]; Alessandrini and Delhommeau, 1996 [1]; Miyata, 1996 [67], Löhner et al., 1998 [64]). A good review of the status of CFD in ship hydrodynamics in the last part of the 20th century can be found in Larsson et al. (1998) [62].
Independently of the flow equations used, the freesurface boundary condition has been solved in different manners. The exact free surface condition is nonlinear and several linearizations have been proposed (Baba and Takekuma, 1975 [6]; Newmann, 1976 [68]; Idelsohn et al., 1999 [52]). Some of them use a fixed domain and others a moving one. An alternative is to solve the full nonlinear free surface equation on a reference surface which does not necessarily coincides with the free surface itself. In this way the updating of the surface mesh is minimized and sometimes is not even necessary.
The solution of the freesurface equation in a bounded domain brings in the necessity of a radiation condition to eliminate spurious waves. A way to introduce this condition was proposed by Dawson (1977) [29] who used a finite difference (FD) formula based in four upwind points to evaluate the first derivatives appearing in the freesurface equation. This method became very popular and this is probably the main reason why a large majority of codes predicting the wave resistance of ships use FD methods on structured meshes (Larsson et al., 1998 [62]).
Indeed the 1990's were a decade of considerable progress in CFD methods for ship hydrodynamics and the most important breakthrough was perhaps the coupled solution of the freesurface equation with the fluid flow equations. Here a number of viscous and inviscid solutions for the surface ship wave problem using finite element (FE) and finite volume (FV) methods with non structured grids were reported (Farmer et al., 1993 [35], Hino et al., 1993 [41]; Luo et al., 1995 [66]; Storti et al., 1998a,b [9091]; García et al., 1998 [36]; García, 1999 [37]; Alessandrini and Delhommeau, 1999 [2]; Idelsohn et al., 1999 [52]; Löhner et al., 1998 [64], 1999 [65]).
The current challenges in CFD research for ship hydrodynamics focus in the development of robust (stable) and computationaly efficient numerical methods able to capture the different scales involved in the analysis of practical ship hydrodynamics situations. Wave resistance coefficients for modern ship design are needed for a wide range of speeds and here the accurate prediction of the wave pattern and the hull pressure distribution at low speed (say below Froude number () ) are still major challenges. Great difficulties also exist in the computation of the viscous resistance which requires very fine grids in the vicinity of the hull, resulting in overall meshes involving (at least) some elements. Other relevant problems are the prediction of the wake details and the propellerhull interaction. Fine meshing and advanced turbulence models are crucial for the realistic solution of these problems. Indeed the use of unstructured meshes is essential for problems involving complex shapes.
A different class of ship hydrodynamic problems require the modelling of breaking waves or the prediction of water inside the hull (green water) due to large amplitude waves typical of sea keeping problems. Here Lagrangian flow methods where the motion of each flow particle is individually tracked using techniques developed for (incompressible) solid mechanics are a promising trend for solving a wide class of ship hydrodynamics problems.
The content of the chapter is structured as follows. In the next section the standard NavierStokes equations for an incompressible viscous flow are presented. The equations are formulated in an arbitrary LagrangianEulerian (ALE) description allowing the independent motion of the mesh nodes from that of the fluid particles. Details of the problems posed by the freesurface wave boundary condition are given. The difficulties encountered in the numerical solution of the fluid flow and the freesurface equations, namely the unstabilities induced by the convective terms and the limits in the approximation introduced by the incompressibility constraint are explained. A new procedure for deriving stabilized numerical methods for this type of problems based on the so called finite calculus (FIC) formulation is presented. The FIC method is based in redefining the standard governing equations in fluid mechanics by expressing the balance laws in a domain of finite size. This introduces additional terms in the differential equations of the infinitessimal theory which are essential to derive stabilized numerical schemes (Oñate, 1998 [70], 2000 [71], 2004 [72]). We present here a stabilized finite element method using equalorder linear interpolation for the velocity and the pressure variables. Both monolithic and fractional step time integration procedures are described. A method for solving the coupled fluidstructure interaction problem induced by the motion of the ship due to the hydrodynamic forces is also presented. A mesh moving algorithm for updating the position of the freesurface nodes during the ship motion is given.
In the last part of the chapter the Lagrangian formulation for fluid flow analysis is presented as a particular case of the ALE form. As above mentioned the Lagrangian description has particular advantages for tracking the displacement of the fluid particles in flows where large motions of the fluid surface occur. One of the advantages of the Lagrangian approach is that the convective terms dissapear in the governing equations of the fluid. In return, the updating of the mesh at almost every time step is now a necessity and efficient mesh generation algorithms must be used (Idelsohn et al., 2002 [53], 2003a,b,c [5456]).
The examples show the efficiency of the ALE and fully Lagrangian formulations to solve a variety of ship hydrodynamics problems.
The NavierStokes equations for an incompressible fluid in a domain can be written in an arbitrary LagrangianEulerian form as
Momentum

(1) 
Mass conservation

(2) 
In Eq.(1) is the velocity along the ith global reference axis, is the velocity of the moving mesh nodes, is the relative velocity between the fluid and the moving mesh nodes, is the (constant) density of the fluid, are body forces, is the time, is the pressure and are the viscous stresses related to the viscosity by the standard expression

(3) 
where is the Kronecker delta and the strain rates are

(4) 
The volumetric strain rate terms in the Eqs.(1) and (3) can be neglected following Eq.(2). We will however retain these terms as they are useful for the derivation of the stabilized formulation in Section 5.
Eqs.(1) and (3) are completed with the boundary conditions

(5a) 

(5b) 
where are the total stresses, are the component of the unit normal vector to the boundary and and are prescribed tractions and velocities on the Neumann and Dirichlet boundaries and , respectively where is the boundary of the analysis domain . The boundary conditions are completed with the initial condition for .
In above equations where is the number of space dimension (i.e. for 3D). Also, throughout the text the summation convention for repeated indexes is assumed unless specified otherwise.
The boundary conditions (5a) on the surface tractions can be written in local normal and tangential axes as

(6) 
where and are the normal and tangential stresses, respectively and and are the prescribed normal and tangential tractions, respectively.
On the free boundary we have to ensure at all times that: 1) the pressure (which approximates the normal stres) equals the atmospheric pressure and the tangential tractions are zero unless specified otherwise, and 2) material particles of the fluid belong to the freesurface.
The condition on the pressure is simply especified as

(7) 
where is the atmospheric pressure (usually given a zero value).
The condition on the tangential tractions is satisfied by setting in Eq.(6). This is automatically accounted for by the natural boundary conditions in the weak form of the momentum equations (Zienkiewicz and Taylor Vol. 3, 2000 [107]).
The condition on the material particles is expressed (for steady state conditions) as

(8) 
i.e. the velocity vector is tangent to the freesurface. An alternative form of Eq.(8) can be found by noting that the normal vector has the following components

(9) 
In Eq.(9) is the freesurface elevation measured in the direction of the vertical coordinate relative to some previously known surface, which we shall refer to as the reference surface (Figure 2). This surface may be horizontal (i.e. the undisturbed water surface) or may be simply a previously computed surface.
Introducing Eq.(9) into (8) gives (for 3D)

(10) 
Eq.(10) is generalized for the transient 3D case as

(11) 
We observe that obeys a pure convection equation with playing the role of a (non linear) source term. The solution of Eq.(11) with standard Galerkin FEM, or centred FD or FV methods will therefore suffer from numerical instabilities and some kind of stabilization is needed in order to obtain a physically meaningful solution. A method to solve Eq.(11), very popular in the context of the potential flow formulation, was introduced by Dawson (1977) [29] using a four point FD upwind operator to evaluate the first derivatives of on regular grids. Dawson's method has been extended by different authors to solve a many ship hydrodynamics problems (Raven, 1996 [84]; Larsson et al., 1998 [62]; Idelsohn et al., 1999 [52]).
Solution of Eq.(11) is strongly coupled with that of the fluid flow equations. The solution of the whole problem is highly non linear due to the pressence of the unknown velocities in Eq.(11) and also to the fact that the freesurface position defining the new boundary conditions is also unknown. A number of iterative schemes have been developed for the solution of the non linear surface wave problem (Idelsohn et al., 1999 [52]). They all basically involve solving Eq.(11) for the new freesurface height , for fixed values of the velocity field computed from the fluid solver in a previous iteration within each time increment. At this stage two procedures are possible, either the position of the freesurface is updated after each iteration and this becomes the new reference surface, or else an equivalent pressure of value , where is the gravity constant, is applied at the current reference surface as a boundary condition in the next flow iteration. The first option might require the regeneration of a new mesh, whereas the second one is less accurate but computationaly cheaper. Hence a compromise between the two alternatives is usually chosen in practice. The iterative process continues until a converged solution is found for the velocity, the pressure and the freesurface height at each time step. Details of the computational process are described in a next section.
An alternative method to treat the freesurface equation is based in the volume of fluid (VOF) technique (Hirt and Nichols, 1981 [43]). In the VOF method the freesurface position is defined as the interface between two fluids interacting with each other, where the efect of one fluid on the other is very small (i.e. the water and the surrounding air). An interface function which takes the values 0 and 1 for each of the two fluids is transported with the fluid velocity using a time dependent advection equation. Early applications of the VOF in the context of the stabilized FEM were reported by Codina et al., 1994 [26] and Tezduyar et al., 1998 [102]. Examples of application of the VOF to ship hydrodynamics problems can be found in Tajima and Yabe, 1999 [93] and Azcueta et al., 1999 [5].
The development of efficient and robust numerical methods for incompressible flow problems has been a subject of intensive research in last decades. Much effort has been spent in developing the so called stabilized numerical methods overcoming the two main sources of instability in incompressible flow analysis, namely those originated by the high values of the convective terms and those induced by the difficulty in satisfying the incompressibility constraint.
The solution of above problems in the context of the finite element method (FEM) has been attempted in a number of ways. The underdiffusive character of the Galerkin FEM for high convection flows (which incidentaly also occurs for centred FD and FV methods) has been corrected by adding some kind of artificial viscosity terms to the standard Galerkin equations. A good review of such approach can be found in Zienkiewicz and Taylor, Vol. 3 2000 [107] and Donea and Huerta, 2003 [31].
A popular way to overcome the problems with the incompressibility constraint is by introducing a pseudocompressibility in the flow and using implicit and explicit algorithms developed for this kind of problems, such as artificial compressibility schemes (Chorin, 1967A [14]) popular way to overcome the problems with the incompressibility constraint is by introducing a pseudocompressibility in the flow and using implicit and explicit algorithms developed for this kind of problems, such as artificial compressibility schemes (Chorin, 1967 [14]; Farmer et al., 1993 [35]; Peraire et al., 1994 [82]; Briley et al., 1995 [10]; Sheng et al., 1996 [86]) and preconditioning techniques (Idelsohn et al., 1995 [51]). Other FEM schemes with good stabilization properties for the convective and incompressibility terms are based in PetrovGalerkin (PG) techniques. The background of PG methods are the noncentred (upwind) schemes for computing the first derivatives of the convective operator in FD and FV methods. More recently a general class of Galerkin FEM has been developed where the standard Galerkin variational form is extended with adequate residualbased terms in order to achieve a stabilized numerical scheme (Codina, 1998 [16], 2000 [17]). Among the many FEM of this kind we can name the Streamline Upwind Petrov Galerkin (SUPG) method (Hughes and Brooks, 1979 [45]; Brooks and Hughes, 1982 [11]; Tezduyar and Hughes, 1983 [97]; Hughes and Tezduyar, 1984 [46]; Hughes and Mallet, 1986 [48]; Idelsohn et al., 1995 [51]; Storti et al., 1995 [88], 1997 [89]; Cruchaga and Oñate, 1997 [27], 1999 [28]), the Galerkin Least Square (GLS) method (Hughes et al., 1989 [50]; Tezduyar, 1991 [94]; Tezduyar et al., 1992a [99]), the TaylorGalerkin method (Donea, 1984 [30]), the Characteristic Galerkin method (Douglas and Russell, 1982 [32]; Pironneau, 1982 [83]; Löhner et al., 1984 [63]) and its variant the characteristic Based Split (CBS) method (Zienkiewicz and Codina, 1995 [108]; Codina et al., 1998 [23]; Codina and Zienkiewicz, 2002 [25]), pressure gradient operator methods (Codina and Blasco, 1997 [22], 2000 [24]) and the Subgrid Scale (SGS) method (Hughes, 1995 [44]; Brezzi et al., 1997 [9]; Codina, 2000 [17], 2002 [21]).
In this work a stabilized FEM for incompressible flows is derived taking as the starting point the modified governing equations of the flow problem formulated via a finite calculus (FIC) approach. The FIC method is based in invoking the balance of fluxes in a domain of finite size. This introduces naturally additional terms in the classical differential equations of infinitessimal fluid mechanics which are a function of the balance domain dimensions. The new terms in the modified governing equations provide naturally the necessary stabilization to the standard Galerkin finite element method.
In the next section, the main concepts of the FIC approch are introduced via a simple 1D convectiondiffusion model problem. Then the FIC formulation of the fluid flow equations and the freesurface wave equations using the finite element method are presented.
We will consider a convectiondiffusion problem in a 1D domain of length . The equation of balance of fluxes in a subdomain of size belonging to (Figure 1) is written as

(12) 
where and are the incoming and outgoing fluxes at points and , respectively. The flux includes both convective and diffusive terms; i.e. , where is the transported variable, is the velocity and is the diffusitivity of the material.
We express now the fluxes and in terms of the flux at an arbitrary point within the balance domain (Figure 1). Expanding and in Taylor series around point up to second order terms gives

(13) 
Substituting Eq.(13) into Eq.(12) gives after simplification

(14) 
where and all derivatives are computed at point .
Standard calculus theory assumes that the domain is of infinitessimal size and the resulting balance equation is simply . We will relax this assumption and allow the balance domain to have a finite size. The new balance equation (14) incorporates now the underlined term which introduces the characteristic length . Obviously, accounting for higher order terms in Eq.(13) would lead to new terms in Eq.(14) involving higher powers of .
Distance in Eq.(14) can be interpreted as a free parameter depending on the location of point (note that for ). However, the fact that Eq.(14) is the exact balance equation (up to second order terms) for any 1D domain of finite size and that the position of point is arbitrary, can be used to derive numerical schemes with enhanced properties, simply by computing the characteristic length parameter using an adequate “optimality” rule.
Consider, for instance, the modified equation (14) applied to the convectiondiffusion problem. Neglecting third order derivatives of , Eq.(14) can be written as

(15) 
We see that the FIC procedure introduces naturally an additional diffusion term into the standard convectiondiffusion equation. This is the basis of the popular “artificial diffusion” method (Hirsch, 1990 [42]). The characteristic length is typically expressed as a function of the cell or element dimensions. The optimal or critical value of can be computed from numerical stability conditions such as obtaining a physically meaningful solution, or even obtaining “exact” nodal values (Zienkiewicz and Taylor, 2000 [107]; Oñate and Manzan, 1999 [76], 2000 [77]; Oñate, 2004 [72]).
Equation (13) can be extended to account for source and time effects. The full FIC equation for the transient convectiondiffusion problem can be written in compact form as

(16) 

(17) 
where is the external source and is a time stabilization parameter (Oñate, 1998 [70]; Oñate and Manzan, 1999 [76]). For consistency a FIC form of the Neumann boundary condition should be used. This is obtained by invoking balance of fluxes in a domain of finite size next to the boundary where the flux is prescribed to a value . The modified FIC boundary condition is

(18) 
The definition of the problem is completed with the standard Dirichlet condition prescribing the value of at the boundary and the initial conditions.
The starting point are the FIC equations for a viscous incompressible fluid. For simplicity we will neglect the time stabilization term, as this is not relevant for the purposes of this work. The equations are written as (Oñate, 1998 [70], 2000 [71]; Oñate et al., 2002 [79])
Momentum

(19) 
Mass balance

(20) 
where

(21) 
and all the terms have been defined in Section 2.1.
The Neumann boundary conditions for the FIC formulation are (Oñate, 1998 [70], 2000 [71])

(22) 
The Dirichlet and initial boundary conditions are the standard ones given in Section 2.1.
The in above equations are characteristic lengths of the domain where balance of momentum and mass is enforced. In Eq.(22) these lengths define the domain where equilibrium of boundary tractions is established. The sign in front the term in Eq.(22) is consistent with the definition of in Eq.(21).
Eqs.(19)(22) are the starting point for deriving stabilized FEM for solving the incompressible NavierStokes equations using equalorder interpolation for all variables.
From Eqs.(19) and (3) it can be obtained

(23) 
where

(24) 
Substituting Eq.(23) into Eq.(20) leads to the following stabilized mass balance equation

(25) 
The 's in Eq.(23) are intrinsic time parameters per unit mass. Similar parameters also appear in other stabilized formulations (Hughes and Mallet, 1986a [48]; Tezduyar, 1991 [94], 2001 [95]; Codina, 2002 [21]). Note that the parameters of Eq.(24) emerge naturally form the FIC formation and take the values of and for the viscous limit (Stokes flow) and the inviscid limit (Euler flow), respectively.
The weighted residual form of the governing equations (Eqs.(19), (22) and (25)) is

(26) 

(27) 
where and are arbitrary weighting functions representing virtual velocity and virtual pressure fields, respectively. Integrating by parts above equations leads to

(28) 

(29) 
The fourth integral in Eq.(28) is computed as a sum of the element contributions to allow for discontinuities in the derivatives of along the element interfaces. As usual . In Eq.(28) we have neglected the volumetric strain rate term, whereas in the derivation of Eq.(29) we have assumed that is negligible on the boundaries.
The convective and pressure gradient projections and are defined as

(30) 

(31) 
We can now express in Eqs.(28) and (29) in terms of and , respectively which become additional variables. The system of integral equations is now augmented in the necessary number of equations by imposing that the residuals vanish (in average sense) for both forms given by Eqs.(30) and (31). The final system of integral equations is:

(32) 
where and are appropriate weighting functions and the and terms are introduced in the last two equations for convenience. As usual .
We will derive next the FIC equation for the water wave surface.
Let us consider a 2D freesurface wave problem. Figure 3 shows a typical freesurface segment line . The average vertical velocity for the segment is defined as

(33a) 
where and are the vertical velocities of the end points of the segment and .
The average vertical velocity can be computed from the wave heights at points and as

(33b) 
where is the time which a material particle takes to travel from point to point at the average velocity and is the freesurface wave height. Equaling Eq.(33a) and (33b) gives

(33c) 
We can now express the vertical velocities and the wave height at points and in terms of values at an arbitrary point as

(34) 
where all the derivatives are computed at the arbitrary point .
Substituting Eqs.(34) into (33c) and noting that , with , and (Figure 3) and that the position of point is arbitrary, gives the FIC equation for the freesurface height (neglecting high order terms) as

(35) 
with

(36) 
where and are space and time stabilization parameters. The standard infinitesimal form of the freesurface wave condition is obtained by making in Eq.(35) giving

(37) 
which coincides with equation (11) for the 2D case.
A simpler FIC expression can be derived from Eq.(35) by retaining the second order space term only. This gives

(38) 
This can be interpreted as the addition of an artificial diffusion term where plays the role of the new balancing diffusion coefficient.
In the following we will use the 3D ALE form of Eq.(35) neglecting the time stabilization term, given by

(39) 
where, as usual, are the relative velocities between the fluid and the moving mesh nodes.
We will choose continuous linear interpolations of the velocities, the pressure, the convection projections and the pressure gradient projections over three node triangles (2D) and four node tetrahedra (3D). The interpolations are written as

(40) 
where (4) for triangles (tetrahedra), denotes nodal variables and are the linear shape functions (Zienkiewicz and Taylor, Vol 1 2000 [107]).
Substituting the approximations (40) into Eqs.(32) and choosing the Galerking form with leads to following system of discretized equations

(41a) 

(41b) 

(41c) 

(41d) 
The element contributions are given by (for 2D problems)

with and .
In above is the standard strain rate matrix and the deviatoric constitutive matrix (for ). For 2D problems

(43) 
Note that the stabilization matrix brings in an additional orthotropic diffusivity of value . Matrices and are dependent on the velocity field. The solution process can be advanced in time in a (quasinearly) implicit iterative manner using the following scheme.
Step 1

(44) 
Step 2

(45) 
Step 3

(46) 
Step 4

(47) 
In above are time integration parameters with and denotes nodal values at the th time step and the ith iteration. etc. Also for the computations in step 1 at the onset of the iterations.
Steps 1, 3 and 4 can be solved explicitely by choosing a lumped (diagonal) form of matrices and . In this manner the main computational cost is the solution of step 2 involving the inverse of a Laplacian matrix. This can be solved very effectively using an iterative method.
For the iterative proces is unavoidable. The iterations follow until convergence is reached. This can be measured using an adequate error norm in terms of the velocity and pressure variables, or the residuals. Indeed some ot the 's can be made equal to zero. Note that for the algorithm is inconditionable unstable. A simple form is obtained making . This elliminates the non linear dependence with the velocity of matrices and during the iterative scheme.
Convergence of above solution scheme is however difficult for some problems. An enhanced version of the algorithm can be obtained by adding the term where to the equation for the computation of the pressure in the second step. The new term acts as a preconditioner of the pressure equation given now by

(48) 
Note that the added term vanishes for the converged solution (i.e. when ).
An alternative to above algorithm is to use the fractional step method described in the nex section.
The pressure can be split from the discretized momentum equations (see Eq.(44)) as

(49) 

(50) 
In above equations is a variable taking values equal to zero or one. For , (first order scheme) and for , (second order scheme) (Codina, 2001 [19]). Note that in both cases the sum of Eqs.(49) and (50) gives the time discretization of the momentum equations with the pressures computed at . The value of from Eq.(50) is substituted now into (41b) to give

(51a) 
The product can be approximated by a laplacian matrix, i.e.

(51b) 
A semiimplicit algorithm can be derived as follows.
Step 1 Compute the nodal fractional velocities explicitely from Eq.(49) with where subscript denotes a diagonal matrix.
Step 2 Compute from Eq.(51a) (using Eq.(51b)) as

(52) 
Step 3 Compute the nodal velocities explicitely from Eq.(50) with
Step 4 Compute explicitely from Eq.(46) using .
Step 5 Compute explicitely from Eq.(47) with .
This algorithm has an additional step than the iterative algorithm of Section 6.1. The advantage is that now Steps 1 and 2 can be fully linearized by choosing . Also the equation for the pressure variables in Step 2 has improved stabilization properties due to the additional laplacian matrix .
Details on the stability properties of the FIC formulation for incompressible fluid flow problems can be found in Oñate et al. (2002) [79].
The solution in time of Eq.(39) can be written in terms of the nodal velocities computed from the flow solution, as

(53) 
Eq.(53) can now be discretized in space using the standard Galerkin method and solved explicitely for the nodal wave heights at . Typically the general algorithm will be as follows:
The mesh updating process can also include the freesurface nodes, although this is not strictly necessary. An hydrostatic adjustement can be implemented once the new freesurface elevation is computed simply by imposing the pressure at the nodes on the reference surface as

(54) 
where is the gravity constant.
Eq.(54) takes into account the changes in the freesurface without the need of updating the reference surface nodes. A higher accuracy in the flow solution can be obtained by updating these nodes after a number of time steps.
The algorithms of previous section can be extended to account for the ship motion due to the hydrodynamic forces. Here the ship hull structure can be modelled as a rigid solid defined by the three translations and the three rotations of its center of gravity. Alternatively, the full deformation of the ship structure can be computed by modelling the actual stiffness of the hull, the decks and the different structural members. Indeed, the first option usually suffices for the hull shape optimization.
In both cases the computation of the ship motion involves solving the dynamic equations of the ship structure written as

(55) 
where and are the displacement and acceleration vectors of the nodes discretizing the ship structure, respectively, and are the mass and stiffness matrices of the structure and is the vector of external nodal forces accounting for the fluid flow loads induced by the pressure and the viscous stresses. Clearly the main driving forces for the motion of the ship is the fluid pressure which acts in the form of a surface traction. Indeed Eq.(55) can be augmented with an appropriate damping term. The form of all the relevant matrices and vectors can be found in standard books on FEM for structural analysis (Zienkiewicz and Taylor, Vol 2 2000 [107]).
Solution of Eq.(55) in time can be performed using implicit or fully explicit time integration algorithms. In both cases the values of the nodal displacements, the velocities and the accelerations at time are found.
A simple coupled fluidshipstructure solution in time using, for instance, the fractional step method of Section 6.2 (for ) involves the following steps.
Step 1 Solve for the fractional velocities using Eq.(49). Here use of is recommended.
Step 2 Compute from Eq.(51a) solving a simultaneous system of equations.
Step 3 Compute explicitely the nodal velocities from Eq.(50) with a diagonal mass matrix.
Step 4 Compute explicitely the projected convective variables from Eq.(46) using .
Step 5 Compute explicitely the projected pressure gradients from Eq.(47) using .
Step 6 Compute explicitely the new position of the freesurface elevation from Eq.(53).
Step 7 Compute the movement of the ship by solving the dynamic equations of motion for the ship structure under the hydrodynamic forces induced by the pressures and the viscous stresses .
Step 8 Update the position of the mesh nodes in the fluid domain at by using the mesh update algorithm described in the next section. The updating process can also include the freesurface nodes although this is not strictly necessary to be done at every time step and the hydrostatic adjustment of the pressure acting on the freesurface (Section 6.3) can be used as an alternative.
Obviously, the use of a value different from zero for the or parameters will require an iterative process between steps 1 and 8 until a converged solution for the variables describing the fluid and ship motions is found.
A common option is to update the position of the mesh nodes only when the iterative process for computing the fluid and ship variables has converged. Clearly the regeneration of the mesh is unavoidable when the distorsion of the elements exceed a certain limit.
Different techniques have been proposed for dealing with mesh updating in fluidstructure interaction problems. The general aim of all methods is to prevent element distortion during mesh deformation (Tezduyar, 2001a [95]; Tezduyar et al., 1992bc [100101]).
Tezduyar et al. (1992c) [101] and Chiandussi et al. (2000) [13] have proposed simple method for moving the mesh nodes based on the iterative solution of a fictious linear elastic problem on the mesh domain. In the method introduced in Tezduyar et al. (1992c) [101], the mesh deformation is handled selectively based on the element sizes and deformation modes, with the objective to increase stiffening of the smaller elements, which are typically located near solid surfaces. In Chiandusi et al. (2000) [13] in order to minimize the mesh deformation the “elastic” properties of each mesh element are appropiately selected so that elements suffering greater movements are stiffer. A simple and effective procedure is to select the Poisson's ratio and compute the “equivalent” Young modulus for each element by

(56) 
where are the principal strains, is an arbitrary value of the Young modulus and is a prescribed uniform strain field. and are constant for all the elements in the mesh.
In summary, the solution scheme proposed by Chiandusi et al. (2000) [13] includes the following two steps.
Step 1. Consider the FE mesh as a linear elastic solid with homogeneous material properties characterized by a prescribed uniform strain field , an arbitrary Young modulus and . Solve a linear elastic problem with imposed displacements at the mesh boundary defined by the actual movement of the boundary nodes. An approximate solution to this linear elastic problem, such as that given by the first iterations of a conjugate gradient solution scheme, suffices for practical purposes.
Step 2. Compute the principal strains in each element. Repeat the (approximate) FE solution of the linear elastic problem with prescribed boundary displacements using the values of of Eq.(56).
The previous algorithm is able to treat the movement of the mesh due to changes in position of fully submerged and semisubmerged bodies such as ships. However if the floating body intersects the freesurface, the changes in the analysis domain can be very important as emersion or inmersion of significant parts of the body can occur within a time step. A possible solution to this problem is to remesh the analysis domain. However, for most problems a mapping of the moving surfaces linked to mesh updating algorithm described above can avoid remeshing. The surface mapping technique used by the authors is based on transforming the 3D curved surfaces into reference planes (Figure 4). This makes it possible to compute within each plane the local (inplane) coordinates of the nodes for the final surface mesh accordingly to the changes in the floating line. The final step is to transform back the local coordinates of the surface mesh in the reference plane to the final curved configuration which incorporates the new floating line (García, 1999 [37]; Oñate and García, 2001 [78]).
The transom stern causes a discontinuity in the domain and the solution of the freesurface equation close to this region is inconsistent with the convective nature of the equation. This leads to instability of the wave height close to the transom region. This instability is found experimentally for low speeds. The flow at a sufficient high speed is physically more stable although it still can not be reproduced by standard numerical techniques (Reed et al., 1990) [85].
A solution to this problem is to apply adequate freesurface boundary conditions at the transom boundary. The obvious condition is to fix both the free surface elevation and its derivative along the corresponding streamline to values given by the transom position and the surface gradient. However, prescribing those values can influence the transition between the transom flux and the lateral flux, resulting in unaccurate wave maps.
The method proposed in García and Oñate (2003) [38] is to extend the freesurface below the ship. In this way the necessary Dirichlet boundary conditions imposed at the inflow domain are enough to define a well posed problem. This method is valid both for the wetted and dry transom cases and it is also applicable to ships with regular stern.
This scheme does not work for partially wetted transoms. This situation can occur for highly unsteady flows where wake vortex induces the freesurface deformation and the flow remains adhered to the transom. To favour the computation of the freesurface, an artificial viscosity term is added to the freesurface equation in the vicinity of the transom in these cases.
The Lagrangian formulation is an effective (and relatively simple) procedure for modelling the flow of fluid particles undergoing severe distorsions such as water jets, high amplitude waves, braking waves, water splashing, filling of cavities, etc. Indeed the Lagrangian formulation is very suitable for treating ship hydrodynamic problems where the ship undergoes large motions. An obvious “a priori” advantage of the Lagrangian formulation is that both the ship and the fluid motion are defined in the same frame of reference.
The Lagrangian fluid flow equations are obtained by noting that the velocity of the mesh nodes and that of the fluid particles are the same. Hence the relative velocities are zero in Eq.(21) and the convective terms vanish in the momentum equations, while the rest of the fluid flow equations remain unchanged. Also, the solution of the free surface equation is not needed in the Lagrangian formulation, as the free surface motion is modelled naturally by computing the displacement of (all) the fluid particles at each time step. In any case, the new position of the free surface nodes must be properly identified as described below.
The FEM algorithms for solving the Lagrangian flow equations are very similar to those for the ALE description presented earlier. We will focus in the second order fractional step algorithm of Section 6.2 (for and ) accounting also for fluidship interaction effects.
Step 1 Compute explicitely a predicted value of the velocities as

(57) 
Note that matrices and of Eq.(48) emanating from the convective terms have been elliminated.
Step 2 Compute from Eq.(52).
Step 3 Compute explicitely from Eq.(50) with .
Step 4 Compute explicitely from Eq.(46).
Step 5 Solve for the motion of the ship structure by integrating Eq.(55).
Step 6 Update the mesh nodes in a Lagrangian manner as

(58) 
Step 7 Generate a new mesh and identify the new free surface nodes.
Further details of above algorithm can be found in Idelsohn et al. (2002 [53], 2003a [54]) and Oñate et al. (2003 [80], 2004 [81]).
The mesh regeneration process can be effectively performed using the extended Delaunay Tesselation described in Idelsohn et al. (2003b,c [5556])). This method allows the fast generation of good quality meshes combining four node tetrahedra (or three node triangles in 2D) with hexahedra and non standard polyhedra such as pentahedra (or quadrilaterals and pentagons in 2D) where linear shape functions are derived using nonSibsonian interpolation rules. The mesh regeneration can take place at each time step (this has been the case for the examples presented in the paper), after a prescribed number of time steps, or when the nodal displacements induce significant element distorsions.
The identification of the freesurface nodes in the Lagrangian analysis can be made using the Alpha Shape method. This is based on the search of all nodes which are on an empty Voronoi sphere with a radius greater than a specified distance defined in terms of the mesh size. For details see Edelsbrunner and Mucke (1994) [34] and Idelsohn et al. (2002 [53], 2003a,b,c [5456]). A known value of the pressure (typically zero pressure or the atmospheric pressure) is prescribed at the free surface nodes.
The conditions of prescribed velocities or pressures at the solid boundaries in the Lagrangian formulation are usually applied on a layer of nodes adjacent to the boundary. These nodes typically remain fixed during the solution process. Contact between water particles and the solid boundaries is accounted for by the incompressibility condition which naturally prevents the water nodes to penetrate into the solid boundaries. This simple way to treat the waterwall contact is another attractive feature of the Lagrangian flow formulation. For details see Idelsohn et al. (2002 [53], 2003a [54]) and Oñate et al. (2003 [80], 2004 [81]).
A simple and yet effective way to analyze the rigid motion of solid bodies in fluids with the Lagrangian flow description is to model the solid as a fluid with a viscosity much higher than that of the surrounding fluid. The fractional step scheme of Section 10 can be readily applied skipping now step 5 to solve for the motion of both fluid domains (the actual fluid and the fictitious fluid modelling the rigid body displacements of the solid). An example of this type is presented in Sections 14.6.4 and 14.6.7.
Indeed this approach can be further extended to account for the elastic deformation of the solid treated now as a viscoelastic fluid. This will however introduce some complexity in the formulation and the full coupled fluidstructure interaction scheme previously described is preferable.
The evaluation of the stabilization parameters is a crucial issue in stabilized methods. Most existing methods use expressions which are direct extensions of the values obtained for the simplest 1D case. It is also usual to accept the so called SUPG assumption, i.e. to admit that vector has the direction of the velocity field. This restriction leads to instabilities when sharp layers transversal to the velocity direction are present. This deficiency is usually corrected by adding a shock capturing or crosswind stabilization term (Hughes and Mallet, 1986 [48]; Tezduyar and Partk, 1986 [98]; Codina, 1993 [15]). Indeed, in the FIC formulation the components of introduce the necessary stabilization along both the streamline and transversal directions to the flow.
Excellent results have been obtained in all problems solved with the ALE formulation using linear tetrahedra with the value of the characteristic length vector defined by

(59) 
where and and are the “streamline” and “cross wind” contributions given by

(60) 

(61) 
where are the vectors defining the element sides ( for tetrahedra).
The form chosen for the characteristic length parameters is similar to that proposed by other authors (see references of previous paragraph and also Donea and Huerta, 2003 [31] and Tezduyar, 2001b [96]).
As for the freesurface equation the following value of the characteristic length vector has been taken

(62) 
The streamline parameter has been obtained by Eq.(60) using the value of the velocity vector over the 3 node triangles discretizing the freesurface and .
The cross wind parameter has been computed by

(63) 
The crosswind terms in eqs.(59) and (62) take into account the gradient of the solution in the stabilization parameters. This is a standard assumption in shockcapturing stabilization procedures.
A more consistent evaluation of based on a diminishing residual technique can be found in Oñate and García (2001) [78].
The detailed discussion on the treatment of turbulent effects in the flow equations falls outside the scope of this chapter. Indeed any of the existing turbulence models is applicable.
In the examples presented next we have chosen a turbulence model based on the Reynolds averaged NavierStokes equations where the deviatoric stresses are computed as sum of the standard viscous contributions and the so called Reynold stresses. Here we have chosen the Boussinesq assumption leading to a modification of the viscosity in the standard NavierStokes equations as sum of the “physical” viscosity and a turbulent viscosity .
One of the simplest and more effective choices for is the Smagorinski LES model giving

(64) 
where is the element size and is a constant (). In the examples analyzed has been taken equal to and for 2D and 3D situations, respectively.
Many other options are possible such as the one and two equations turbulence models (i.e. the model and the and models) and the algebraic stress models. For further details the reader is refered to specialized publications (Celik el al., 1982 [12]; Wilcox, 1994 [105]).
The examples chosen show the applicability of the ALE and Lagrangian formulations presented to solve ship hydrodynamics problems. The fractional step algorithm of Section 6.2 with and has been used. The first ALE example is the flow past a submerged NACA 0012 profile. The next ALE examples include the analysis of a Wigley hull, a scale model of a commercial ship and two American Cup racing sail boats. Numerical results obtained with linear tetrahedral meshes are compared with experimental data in all cases.
The last series of examples show applications of the Lagrangian formulation to simple schematic 2D ship hydrodynamics situations.
A submerged NACA0012 profile at angle of attack is studied using a 3D finite element model. This configuration was tested experimentally by Duncan (1983) [33] and modelled numerically using the Euler equations by several authors (Hino et al., 1993 [41]; Löhner et al., 1998 [64], 1999 [65]; Idelsohn et al., 1999 [52]). The submerged depth of the airfoil is equal to the chord length L. The Froude number for all the cases tested was set to where is the incoming flow velocity at infinity.
The stationary free surface and the pressure distribution in the domain are shown in Figure 5. The nondimensional wave heights compare well with the experimental results.
We study here the well known Wigley Hull, given by the analytical formula where and are the beam and the draft of the ship hull at still water.
The same configuration was tested experimentally (Noblesse and McCarthy, 1983 [69]; Wigley, 1983 [104] ) and modelled numerically by several authors (Farmer et al., 1993 [35]; Storti et al., 1998 [91]; Idelsohn et al., 1999 [52]; Löhner et al., 1999 [65]). We use an unstructured 3D finite element mesh of linear tetrahedra, with a reference surface of triangles, partially represented in Figure 6. A Smagorinsky turbulence model was chosen.
Three different viscous cases were studied for . In the first case the volume mesh was considered fixed, not allowing freesurface nor ship movements. Next, the volume mesh was updated due to freesurface movement, while keeping the model fixed. The third case corresponds to the analysis of a real free model including the mesh updating due to freesurface displacement and ship movement (sinkage and trim).
Table 1 shows the obtained total resistance coefficient in the three cases studied compared with experimental data.
Experimental  Numerical  
Test 1  
Test 2  
Test 3 
Numerical values obtained for sinkage and trim were and , respectively, while experiments gave and .
Figure 6a shows the pressure distribution obtained near the Wigley hull for the free model. Some streamlines have also been plotted. The mesh deformation in this case is shown in Figure 6a
Comparison of the obtained body wave profile with experimental data for the free and fixed models is shown in Figure 6b.Significant differences are found close to stern for the fixed model. The freesurface contours for the truly free ship motion are shown in Figure 6c.
The example is the analysis of the KVLCC2 benchmark model. Here a partially wetted tramsom stern is expected due to the low Froude number of the test. Figure 7 shows the NURBS geometry obtained from the Hydrodynamic Performance Research team of Korea (KRISO). Numerical results are compared with the experimental data available from KRISO (2000).
The smallest element size used was m and the largest m. The surface mesh chosen is also shown in Figure 7. The 3D mesh included tetrahedra. The tramsom stern flow model described was used.
Test 1. Wave pattern calculation. The main characteristics of the analysis are listed below:
The turbulence model chosen was the model. Figures 8 and 9 show the wave profiles on the hull and in a cut at obtained for Test 1, compared to experimental data. The results are quantitatively good.
Test 2. Wake analysis at different planes. Several turbulence models were used (Smagorinsky, and model) in order to verify the quality of the results. Here, only the results from the model are shown. The velocity maps obtained even for the simplest model were qualitatively good. The main characteristics of this analysis are:
Figures 10 and 11 present some results for Test 2. Figure 10 shows the contours of the axial (X) component of the velocity on a plane at from the orthogonal aft. Figure 11 shows the map of the kinetic energy at this plane. Experimental results are also plotted for comparison. Further results are reported in García and Oñate (2003)[38].
The Rioja de España boat representing Spain in the American Cup's edition of 1995 was analyzed. Figure 12 shows the geometry of the boat based on standard NURBS surfaces. Numerical results are compared with an extrapolation of experimental data obtained in CEHIPAR basin (Spain) using a 1/3.5 scale model. Resistance tests were performed with the model free to sink and trim. Experimental data include drag, lift, sinkage, trim angles and wave profiles at (real scale) from symmetry plane. Every model was towed at equivalent velocities of , , , , and knots.
Numerical analysis were carried out at real scale. Characteristics of unstructured grids of four node linear tetrahedra used are shown in Table 2 together with the parameters of some of the model studied.
Test  Geometry  Heel  Drift  Symmetry  Elements  nodes 
E0D0  Hull, bulb and keel  Yes  700 000  175 000  
E15D2  Hull, bulb, keel and rudder  No  1 500 000  380 000  
E15D4  Hull, bulb, keel and rudder  No  1 500 000  380 000  
E25D2  Hull, bulb, keel and rudder  No  1 500 000  380 000 
All grids used were generated with the same quality criteria and using element sizes from to . Some details of the grid used in the E0D0 case are shown in Figures 13 and 14. A turbulence model in combination with an extended law of the wall was chosen for the analysis.
A time step of was used and this sufficed to achieve a steady state solution in all cases. Additional calculations were carried out with and , in order to verify the influence of the time increment in the accuracy of the results. No significant influence was detected for the selected results.
Figure 16 shows a comparison of simulated and experimental wave profiles. The origin of the x axis has been taken at the fore perpendicular and the x+ sense is afore. In the non symmetric case, the measurements were performed at the opposite side of the heeled board. The ratio of maximum amplitudes of fluctuations (noise) and waves in the experimental measurements varies from to
Figures 16 and 17 show pressure contours on the bulb and keel for different cases, corresponding to a velocity of . Figures 18 to 19 show some of the wave patterns obtained for a velocity of . Figures 21 and 22 show some perspective views, including pressure and velocity contours, streamlines and cuts for some cases analyzed. Finally Figures 23 and 24 show resistance graphs where numerical results are compared with values extrapolated from experimental data. For further details see García et al. (2002) [39].
The finite element formulation presented was also applied to study the racing sail boat Bravo España participating in the 1999 edition of the American Cup. The mesh of linear tetrahedra used is shown in Figure 25. Results presented in Figures 25–27 correspond to a non symmetrical case including appendages. Good comparison between experimental data and numerical results was again obtained Further details can be found in García and Oñate (2002) [38].
A number of simple problems have been solved in order to test the capabilities of the Lagrangian flow formulation to solve ship hydrodynamics problem. All examples have been analyzed with the algorithm described in Section 10. We note again that the regeneration of the mesh has been performed at each time step.
The first example is a very schematic representation of a ship when is hit by a big wave produced by the collapse of a water recipient (Fig. 28). The ship can not move and initially the freesurface near the ship is horizontal. Fixed nodes represent the ship as well as the domain walls. The example tests the suitability of the Lagrangian flow formulation to solve waterwall contact situations even in the presence of curved walls. Note the breaking and splashing of the waves under the ship prow and the rebound of the incoming wave. It is also interesting to see the different waterwall contact situations at the internal and external ship surfaces and the moving freesurface at the back of the ship.
In the next example (Fig. 29) the same ship of the previous example moves now horizontally at a fixed velocity in a water reservoir. The freesurface, which was initially horizontal, takes the correct position at the ship prow and stern. Again, the complex waterwall contact problem is naturally solved in the curved prow region.
The example shown in Figure 30 is the analysis of a rotating water mill semi submerged in water. The blades of the mill are treated as a rigid body with an imposed rotating velocity, while the water is initially in a stationary flat position. Fluid structure interactions with freesurfaces and water fragmentation are well reproduced in this example.
The next example shows an initially stationary recipient with a floating piece of wood where a wave is produced on the left side. The wood has been simulated by a liquid of higher viscosity as described in Section 11. The wave intercepts the wood piece producing a breaking wave and displacing the floating wood as shown in Figure 31.
In the example of Figure 32 the motion of a fictitious rigid ship hit by an incoming wave is analyzed. The dynamic motion of the ship is induced by the resultant of the pressure and the viscous forces acting on the ship boundaries. We note that the horizontal displacement of the gravity center of the ship was fixed to zero. In this way, the ship moves only vertically although it can freely rotate. The position of the ship boundary at each time step defines a moving boundary condition for the free surface particles in the fluid.
Figure 32 shows instants of the motion of the ship and the surrounding fluid. It is interesting to see the breaking of a wave on the ship prow as well as on the stern when the wave goes back. Note that some water drops slip over the ship deck.
Figure 33 shows the analysis of a similar problem using precisely the same approach. The section of the ship analyzed corresponds now to that of a real container ship. Differently to the previous case the rigid ship is free to move laterally due to the sea wave forces. The objective of the study was to asses the influence of the stabilizers in the ship roll. The figures show clearly how the lagrangian method predicts the ship and wave motions in a realistic manner.
Figure 34 represents a transversal section of a tankship and a wave approaching to it. The tankship is modelled as a rigid body and it carries a liquid inside which moves freely within the ship domain.
Initially () the ship is released from a fixed position and the equilibrium configuration found is consistent with Arquimides principle. During the following time steps the external wave hits the ship and both the internal and the external fluids interact with the ship boundaries. At times and breaking waves and some water drops slipping along the ship deck can be observed. Figure 35 shows the pressure pattern at two time steps.
The last examples show the capacity of the Lagrangian formulation to analyze the motion of solids within water. In the first example a solid cube is initially free and falls down within a water recipient. In this example, the rigid solid is modeled first as a fictitious fluid with a higher viscosity, similarly as for the floating solid of Section 14.6.4. The results of this analysis are shown in Figure 36. Note that the method reproduces very well the interaction of the cube with the free surface as well as the overall sinking process. A small deformation of the cube is produced. This can be reduced by increasing further the fictitious viscosity of the cube particles.
The same problem is analyzed again considering now the cube as a rigid solid subjected to pressure and viscous forces acting in its boundaries. The resultant of the fluid forces and the weight of the cube are applied to the center of the cube. These forces govern the displacement of the cube which is computed by solving the dynamic equations of motion as described in the fractional step algorithm of Section 10, similarly as for the rigid ships of the three previous examples. Similarly as in those cases the moving cube contours define a boundary condition for the fluid particles at each time step.
Initially the solid falls down freely due to the gravity forces (Figure 37). Once in contact with the water surface, the alphashape method recognizes the different boundary contours which are shown with a thick line in the figure. The pressure and viscous forces are evaluated in all the domain and in particular on the cube contours. The fluid forces introduce a negative acceleration in the vertical motion until, once the cube is completely inside the water, the vertical velocity becomes zero. Then, Arquimides forces bring the cube up to the freesurface. It is interesting to observe that there is a rotation of the cube. The reason is that the center of the floating forces is higher in the rotated position than in the initial ones.
Figure 38 shows a repetition of the same problem showing now all the finite elements in the mesh discretizing the fluid. We recall that in all the Lagrangian problems here described the mesh in the fluid domain is regenerated at each time step combining linear triangles and quadrilateral elements as described in Section 10. Note that some fluid particles separate from the fluid domain. These particles are treated as free boundary points with zero pressure and hence fall down due to gravity.
It is interesting to see that the final position of the cube is different from that of Figure 37. This is due to the unstable character of the cube motion. A small difference in the numerical computations (for instance in the mesh generation process) shifts the movement of the cube towards the right or the left. Note however that a final rotated equilibrium position is found in both cases.
Further examples of the Lagrangian flow formulation can be found in Oñate et al. (2003 [80], 2004 [81]) and Idelsohn et al. (2002 [53], 2003a [54]).
An overview of different finite element schemes for solving ship hydrodynamics problems has been presented. The necessary stabilization in the numerical computation has been provided by the finite calculus formulation of the fluid flow equations accounting for surface wave effects.
Stabilized finite algorithms for solving a variety of ship hydrodynamics situations using ALE and fully Lagrangian descriptions have been presented. Both monolithic and fractional step algorithms have been derived. The interaction of the motion of the ship structure with the fluid equations can be accounted for in a straight forward manner within the general flow solution schemes. The ALE formulation is particularly adequate for ship hydrodynamics problems involving freesurface waves of moderate amplitude. The Lagrangian flow formulation allows to solve in a simple and effective manner ship hydrodynamics problems involving large motions of the freesurface and complex fluidstructure interaction situations.
The authors thank Prof. R. Löhner, Prof. P. Bergan, Dr. C. Sacco and Dr. H. Sierra for many useful discussions.
Thanks are also given to Dr. F. del Pin and Mr. N. Calvo for their help and involvement in solving the problems with the Lagrangian formulation.
The authors are grateful to Copa America Desafio Español SA and CEHIPAR for providing the geometry and experimental data for the racing boats analyzed and to Det Norske Veritas for providing the geometry of the ship with stabilizers for the analysis shown in Figure 33.
The Spanish company COMPASS Ingeniería y Sistemas SA provided the computer code Tdyn for the numerical simulation of the examples solved with the ALE formulation (Tdyn, 2004 [92]). This support is gratefully acknowledged.
[1] Alessandrini B and Delhommeau G. A multigrid velocitypressurefreesurface elevation fully coupled solver for calculation of turbulent incompressible flow around a hull. In Proc. of the 21st Symp. on Naval Hydrodynamics, Trondheim, Norway, 1996, 40–55.
[2] Alessandrini B and Delhommeau G. A Fully Coupled NavierStokes Solver for Calculation of Turbulent Incompressible FreeSurface Flow Past a Ship Hull. J. Numer. Meth. Fluids 1999; 29:125–142.
[3] Aliabadi S and Shujaee S. freesurface flow simulations using parallel finite element method. Simulation 2001; 76(5):257–262.
[4] Aliabadi S, Abedi J and Zellars B. Parallel finite element simulation of mooring forces on floating objects. Submitted to Int. J. Num. Meth. Fluids 2002.
[5] Azcueta R, Muzaferija S and Peric M. Computation of breaking bow waves for a very full hull ship. 7th Int. Conf. on Num. Ship Hydro. 1999, Nantes, France, 6.21, 11.
[6] Baba E and Takekuma K. A study on freesurface flow around bow of slowly moving full forms. J. Soc. Naval Architects Japan 1975; 137.
[7] Bai KJ and McCarthy JH (eds). Proceedings of the Workshop on Ship WaveResistance Computations, Bethesda, MD, USA, 1979.
[8] Beck, RF, Cao, Y and Lee TH. Fully nonlinear water wave computations using the desingularized method. Proceedings Sixth International Conference on Numerical Ship Hydrodynamics, The University of Iowa, Iowa City, Aug. 1993.
[9] Brezzi F, Franca LP, Hughes TJR and Russo A. . Comput. Methods Appl. Mech. Engrg. 1997; 145:329–339.
[10] Briley WR, Neerarambam SS and Whitfield DL. Multigrid algorithm for threedimensional incompressible highReynolds number turbulent flows. AIAA Journal 1995; 33 (1):2073–2079.
[11] Brooks A and Hughes TJR. Streamline upwind/PetrovGalerkin formulation for convection dominated flows with particular emphasis on the incompressible NavierStokes equations. Comput. Methods Appl. Mech. Engrg 1982; 32:199–259.
[12] Celik I; Rodi W and Hossain MS. Modelling of freesurface proximity effects on turbulence. Proc. Refined Modelling of Flows, Paris, 1982.
[13] Chiandusi G, Bugeda G and Oñate E. A simple method for update of finite element meshes. Commun, Numer. Meth. Engng. 2000; 16:1–9.
[14] Chorin AJ. A numerical solution for solving incompressible viscous flow problems. J. Comp. Phys. 1967; 2:12–26.
[15] Codina R. A discontinuitycapturing crosswind dissipation for the finite element solution of the convectiondiffusion equation. Comput. Methods Appl. Mech. Engrg. 1993; 110:325–342.
[16] Codina R. Comparison of some finite element methods for solving the diffusionconvectionreaction equation. Comput. Methods Appl. Mech. Engrg. 1998; 156:185–210.
[17] Codina R. On stabilized finite element methods for linear systems of convectiondiffusionreaction equation. Comput. Methods Appl. Mech. Engrg. 2000; 188:61–83.
[18] Codina R. Stabilization of incompressibility and convection through orthogonal subscales in finite element methods. Comput. Methods Appl. Mech. Engrg. 2000; 191:4295–4321.
[19] Codina R. A stabilized finite element method for generalized stationary incompressible flows. Computer Methods in Applied Mechanics and Engineering 2001; 190:26812706.
[20] Codina R. Pressure stability in fractional step finite element methods for incompressible flows. Journal of Computational Physics 2001; 170:112–140.
[21] Codina R. Stabilized finite element approximation of transient incompressible flows using orthogonal subscales. Comput. Methods Appl. Mech. Engrg. 2002; 191:4295–4321.
[22] Codina R and Blasco J. A finite element formulation for the Stokes problem allowing equal velocity  pressure interpolation. Comput. Methods Appl. Mech. Engrg. 1997; 143:373–391.
[23] Codina R, Vazquez M, Zienkiewicz OC. A general algorithm for compressible and incompressible flow  Part III. The semiimplicit form. Int. J. Num. Meth. in Fluids 1998; 27:13–32.
[24] Codina R and Blasco J. Stabilized finite element method for the transient NavierStokes equations based on a pressure gradient operator. Comput. Methods in Appl. Mech. Engrg. 2000; 182:277–301.
[25] Codina R and Zienkiewicz OC. CBS versus GLS stabilization of the incompressible NavierStokes equations and the role of the time step as stabilization parameter. Communications in Numerical Methods in Engineering 2002; 18 (2):99–112.
[26] Codina R, Schafer U and Oñate E. Mould filling simulation using finite elements. Int. J. Num. Meth. Heat Fluid Flows 1994; 4:291–310.
[27] Cruchaga MA and Oñate E. A finite element formulation for incompressible flow problems using a generalized streamline operator. Comput. Methods in Appl. Mech. Engrg. 1997; 143:49–67.
[28] Cruchaga MA and Oñate E. A generalized streamline finite element approach for the analysis of incompressible flow problems including moving surfaces. Comput. Methods in Appl. Mech. Engrg. 1999; 173:241–255.
[29] Dawson CW. A practical computer method for solving shipwave problems. Proc. 2nd Int. Conf. Num. Ship Hydrodynamics, USA, 1977.
[30] Donea J. A TaylorGalerkin method for convective transport problems. Int. J. Num. Meth. Engng. 1984; 20:101–119.
[31] Donea J and Huerta A. Finite element method for flow problems, J. Wiley, 2003.
[32] Douglas J, Russell TF. Numerical methods for convection dominated diffusion problems based on combining the method of characteristics with finite element or finite difference procedures. SIAM J. Numer. Anal. 1982; 19:871.
[33] Duncan JH. The breaking and nonbreaking wave resistance of a twodimensional hydrofoil. J. Fluid Mech. 1983; 126.
[34] Edelsbrunner H and Mucke EP. Three dimensional alphashape. ACM Transaction on Graphics 1994; 3:43–72.
[35] Farmer JR, Martinelli L and Jameson A. A fast multigrid method for solving incompressible hydrodynamic problems with freesurfaces. AIAA J. 1993; 32 (6):1175–1182.
[36] García J, Oñate E, Sierra H, Sacco C and Idelsohn SR. A stabilized numerical method for analysis of ship hydrodynamics. ECCOMAS CFD98, Papaliou K et al. (Eds), J. Wiley, 1998.
[37] García J. A finite element method for analysis of naval structures (in Spanish). Ph.D. Thesis, Univ. Politecnica de Catalunya, Barcelona, Spain, December, 1999.
[38] García J and Oñate E. An unstructured finite element solver for ship hydrodynamic problems. J. Appl. Mech. 2003; 70:18–26.
[39] García J, LucoSalman R, Salas M, LópezRodríguez M and Oñate E. An advanced finite element method for fluiddynamic analysis of America's Cup boat. High Performance Yatch Design Conference, Auckland, 4–6, December, 2002.
[40] Hansbo P and Szepessy A. A velocitypressure streamline diffusion finite element method for the incompressible NavierStokes equations. Comput. Methods Appl. Mech. Engrg. 1990; 84:175–192.
[41] Hino T, Martinelli L and Jameson A. A finite volumen method with unstructured grid for freesurface flow. In Proc. of the 6th Int. Conf.Num. Ship Hydrodynamics, Iowa City, Iowa, 1993, 173–194.
[42] Hirsch C. Numerical computation of internal and external flow, J. Wiley, Vol. 1 1988, Vol. 2, 1990.
[43] Hirt CW and Nichols BD. Volume of fluid (VOF) method for the dynamics of free boundaries. J. Comput. Physics 1981; 39:201–225.
[44] Hughes TJR. Multiscale phenomena: Green functions, subgrid scale models, bubbles and the origins of stabilized methods. Comput. Methods Appl. Mech. Engrg 1995; 127:387–401.
[45] Hughes TJR and Brooks, AN. A multidimensional upwind scheme with no crosswind diffusion. Finite Element Methods for Convection Dominated Flows 1979; AMDVol. 34, ASME, New York.
[46] Hughes TJR and Tezduyar, TE. Finite element methods for firstorder hyperbolic systems with particular emphasis on the compressible Euler equations. Computer Methods in Applied Mechanics and Enginerring 1984; 45:217–284.
[47] Hughes TJR, Franca LP and Balestra M. A new finite element formulation for computational fluid dynamics. V Circumventing the BabuskaBrezzi condition: A stable PetrovGalerkin formulation of the Stokes problem accomodating equalorder interpolations. Comput. Methods Appl. Mech. Engrg. 1986; 59:85–89.
[48] Hughes TJR and Mallet M. A new finite element formulations for computational fluid dynamics: III. The generalized streamline operator for multidimensional advectivediffusive systems. Comput Methods Appl. Mech. Engrg. 1986; 58:305–328.
[49] Hughes TJR and Mallet M. A new finite element formulations for computational fluid dynamics: IV. A discontinuity capturing operator for multidimensional advectivediffusive system. Comput. Methods Appl. Mech. Engrg. 1986; 58:329–336.
[50] Hughes TJR, Franca LP and Hulbert GM. A new finite element formulation for computational fluid dynamics: VIII. The Galerkin/leastsquares method for advectivediffusive equations. Comput. Methods Appl. Mech. Engrg. 1989; 73:173–189.
[51] Idelsohn SR, Storti M and Nigro N. Stability analysis of mixed finite element formulation with special mention to stabilized equalorder interpolations. Int. J. for Num. Meth. in Fluids 1995; 20:10031022.
[52] Idelsohn SR, Oñate E and Sacco C. Finite element solution of freesurface shipwave problem. Int. J. Num. Meth. Engng. 1999; 45:503–508.
[53] Idelsohn SR, Oñate E, Del Pin F and Calvo N. Lagrangian formulation: the only way to solve some freesurface fluid mechanics problems. Fith World Congress on Computational Mechanics, Mang HA, Rammerstorfer FG and Eberhardsteiner J. (eds), July 7–12, Viena, Austria, 2002, Web: wccm.tuwien.ac.at.
[54] Idelsohn SR, Oñate E and Del Pin F. A Lagrangian meshless finite element method applied to fluidstructure interaction problems. Computer and Structures 2003; 81:655–671.
[55] Idelsohn SR, Oñate E, Calvo N and Del Pin F. The meshless finite element method. Int. J. Num. Meth. Engng 2003; 58(6):893–912.
[56] Idelsohn SR, Calvo N and Oñate E. Polyhedrization of an arbitrary point set. Comput. Method Appl. Mech. Engng. 2003; 192(22–24):2649–2668.
[57] Janson C and Larsson L. A method for the optimization of ship hulls from a resistance point of view. Proc. of the 21th Symp. on Naval Hidrodynamics, Trondheim, Norway, 1996.
[58] Jenson G and Soding H. Ship wave resistance computation. Finite Approximations in Fluid Mechanics 1989; II, Vol.25.
[59] Kim YH and Lucas T. Nonlinear ship waves. Proc. 18th Symp. on Naval Hidrodynamics, MI, USA, 1990.
[60] KRISO (Korea Research Institute of Ships and Ocean Engineering).http://www.iihr.uiowa.edu/gothenburg2000/KVLCC/tanker.html.
[61] Larsson L. Procedings of the 1980 SSPAITTC Workshop on Ship Boundary Layers. SSPA Publication No. 90, 1981.
[62] Larsson L, Regnströn B, Broberg L, Li DQ and Janson CE. Failures, fantasies and feats in the theoretical/numerical prediction of ship performance. 22d Symposium on Naval Hydrodynamics, Washington DC, 10–14, August, 1998.
[63] Löhner R, Morgan K and Zienkiewicz OC. The solution of nonlinear hyperbolic equation systems by the finite element method. Int. J. Num. Meth. in Fluids 1984; 4:1043.
[64] Löhner R, Yang C and Oñate E. Viscous free surface hydrodynamics using unstructured grids. 22nd Simposium on Naval Hydrodynamics, Washington CD, September 1998.
[65] Löhner R, Yang C, Oñate E and Idelsohn SR. An unstructured gridbased parallel freesurface solver. Appl. Num. Math. 1999; 31:271–293.
[66] Luo H, Baum JD and Löhner R. A finite volume scheme for hydrodynamic free boundary problems on unstructured grids. AIAA950668, 1995.
[67] Miyata H. Timemarching CFD simulation for moving boundary problems. In Proc. of the 21st Symp. on Naval Hydrodynamics, Trondheim, Norway, 1996, 1–21.
[68] Newman JN. Linearized wave resistance theory. Int. Seminar on Wave Resistance, Tokyo/Osaka, J. Soc. Naval Architects Japan, 1976.
[69] Noblesse F and McCarthy JH (eds). Proceedings of the Second DTNSRDC Workshop on Ship WaveResistance Computations, Bethesda, MD, USA, 1983.
[70] Oñate E. Derivation of stabilized equations for advectivediffusive transport and fluid flow problems. Comput. Methods Appl. Mech. Engrg. 1998; 151 (12):233–267.
[71] Oñate E. A stabilized finite element method for incompressible viscous flows using a finite increment calculus formulation. Comput. Methods Appl. Mech. Engrg. 2000; 182 (1–2):355–370.
[72] Oñate E. Possibilities of finite calculus in computational mechanics. To be published in Int. J. Num. Meth. Engng 2004; 59(54).
[73] Oñate E and Idelsohn SR. A mesh free finite point method for advectivediffusive transport and fluid flow problems. Computational Mechanics 1988; 21:283–292.
[74] Oñate E, García J and Idelsohn SR. Computation of the stabilization parameter for the finite element solution of advectivediffusive problems. Int. J. Num. Meth. Fluids 1997; 25:1385–1407.
[75] Oñate E, García J and Idelsohn SR. An alphaadaptive approach for stabilized finite element solution of advectivediffusive problems with sharp gradients. New Adv. in Adaptive Comput. Methods in Mech., Ladeveze P and Oden JT (eds), Elsevier, 1998.
[76] Oñate E and Manzan M. A general procedure for deriving stabilized spacetime finite element methods for advectivediffusive problems. Int. J. Num. Meth. Fluids 1999; 31:203–221.
[77] Oñate E and Manzan M. Stabilization techniques for finite element analysis of convection diffusion problems. In Computational Analysis of Heat Transfer, Comini G and Sunden B (eds), WIT Press: Southampton, 2000.
[78] Oñate E and García J. A finite element method for fluidstructure interaction with surface waves using a finite calculus formulation. Comput. Methods Appl. Mech. Engrg. 2001; 191:635–660.
[79] Oñate E, García J, Bugeda G and Idelsohn SR. A general stabilized formulation for incompressible fluid flow using finite calculus and the finite element method. In Towards a New Fluid Dynamics with its Challenges in Aeronautics, Periaux J, Chamption D, Pironneau O and Thomas Ph. (eds), CIMNE, Barcelona, Spain 2002.
[80] Oñate E, Idelsohn SR and del Pin, F. Fully Lagrangian formulation for fluidstructure interaction problems using the finite element method and finite calculus. Computational Mechanics  Theory and Practice, K.M. Mathisen, T. Kvamsdal and K.M. Okstad (eds.), CIMNE, Barcelona, Spain 2003.
[81] Oñate E, Idelsohn SR, del Pin F, Calvo N. and Aubry R. The particle finite element method. An overview. To be published in Int. J. Computational Methods, 2004.
[82] Peraire J, Morgan K and Peiro J. The simulation of 3D incompressible flows using unstructured grids. In Frontiers of Computational Fluid Dynamics, Caughey DA and Hafez MM. (eds), Chapter 16, J. Wiley, 1994.
[83] Pironneau O. On the transportdiffusion algorithm and its applications to the NavierStokes equations. Numer. Math. 1982; 38:309.
[84] Raven H. A solution method for the nonlinear ship resistance problem. Doctoral Thesis, Maritime Research Institute, Netherlands, 1996.
[85] Reed AM, Telste JG, Scragg CA and Liepmann D. Analysis of transon stern flows. Proc. 17th Symp. on Naval Hidrodynamics, The Hague, The Netherlands, 1990.
[86] Sheng C, Taylor LK and Whitfield DL. Implicit lowerupper/approximatefactorization schemes for incompressible flows. Journal of Computational Physics 1996; 128 (1):32–42.
[87] Soding H. Advances in panel methods. Proc. of the 21th Symp. on Naval Hidrodynamics, Trondheim, Norway, 1996.
[88] Storti M, Nigro N and Idelsohn SR. Steady state incompressible flows using explicit schemes with an optimal local preconditioning. Comput. Meth. Appl. Mech. and Engrg. 1995; 124:231–252.
[89] Storti M, Nigro N, Idelsohn SR. equalorder interpolations. A unified approach to stabilize the incompressible and convective effects. Comput. Meth. Appl. Mech. and Engrg. 1997; 143:317–331.
[90] Storti M, d'Elia J and Idelsohn SR. Algebraic discrete nonlocal (DNL) absorbing boundary condition for the ship wave resistance problem. J. Computational Physics 1998a; 146(2): 570–602.
[91] Storti M, d'Elia J and Idelsohn SR. Computing ship wave resistance form wave amplitude with the DNL absorbing boundary condition. Comun. in Numer. Meth. in Engng 1998b; 14:997–1012.
[92] Tdyn. A finite element code for fluiddynamic analysis, COMPASS Ingeniería y Sistemas SA, Barcelona, Spain, www.compassis.com, (2004).
[93] Tajima M and Yabe T. Simulation on slamming of a vessel by CIP method. J. Phys. Soc. Japan 1999; 68:2576–2584.
[94] Tezduyar TE. Stabilized finite element formulations for incompressible flow computations. Advances in Applied Mechanics 1991; 28:1–44.
[95] Tezduyar TE. Finite element methods for flow problems with moving boundaries and interfaces. Arch. Comput. Meth. Engng. 2001; 8:83–130.
[96] Tezduyar TE. Adaptive determination of the finite element stabilization parameters. Proceeding of the ECCOMAS Computational Fluid Dynamics Conference 2001, Swansea, Wales, United Kingdom, CDROM (2001).
[97] Tezduyar TE and Hughes, TJR. Finite element formulations for convection dominated flows with particular emphasis on the compressible Euler equations. AIAA Paper 830125, Proceedings of AIAA 21st Aerospace Sciences Metting, Reno, Nevada, 1983.
[98] Tezduyar TE and Park YJ. Discontinuity capturing finite element formulations for nonlinear convectiondiffusionreaction problems. Comput. Meth. Appl. Mech. and Engrg. 1986; 59:307–325.
[99] Tezduyar TE, Mittal S, Ray SE and Shih R. Incompressible flow computations with stabilized bilinear and linear equalorder interpolation velocity–pressure elements. Comput. Methods Appl. Mech. Engrg. 1992; 95:221–242.
[100] Tezduyar TE, Behr M and Liou J. A new strategy for finite element computations involving moving boundaries and interfaces  the deformingspatialdomain/spacetime procedure: I. The concept and the preliminary tests. Comput. Methods in Appl. Mech. and Engrg. 1992; 94:339–351.
[101] Tezduyar TE, Behr M, Mittal S and Johnson AA. Computation of unsteady incompressible flows with the stabilized finite element methods. Spacetime formulations, iterative strategies and massively parallel implementations. New Methods in Transient Analysis (eds. P. Smolinki, W.K. Liu, G. Hulbert and K. Tamma), AMDVol. 143, ASME, New York (1992) 7–24.
[102] Tezduyar TE, Aliabadi S and Behr M. Parallel finite element computing methods for unsteady flows with interfaces. Computational Fluid Dynamics Review 1998 (eds. M. Hafez and K. Oshima), World Scientific (1998) 643–667.
[103] Wehausen JV. The wave resistance of ships. Advances in Applied Mechanics, 1973.
[104] Wigley. Comparative experiment on Wigley parabolic models in Japan. 17th ITTC Resistance Committee Report, 2nd ed., 1983.
[105] Wilcox DC. Turbulence modeling for CFD. DCW Industries Inc., 1994.
[106] Xia F. Numerical calculation of ship flows with special emphasis on the freesurface potential flow. Doctoral Thesis, Chalmers University of Technology, Sweden, 1986.
[107] Zienkiewicz OC and Taylor RC. The finite element method, 5th Edition, 3 Volumes, Butterworth–Heinemann, 2000.
[108] Zienkiewicz OC and Codina R. A general algorithm for compressible and incompressible flow. Part I: The split characteristic based scheme. Int. J. Num. Meth. in Fluids 1995; 20:869–85.
[109] Nakos DE and Sclavounos PD. On steady and unsteady ship wave patterns. J. of Fluid Mechanics 1990; 215:256–288.
Figure 1. Equilibrium of fluxes in a balance domain of finite size
Figure 2. Reference surface for the wave height
Figure 3. Definition of freesurface parameters
Figure 4. Changes on the fluid interface in a floating body
Figure 5. Submerged NACA0012 profile. a) Detail of the mesh of 70000 linear tetrahedra chosen. b) Pressure contours. c) Stationary wave profile
Figure 6. Wigley hull. a) Pressure distribution and mesh deformation of the wigley hull (free model). b) Numerical and experimental body wave profiles. c) freesurface contours for the truly free ship motion
Figure 7. KVLCC2 model. Geometrical definition based on NURBS surfaces and surface mesh used in the analysis
Figure 8. KVLCC2 model. Wave profile on the hull. Thick line shows numerical results
Figure 9. KVLCC2 model. Wave profile on a cut at y/L=0.0964. Thick line shows numerical results
Figure 10. KVLCC2 model. Map of the X component of the velocity on a plane at 2.71 m from the orthogonal aft. Experimental data shown in the right
Figure 11. KVLCC2 model. Map of the eddy kinetic energy () on a plane at 2.71 m from the orthogonal aft. Experimental data shown in the right
Figure 12. NURBSbased geometry used in the analysis of the Rioja de España America Cup boat
Figure 13. Detail of the final (boundary) mesh used in the E0D0 case. The mesh has been adapted taking into account sinkage, trim and freesurface deformation
Figure 14. Detail of the mesh used in the analysis of the E0D0 case, around keelbulb union
Figure 15. Wave elevation profile for 10kn (left: E0D0, right: E15D2)
Figure 16. E0D0 8kn. Pressure contours on bulb
Figure 17. E25D2 8kn. Pressure contours on bulb
Figure 18. E15D2 9kn Wave map
Figure 19. E25D2 9kn Wave map
Figure 20. E15D2 7.5 kn. Pressure map on appendages and streamlines. Perspective view
Figure 21. E15D4 7.5 kn. Velocity modulus contours. Perspective view
Figure 22. E25D2 7.5 kn. Pressure contours on appendages and cuts. Perspective view
Figure 23. E0D0. Resistance graph. Comparison with results extrapolated from experimental data
Figure 24. E15D4. Resistance graph. Comparison with results extrapolated from experimental data
Figure 25. American Cup racing sail boat. a) Mesh used in the analysis. b) Velocity contours
Figure 26. . Streamlines
Figure 27. . Resistance test. Comparison of numerical results with experimental data
Figure 28. Fixed ship hit by incoming wave
Figure 29. Moving ship with fixed velocity
Figure 30. Rotating water mill
Figure 31. Floating wood piece hit by a wave
Figure 32. Motion of a rigid ship hit by an income wave. The ship is modelled as a rigid solid restrained to move in the vertical direction
Figure 33. Ship with stabilizers hit by a lateral wave
Figure 34. Tank ship carrying an internal liquid hit by wave. Ship and fluids motion at different time steps
Figure 35. Tank ship under lateral wave. Pressure distribution at two time steps.
Figure 36. Solid cube falling into a recipient with water. The cube is modelled as a very viscous fluid
Figure 37. Cube falling into a recipient with water. The cube is modelled as a rigid solid. Motion of the cube and free surface positions at different time steps
Figure 38. Cube falling in a water recipient. The cube is modeled as a rigid solid. The finite element meshes generated at the selected instants are shown
Figure 1: Equilibrium of fluxes in a balance domain of finite size 
Figure 2: Reference surface for the wave height 
Figure 3: Definition of freesurface parameters 
Figure 4: Changes on the fluid interface in a floating body 
a)  
b)  
c)  
Figure 5: Submerged NACA0012 profile. a) Detail of the mesh of 70000 linear tetrahedra chosen. b) Pressure contours. c) Stationary wave profile 
Figure 7: KVLCC2 model. Geometrical definition based on NURBS surfaces and surface mesh used in the analysis 
Figure 8: KVLCC2 model. Wave profile on the hull. Thick line shows numerical results 
Figure 9: KVLCC2 model. Wave profile on a cut at y/L=0.0964. Thick line shows numerical results 
Figure 10: KVLCC2 model. Map of the X component of the velocity on a plane at 2.71 m from the orthogonal aft. Experimental data shown in the right 
Figure 11: KVLCC2 model. Map of the eddy kinetic energy () on a plane at 2.71 m from the orthogonal aft. Experimental data shown in the right 
Figure 12: NURBSbased geometry used in the analysis of the Rioja de España America Cup boat 
Figure 13: Detail of the final (boundary) mesh used in the E0D0 case. The mesh has been adapted taking into account sinkage, trim and freesurface deformation 
Figure 14: Detail of the mesh used in the analysis of the E0D0 case, around keelbulb union 
Figure 15: Wave elevation profile for 10kn (left: E0D0, right: E15D2) 
Figure 16: E0D0 8kn. Pressure contours on bulb 
Figure 17: E25D2 8kn. Pressure contours on bulb 
Figure 18: E15D2 9kn Wave map 
Figure 19: E25D2 9kn Wave map 
Figure 20: E15D2 7.5 kn. Pressure map on appendages and streamlines. Perspective view 
Figure 21:E15D4 7.5 kn. Velocity modulus contours. Perspective view 
Figure 22:E25D2 7.5 kn. Pressure contours on appendages and cuts. Perspective view 
Figure 23:E0D0. Resistance graph. Comparison with results extrapolated from experimental data 
Figure 24:E15D4. Resistance graph. Comparison with results extrapolated from experimental data 
a)  b) 
Figure 25:American Cup Bravo España racing sail boat. a) Mesh used in the analysis. b) Velocity contours 
Figure 26:Bravo España. Streamlines 
Figure 27:Bravo España. Resistance test. Comparison of numerical results with experimental data 
Figure 28: Fixed ship hit by incoming wave 
Figure 29: Moving ship with fixed velocity 
Figure 30:Rotating water mill. 
Figure 31:Floating wood piece hit by a wave 
Figure 32:Motion of a rigid ship hit by an incoming wave. The ship is modelled as a rigid solid restrained to move in the vertical direction. 
Figure 33:Ship with stabilizers hit by a lateral wave 
Figure 34:Tanks hip carrying an internal liquid hit by wave. Ship and fluids motion at different time steps 
Figure 34:Cont. 
Figure 35:Tanks hip under lateral wave. Pressure distribution at two time steps 
Figure 36:Solid cube falling into a recipient with water. The cube is modelled as a very viscous fluid. 
Figure 37:Cube falling into a recipient with water. The cube is modelled as a rigid solid. Motion of the cube and free surface positions at different time steps. 
Figure 38:Cube falling in a water recipient. The cube is modeled as a rigid solid. The finite element meshes generated at the selected instants are shown. 
Figure 38:cont 