A computational model for the evaluation of the spray generation of a Wave Adaptive Modular Vessel

==ABSTRACT==

This paper presents part of the work done within the project ‘Advanced Numerical Simulation and Performance Evaluation of WAM-V ® in Spray Generating Conditions’ developed by the International Center for Numerical Methods in Engineering (CIMNE) under Navy Grant N62909-12-1-7101 issued by the Office of Naval Research Global.

One of the primary goals of that project was the development of a computational model for simulation of the Wave Adaptive Modular Vessel (WAM-V®) under spray generating conditions.

For this purpose, a Semi-Lagrangian Particle Finite Element Method (SL-PFEM) has been applied. This is the latest development within the framework of the so-called Particle Finite Element Method (PFEM), using the X-IVAS (eXplicit Integration along the Velocity and Acceleration Streamlines) scheme.

In this paper we demonstrate the applicability of the SL-PFEM using the X-IVAS scheme for the simulation of the Wave Adaptive Modular Vehicle under spray generating conditions.

==INTRODUCTION==

A Wave Adaptive Modular Vessel (WAM-V®) is a new class of ship that uses inflatable flexible hulls to conform to the surface of the water. It is similar in design to a catamaran, in that it has a twin hull design and no keel. However, the superstructure is not rigidly attached to the hulls; it uses shock absorbers and ball joints to articulate the vessel, which allows WAM-V to conform to the surface of the water while mitigating the stresses transmitted to the structure. Moreover, the inflatable hulls help to absorb the high frequency wave-loads. These features allow WAM-V to travel efficiently with low wave resistance in rough seas, by surfing on top of the waves rather than cut through them.

The objective of the WAM-V is to be a lightweight watercraft capable of moving fast and efficiently on the surface of the sea. WAM-Vs are designed to allow for a variety of applications for either manned or unmanned operations and can be built in different lengths to match specific services.

This paper presents part of the work done in the project ‘Advanced Numerical Simulation and Performance Evaluation of WAM-V ® in Spray Generating Conditions’ developed by the International Center for Numerical Methods in Engineering (CIMNE) under Navy Grant N62909-12-1-7101 issued by the Office of Naval Research Global. The scope of that project included the performance analysis of the WAM-V in waves, taking into account the flexibility of the ship hulls, using fluid-structure interaction computational models (see Figure 1). However, the focus of this paper is one of the primary concerns of that project; the development of a computational model for simulation of the WAM-V under spray generating conditions. In this regards, the final goal was to develop and demonstrate a computational engineering solver that could be used to design strategies to reduce the spray generation of the vessel.

{| style="width: 100%;border-collapse: collapse;" 
|-
|  style="vertical-align: top;"|
[[Image:draft_García-Espinosa_900443996-image1-c.jpeg|center|600px]]

|-
|  style="text-align: center;vertical-align: top;"|<span id='_Ref451935418'></span>

<span style="text-align: center; font-size: 75%;">'''Figure 1.''' Snapshot of a fluid-structure interaction analysis of the WAM-V in irregular sea (colormap shows free surface elevation).</span>
|}

When the rest of the components of drag are significantly reduced, the viscous components and any other source of energy dissipation induced by the movement of the vessel, become increasingly important. Because of the inflatable nature of the hulls of the WAM-V, there is little room for hydrodynamic shape optimization. Furthermore, as a consequence of the shape of the hulls, it is likely that they will generate spray when touching the sea surface. Therefore, spray might become an important source of energy dissipation in these little optimized hull shapes. In addition, excessive spray generation can increase the difficulties associated with operating the ship in certain cases (e.g. the possibility of the spray reaching the deck of the craft is a design issue depending on the particular operations the vessel is set to perform).

Therefore it becomes obvious the need to characterize and reduce the spray generation, in order to increase the range of operation of this class of vessels. Furthermore, there is a need to understand the dynamics of the vessel and the hulls in different sea states and the generation of spray when sailing in a seaway.

The Particle Finite Element Method (PFEM, Idelsohn et al., 2004) is a versatile framework for the analysis of fluid-structure interaction problems. The PFEM combines Lagrangian particle-based techniques with the advantage of the integral formulation of the Finite Element Method (FEM).

It has been shown (Idelsohn et al., 2004; Becker, 2015) to successfully simulate a wide variety of complex engineering problems, e.g. free-surface/multi-fluid flows with violent interface motions, multi-fluid mixing and buoyancy-driven segregation problems etc.

The latest development within the framework of the PFEM is the X-IVAS (eXplicit Integration along the Velocity and Acceleration Streamlines) scheme (Idelsohn et al., 2012). It is a semi-implicit scheme built over a Semi-Lagrangian (SL) formulation of the PFEM.

In this paper we present the application of the SL-PFEM using the X-IVAS scheme for the simulation of the Wave Adaptive Modular Vehicle under spray generating conditions.

==SEMI-LAGRANGIAN PARTICLE FINITE ELEMENT METHOD==

''Notation'': Vectors are written using bold italic font and matrices are written using bold upright font.

The independent variables in Lagrangian kinematics are <math display="inline">\lbrace \lambda ,t\rbrace</math> , where <math display="inline">\lambda</math>  represents a label to identify particles and <math display="inline">t</math> represents the time elapsed after labeling.

The primary dependent variable is the fluid particle trajectory denoted as <math display="inline">\mathit{\boldsymbol{X}}\left( \lambda ,t\right)</math> . The independent variables in Eulerian kinematics are <math display="inline">\left( \mathit{\boldsymbol{x}},t\right)</math> , where <math display="inline">\mathit{\boldsymbol{x}}</math> denotes the spatial coordinates. The primary dependent variable is the fluid velocity <math display="inline">\mathit{\boldsymbol{u}}\left( \mathit{\boldsymbol{x}},t\right)</math> .

Consider the Eulerian description of the incompressible Navier-Stokes equations.

{| style="width: 100%;border-collapse: collapse;width: 100%;text-align: center;" class="formulaSCP" 
|-
| 
{| style="text-align: center; margin:auto;width: 100%;"
|-
| <math display="inline">{\partial }_{t}\mathit{\boldsymbol{u+}}\left( \mathit{\boldsymbol{u}}\cdot \nabla \right) \mathit{\boldsymbol{u-}}\nu \Delta \mathit{\boldsymbol{u}}+\nabla (p/\rho )=</math>'''f'''
|}
|  style="text-align: right;width: 5px;text-align: right;white-space: nowrap;"|<span id='_Ref451870648'></span>(1)
|-
| 
{| style="text-align: center; margin:auto;width: 100%;"
|-
| <math>\nabla \cdot \mathit{\boldsymbol{u}}=0</math>
|}
|  style="text-align: right;vertical-align: top;width: 5px;text-align: right;white-space: nowrap;"|<span id='_Ref451870658'></span>(2)
|}

where <math display="inline">\nu</math>  is the kinematic viscosity and <math display="inline">p\left( \mathit{\boldsymbol{x}},t\right)</math> , <math display="inline">\mathit{\boldsymbol{f}}\left( \mathit{\boldsymbol{x}},t\right)</math>  are the pressure and the external acceleration fields, respectively.

The effective acceleration field <math display="inline">\mathit{\boldsymbol{a}}\left( \mathit{\boldsymbol{x}},t\right) \,</math> in the fluid domain is obtained from the momentum balance equation of the flow.

{| style="width: 100%;border-collapse: collapse;width: 100%;text-align: center;" class="formulaSCP" 
|-
| 
{| style="text-align: center; margin:auto;width: 100%;"
|-
| <math display="inline">\mathit{\boldsymbol{a}}={\partial }_{t}\mathit{\boldsymbol{u+}}\left( \mathit{\boldsymbol{u}}\cdot \nabla \right) \mathit{\boldsymbol{u}}=</math><math>\nu \Delta \mathit{\boldsymbol{u}}-\nabla \left( \frac{p}{\rho }\right) +</math>'''f'''
|}
|  style="text-align: right;width: 5px;text-align: right;white-space: nowrap;"|<span id='_Ref451870662'></span>(3)
|}

Note that the functional dependence on the independent variables is suppressed in equations <span id='cite-_Ref451870648'></span>[[#_Ref451870648|(1)]], <span id='cite-_Ref451870658'></span>[[#_Ref451870658|(2)]] and <span id='cite-_Ref451870662'></span>[[#_Ref451870662|(3)]] for brevity.

The fundamental principle of kinematics relates the Eulerian description of the flow with the Lagrangian description as follows.

{| style="width: 100%;border-collapse: collapse;width: 100%;text-align: center;" class="formulaSCP" 
|-
| 
{| style="text-align: center; margin:auto;width: 100%;"
|-
| <math>\mathit{\boldsymbol{U}}\left( \lambda ,t\right) :=\frac{d\mathit{\boldsymbol{X}}\left( \lambda ,t\right) }{dt}=</math><math>\mathit{\boldsymbol{u}}\left( \mathit{\boldsymbol{X}}\left( \lambda ,t\right) ,t\right)</math> 
|}
|  style="text-align: right;width: 5px;text-align: right;white-space: nowrap;"|(4)
|-
| 
{| style="text-align: center; margin:auto;width: 100%;"
|-
| <math>\frac{d\mathit{\boldsymbol{U}}\left( \lambda ,t\right) }{dt}=\frac{{d}^{2}\mathit{\boldsymbol{X}}\left( \lambda ,t\right) }{d{t}^{2}}=</math><math>\mathit{\boldsymbol{a}}\left( \mathit{\boldsymbol{X}}\left( \lambda ,t\right) ,t\right)</math> 
|}
|  style="text-align: right;width: 5px;text-align: right;white-space: nowrap;"|(5)
|}

The basic idea of the X-IVAS scheme is to update the fluid particle position and velocity within a time-step <math display="inline">{t}^{n}\, \leq t\, \leq \, {t}^{n+1}</math> using

{| style="width: 100%;border-collapse: collapse;width: 100%;text-align: center;" class="formulaSCP" 
|-
| 
{| style="text-align: center; margin:auto;width: 100%;"
|-
| <math>\frac{d{\mathit{\boldsymbol{X}}}^{h}\left( \lambda ,t\right) }{dt}={\mathit{\boldsymbol{u}}}^{h}\left( {\mathit{\boldsymbol{X}}}^{\mathit{\boldsymbol{h}}}\left( \lambda ,t\right) ,{t}^{n}\right) =</math><math>\boldsymbol{A}\mathit{\boldsymbol{\, }}{\mathit{\boldsymbol{X}}}^{h}\left( \lambda ,t\right) +</math><math>\mathit{\boldsymbol{b}}</math>
|}
|  style="text-align: right;width: 5px;text-align: right;white-space: nowrap;"|<span id='_Ref451872392'></span>(6)
|-
| 
{| style="text-align: center; margin:auto;width: 100%;"
|-
| <math>\frac{d{\mathit{\boldsymbol{U}}}^{h}\left( \lambda ,t\right) }{dt}={\mathit{\boldsymbol{a}}}^{h}\left( {\mathit{\boldsymbol{X}}}^{h}\left( \lambda ,t\right) ,{t}^{n}\right) =</math><math>\boldsymbol{C}\mathit{\boldsymbol{\, }}{\mathit{\boldsymbol{X}}}^{h}\left( \lambda ,t\right) +</math><math>\mathit{\boldsymbol{d}}</math>
|}
|  style="text-align: right;width: 5px;text-align: right;white-space: nowrap;"|<span id='_Ref451872467'></span>(7)
|}

where <math display="inline">{\mathit{\boldsymbol{u}}}^{h}\left( \mathit{\boldsymbol{x}},t\right)</math>  and <math display="inline">{\mathit{\boldsymbol{a}}}^{h}\left( \mathit{\boldsymbol{x}},t\right)</math>  denote spatially continuous piecewise linear approximations of the velocity and acceleration defined on a background simplicial mesh. The matrices <math display="inline">\boldsymbol{A}</math>, <math display="inline">\boldsymbol{C}</math> and the vectors <math display="inline">\mathit{\boldsymbol{b}}</math>, <math display="inline">\mathit{\boldsymbol{d}}</math> are spatially piecewise constant and depend on the time <math display="inline">{t}^{n}</math>. The particle trajectory and its velocity computed in this manner are denoted as <math display="inline">{\mathit{\boldsymbol{X}}}^{\mathit{\boldsymbol{h}}}\left( \lambda ,t\right)</math>  and <math display="inline">{\mathit{\boldsymbol{U}}}^{\mathit{\boldsymbol{h}}}\left( \lambda ,t\right)</math> , respectively.

Nielson and Jung (1999) presented formulas in 2D and 3D to compute the closed-form analytical solution of tangent curves for piecewise linear vector fields defined over simplicial meshes.

Thus, the Nielson--Jung formulas can be used to compute the analytical solution of <span id='cite-_Ref451872392'></span>[[#_Ref451872392|(6)]]. Idelsohn et al. (2012) presented a procedure to compute the analytical solution of <span id='cite-_Ref451872392'></span>[[#_Ref451872392|(6)]] and <span id='cite-_Ref451872467'></span>[[#_Ref451872467|(7)]] in 2D. However the Nielson--Jung formulas and the calculation procedure described by Idelsohn et al. to compute the analytical solution are not numerically stable; loss of significance occurs due to subtractive cancellations near removable singularities. Recently, Nadukandi (2015) presented numerically stable formulas in 2D and 3D for the closed-form analytical solution of <span id='cite-_Ref451872392'></span>[[#_Ref451872392|(6)]] and <span id='cite-_Ref451872467'></span>[[#_Ref451872467|(7)]].

In the following, we briefly describe the algorithm to implement the SL-PFEM using the X-IVAS scheme.

First the Lagrangian advection of the particles: <math display="inline">{\mathit{\boldsymbol{X}}}^{\mathit{\boldsymbol{h}}}\left( \lambda ,{t}^{n}\right) \, \, \rightarrow</math>  <math display="inline">{\mathit{\boldsymbol{X}}}^{\mathit{\boldsymbol{h}}}\left( \lambda ,{t}^{n+1}\right)</math>  and <math display="inline">{\mathit{\boldsymbol{U}}}^{\mathit{\boldsymbol{h}}}\left( \lambda ,{t}^{n}\right) \,</math>  <math display="inline">\rightarrow \, {\mathit{\boldsymbol{U}}}^{\mathit{\boldsymbol{h}}}\left( \lambda ,{t}^{n+1}\right)</math>  are done solving the following equations

{| style="width: 100%;border-collapse: collapse;width: 100%;text-align: center;" class="formulaSCP" 
|-
| 
{| style="text-align: center; margin:auto;width: 100%;"
|-
| <math>{\mathit{\boldsymbol{X}}}^{h}\left( \lambda ,{t}^{n+1}\right) =={\mathit{\boldsymbol{X}}}^{h}\left( \lambda ,{t}^{n}\right) +</math><math>\int_{{t}^{n}}^{{t}^{n+1}}{\mathit{\boldsymbol{u}}}^{h}\left( {\mathit{\boldsymbol{X}}}^{h}\left( \lambda ,\tau \right) ,{t}^{n}\right) d\tau</math> 
|}
|  style="text-align: right;width: 5px;text-align: right;white-space: nowrap;"|(8)
|}

{| style="width: 100%;border-collapse: collapse;width: 100%;text-align: center;" class="formulaSCP" 
|-
| 
{| style="text-align: center; margin:auto;width: 100%;"
|-
| <math>{\hat{\mathit{\boldsymbol{U}}}}^{h}\left( \lambda ,{t}^{n+1}\right) ={\mathit{\boldsymbol{U}}}^{h}\left( \lambda ,{t}^{n}\right)</math>

<math>-\gamma \int_{{t}^{n}}^{{t}^{n+1}}\nabla \left( \frac{p\left( {\mathit{\boldsymbol{X}}}^{h}\left( \lambda ,\tau \right) ,{t}^{n}\right) }{\rho }\right) d\tau +</math><math>\theta \int_{{t}^{n}}^{{t}^{n+1}}\left( \nu \Delta {\mathit{\boldsymbol{u}}}^{h}\left( {\mathit{\boldsymbol{X}}}^{h}\left( \lambda ,\tau \right) ,{t}^{n}\right) +{\mathit{\boldsymbol{f}}}^{h}\left( {\mathit{\boldsymbol{X}}}^{h}\left( \lambda ,\tau \right) ,{t}^{n}\right) \right) d\tau</math> 
|}
|  style="text-align: right;width: 5px;text-align: right;white-space: nowrap;"|(9)
|}

where <math display="inline">{\hat{\mathit{\boldsymbol{U}}}}^{h}\left( \lambda ,{t}^{n+1}\right)</math>  is an estimate of the velocity of the particle, and <math display="inline">\, \gamma ,\theta \, \in \, \left[ 0,1\right] \,</math>  depend on the integration scheme to be used.

Then the data advected with particles is projected onto a background finite element (FE) mesh. These data include the particle velocities and identities (in multi-fluid flows) among other problem dependent information. Then, the interface between multiple fluids is reconstructed (Becker, 2015) on the FE mesh using the advected particle identities. Appropriate enrichments are determined (Becker, 2015) for the pressure FE shape functions about the interface. The procedure to project the information stored in the particles uses the standard piecewise linear FEM shape-function, <math display="inline">{\, N}_{x}</math>, of a simplicial mesh. Using the shape functions and an integer <math display="inline">\, \mu</math>  usually chosen from the set <math display="inline">\, \, \left\{ 1,2,4\right\}</math> , the implicit approximation <math display="inline">{f}^{h}\left( x\right)</math>  of the projected variables is constructed from samples of the given function <math display="inline">\, f\left( x\right)</math>  at the data point <math display="inline">\left\{ {D}_{i}\right\}</math>  as follows.

{| style="width: 100%;border-collapse: collapse;width: 100%;text-align: center;" class="formulaSCP" 
|-
| 
{| style="text-align: center; margin:auto;width: 100%;"
|-
| <math>W\left( x,{D}_{i}\right) :=\frac{{\left[ {N}_{x}\left( {D}_{i}\right) \right] }^{\mu }}{\sum _{i}^{}{\left[ {N}_{x}\left( {D}_{i}\right) \right] }^{\mu }}\quad</math>

<math>{f}^{h}\left( x\right) =\sum _{i}^{}W(x,\, {D}_{i})f({D}_{i})</math>
|}
|  style="text-align: right;width: 5px;text-align: right;white-space: nowrap;"|<span id='_Ref452495147'></span><span id='_Ref451883182'></span>(10)
|}

Using the particle velocities projected onto the FE mesh as the solution at the start of the time interval <math display="inline">{t}^{n}\, \leq t\, \leq \, {t}^{n+1}</math>, the Stokes problem on the background FE mesh is solved. Using the backward Euler time integration, and a fractional step method (García-Espinosa and Oñate, 2003), the semi-discrete Stokes system to be solved is

{| style="width: 100%;border-collapse: collapse;width: 100%;text-align: center;" class="formulaSCP" 
|-
| 
{| style="text-align: center; margin:auto;width: 100%;"
|-
| <math>{\hat{\mathit{\boldsymbol{u}}}}^{h}\left( \mathit{\boldsymbol{x}},{t}^{n+1}\right) =</math><math>{\wp }^{h}\left( {\hat{\mathit{\boldsymbol{U}}}}^{h}\left( \lambda ,{t}^{n+1}\right) \right)</math> 
|}
|  style="text-align: right;vertical-align: top;width: 5px;text-align: right;white-space: nowrap;"|(11)
|}

{| style="width: 100%;border-collapse: collapse;width: 100%;text-align: center;" class="formulaSCP" 
|-
| 
{| style="text-align: center; margin:auto;width: 100%;"
|-
| <math>\Delta {p}^{h}\left( \mathit{\boldsymbol{x}},{t}^{n+1}\right) =\rho \frac{\nabla \cdot {\hat{\mathit{\boldsymbol{u}}}}^{h}\left( \mathit{\boldsymbol{x}},{t}^{n+1}\right) }{\Delta t}</math>

<math>+\rho \nabla \cdot \left( \nu \Delta {\mathit{\boldsymbol{u}}}^{h}\left( \mathit{\boldsymbol{x}},{t}^{n+1}\right) +\right. </math><math>\left. {\mathit{\boldsymbol{f}}}^{h}\left( \mathit{\boldsymbol{x}},{t}^{n+1}\right) \right) \,</math>

<math>-\gamma \Delta {p}^{h}\left( \mathit{\boldsymbol{x}},{t}^{n}\right)</math>

<math>-\theta \rho \nabla \cdot \left( \nu \Delta {\mathit{\boldsymbol{u}}}^{h}\left( \mathit{\boldsymbol{x}},{t}^{n}\right) +\right. </math><math>\left. {\mathit{\boldsymbol{f}}}^{h}\left( \mathit{\boldsymbol{x}},{t}^{n}\right) \right)</math> 
|}
|  style="text-align: right;vertical-align: top;width: 5px;text-align: right;white-space: nowrap;"|(12)
|}

{| style="width: 100%;border-collapse: collapse;width: 100%;text-align: center;" class="formulaSCP" 
|-
| 
{| style="text-align: center; margin:auto;width: 100%;"
|-
| <math>\frac{{\mathit{\boldsymbol{u}}}^{h}\left( \mathit{\boldsymbol{x}},{t}^{n+1}\right) \mathit{\boldsymbol{-}}{\hat{\mathit{\boldsymbol{u}}}}^{h}\left( \mathit{\boldsymbol{x}},{t}^{n+1}\right) }{\Delta t}=</math><math>-\nabla \left( \frac{{p}^{h}\left( \mathit{\boldsymbol{x}},{t}^{n+1}\right) }{\rho }\right) +</math><math>\left( \nu \Delta {\mathit{\boldsymbol{u}}}^{h}\left( \mathit{\boldsymbol{x}},{t}^{n+1}\right) +\right. </math><math>\left. {\mathit{\boldsymbol{f}}}^{h}\left( \mathit{\boldsymbol{x}},{t}^{n+1}\right) \right) -</math><math>\gamma \nabla \left( \frac{{p}^{h}\left( \mathit{\boldsymbol{x}},{t}^{n}\right) }{\rho }\right) -</math><math>\theta \left( \nu \Delta {\mathit{\boldsymbol{u}}}^{h}\left( \mathit{\boldsymbol{x}},{t}^{n}\right) +\right. </math><math>\left. {\mathit{\boldsymbol{f}}}^{h}\left( \mathit{\boldsymbol{x}},{t}^{n}\right) \right)</math> 
|}
|  style="text-align: right;vertical-align: top;width: 5px;text-align: right;white-space: nowrap;"|(13)
|}

where <math display="inline">{\wp }^{h}</math> is a projection operator from the particles to the FE mesh, as defined in eq. <span id='cite-_Ref452495147'></span>[[#_Ref452495147|(10]]).

We refer to earlier papers (Idelsohn, 2013, 2012, 2015) on the SL-PFEM for several alternate time integration strategies for the Stokes system to be solved on the background mesh.

Finally the particle velocities are updated by the increment <math display="inline">{\boldsymbol{u}}^{h}\left( \boldsymbol{x},{t}^{n+1}\right) -</math><math>{\hat{\boldsymbol{u}}}^{h}\left( \boldsymbol{x},{t}^{n+1}\right)</math>  evaluated at the particle positions.

{| style="width: 100%;border-collapse: collapse;width: 100%;text-align: center;" class="formulaSCP" 
|-
| 
{| style="text-align: center; margin:auto;width: 100%;"
|-
| <math>{\mathit{\boldsymbol{U}}}^{h}\left( \lambda ,{t}^{n+1}\right) ={\hat{\mathit{\boldsymbol{U}}}}^{h}\left( \lambda ,{t}^{n}\right) +</math><math>{\mathit{\boldsymbol{u}}}^{h}\left( {\mathit{\boldsymbol{X}}}^{h}\left( \lambda ,{t}^{n+1}\right) ,{t}^{n+1}\right) -</math><math>{\hat{\mathit{\boldsymbol{u}}}}^{h}\left( {\mathit{\boldsymbol{X}}}^{h}\left( \lambda ,{t}^{n+1}\right) ,{t}^{n+1}\right)</math> 
|}
|  style="text-align: right;width: 5px;text-align: right;white-space: nowrap;"|(14)
|}

<span id='_Toc426546668'></span>

==MODELING OF SPRAY GENERATION==

Two main reasons justify the relevance of the spray generation phenomenon for the WAM-V craft. First, spray might become a relevant source of energy dissipation in the little optimized hull shapes of the craft. Second, excessive spray generation can increase the difficulties associated with operating the ship in certain cases.

As can be seen in Figure 2, spray generation in the current 33 feet WAM-V design is massive in certain conditions. Furthermore, the tubular shape of the hulls does not help to deflect the spray that can easily reach the deck, and therefore deflectors have to be installed.

{| style="width: 100%;border-collapse: collapse;" 
|-
|  style="vertical-align: top;"|
[[Image:draft_García-Espinosa_900443996-image2.png|center|600px]]

|-
|  style="text-align: center;vertical-align: top;"|<span id='_Ref451874301'></span>

<span style="text-align: center; font-size: 75%;">'''Figure 2.''' Spray generated by the WAM-V hull at 21 knots (Peterson, 2014)</span>
|}

As presented above, in the developed model, the fluid data at any given time is available as discrete samples at the spatial locations occupied by the particles. At any particular instant of time, all elements of the auxiliary background mesh are labeled as water-element, air-element or interface-element. Each particle has an identity which is either air (label: <math display="inline">\, +</math><math>1</math>) or water (label: <math display="inline">\, -1</math>). Within each element and for a given time-step, the particles transfer their identities to the element nodes using equations <span id='cite-_Ref451883182'></span>[[#_Ref451883182|(10)]]. After the assembly of the identities, each node has a value between <math display="inline">\, -</math><math>1</math> and <math display="inline">\, +1</math>. Then, a continuous piecewise linear approximation of the otherwise discrete identity data is obtained on the background mesh. The water-air interface is defined as the piecewise planar surface where the aforesaid approximate identity takes a value 0.

In due course, situations may arise where we will find water particles on the air-side of the interface. These particles compose the spray generated in the simulations (see Figure 3). Should the pockets of water particles be large enough then it creates a situation where there are one or more water-elements surrounded by air-elements. These islands of water-elements are seen as water splash in the simulations which represent violent separation and/or merger of the interface. Naturally, such representation of water spray and its intensity depends on the number of particles chosen in the simulation. Nevertheless, the number of particles that compose the spray is not a representation of the mass of water in the spray. Recall that particles represent material points that carry with them only the ''intrinsic'' properties of the flow. So a smaller number of particles just mean that the spray representation is sampled at a coarser level of detail.

{| style="width: 100%;border-collapse: collapse;" 
|-
|  style="vertical-align: top;"|
[[Image:draft_García-Espinosa_900443996-image3.png|center|600px]]

|-
|  style="text-align: center;vertical-align: top;"|<span id='_Ref451874986'></span>

<span style="text-align: center; font-size: 75%;">'''Figure 3.''' Particles that compose the spray in an air-element. The water particles are shown as blue circles; the air particles are seen as white circles.</span>
|}

It is important to emphasize that other possible physical conditions which may generate spray are not considered in our model. For instance, the viscous action of air motion may separate water particles from the interface or decompose existing water splash into water spray. Additionally, this phenomenon can happen at multiple scales wherein entities which can be classified as water spray at a coarse scale may be classified as water splash at a fine scale (which in turn can be decomposed into spray). This cascade will continue until surface tension forces comes to prominence and protect the integrity of the water droplets. Reproducing such physical conditions is out of the scope of this work.

==ANALYSIS OF THE WAM-V USING THE SL-PFEM==

As stated above, one of the primary concerns of this project is the simulation of the WAM-V under spray generating conditions. This section presents an application example of the SL-PFEM solver, including the modelling of the spray generation of the WAM-V craft, following the procedure presented in the previous sections. In order to reduce the computational effort required, only a symmetric configuration (with forward monochromatic seas) has been studied. Therefore, we just simulate the action of only one hull of the craft, neglecting the possible interaction of the deck. Anyhow, these simulations will allow us to obtain a qualitative understanding of the physical conditions leading to spray generation.

{| style="width: 100%;border-collapse: collapse;" 
|-
|  style="vertical-align: top;"|
[[Image:draft_García-Espinosa_900443996-image4-c.png|center|600px]]

|-
|  style="text-align: center;vertical-align: top;"|<span id='_Ref451883427'></span>

<span style="text-align: center; font-size: 75%;">'''Figure 4.''' The domain dimensions, location of the hull and the waterline</span>
|}

The problem domain is a 3D box with straight walls (see Figure 4). The domain dimensions are: 9 m along the x-axis, 3 m along the y-axis and 2 m along the z-axis. The coordinates (0.0,0.0,-1.0) and (9.0,3.0,1.0) represent two diagonally opposite corners of the domain. The face with coordinates (0.0,0.0,-1.0), (9.0,0.0,-1.0), (9.0,0.0,1.0) and (0.0,0.0,1.0) represents the plane of symmetry of the WAM-V. The CAD geometry of the problem domain is obtained by subtracting the volume occupied by the catamaran-type hull of the WAM-V from the volume of the containing 3D box. The geometry of the hull used in the simulations corresponds to the 14 feet (4.27 meters) USV configuration of the WAM-V. The draft of the hull, i.e. the displacement below the steady waterline is taken as  4 inches.

The domain is discretized by a mesh of 3,080,211 three-node tetrahedral elements. On an average twenty particles (material points that transport intrinsic properties of the fluid) per element were used in the CFD simulations, resulting in approximately 60 million particles.

The auxiliary background mesh is deformed every time step such that its internal boundary always conforms to that of the WAM-V hull. The mesh deformation scheme used is based on a Laplacian solver which is commonly used in the implementation of the Arbitrary Lagrangian–Eulerian formulations (Oñate et al., 2004). This enables us to impose the no-slip velocity boundary conditions at the internal boundary in a straight-forward manner.

Since the interest of this study is focused on the near field, the solution scheme used in this case assumes a flat incident free surface field. Therefore, the water wave motion is generated by imposing the solution of a first-order Stokes wave in a narrow strip of water at the inlet and on the walls of the outlet. The analytical solution of the first-order Stokes wave used here is given by

{| style="width: 100%;border-collapse: collapse;width: 100%;text-align: center;" class="formulaSCP" 
|-
| 
{| style="text-align: center; margin:auto;width: 100%;"
|-
| <math>a=0.05\, m,\, \, k=\frac{4\pi }{3},\, \,</math>

<math>g=9.8\, m{/s}^{2},\, \, \omega =\sqrt{gk}</math>

<math>u\left( x,y,t\right) ==a\mathrm{exp}\,\left( kz\right) \mathrm{cos}\,\left( k\left( x-Ut\right) -\omega t\right) +</math><math>U</math>

<math>v\left( x,y,t\right) =a\mathrm{exp}\,(kz)\mathrm{sin}\,\left( k\left( x-Ut\right) -\omega t\right)</math> 
|}
|  style="text-align: right;width: 5px;text-align: right;white-space: nowrap;"|(15)
|}

where <math display="inline">a</math> is the amplitude of the wave, <math display="inline">k</math> is the angular wavenumber, <math display="inline">g</math> is the acceleration due to gravity and <math display="inline">\omega</math>  is the angular frequency. Further ''U'' is the velocity with which the WAM-V moves relative to water and <math display="inline">u\left( x,y,t\right) ,\, v(x,y,t)</math> represent the spatial velocity components of the water in an inertial reference frame that moves with the WAM-V. The narrow strip at the inlet where the wave velocity is imposed has a width of  0.2m. This periodic velocity condition causes a disturbance which is propagated in the rest of the domain and whose motion is governed by the Navier-Stokes equations. In this case, only the near field of the problem is analyzed, and therefore the standard formulation of the fluid dynamics solver (i.e. not the diffraction-radiation splitting formulation) has been used.

The total physical time of simulation was chosen as 4s. The time evolution was computed using 800 time steps of  0.005s each.

The results of three representative cases are here presented. These correspond to three different speeds of the WAM-V, viz. ''U''=15 knots, ''U''=20 knots and ''U''=25 knots, respectively.


{| style="width: 100%;border-collapse: collapse;" 
|-
| <span style="text-align: center; font-size: 75%;">
[[Image:draft_García-Espinosa_900443996-image5.png|center|312px]]
</span>

<span style="text-align: center; font-size: 75%;">a. Time = 3.47s</span>
| <span style="text-align: center; font-size: 75%;">
[[Image:draft_García-Espinosa_900443996-image6.png|center|312px]]
</span>

<span style="text-align: center; font-size: 75%;">b. Time = 3.52s</span>
|-
| <span style="text-align: center; font-size: 75%;">
[[Image:draft_García-Espinosa_900443996-image7.png|center|312px]]
</span>

<span style="text-align: center; font-size: 75%;">c. Time = 3.57s</span>
| <span style="text-align: center; font-size: 75%;">
[[Image:draft_García-Espinosa_900443996-image8.png|center|312px]]
</span>

<span style="text-align: center; font-size: 75%;">d. Time = 3.62s</span>
|-
| <span style="text-align: center; font-size: 75%;">
[[Image:draft_García-Espinosa_900443996-image9.png|center|312px]]
</span>

<span style="text-align: center; font-size: 75%;">e. Time = 3.67s</span>
| <span style="text-align: center; font-size: 75%;">
[[Image:draft_García-Espinosa_900443996-image10.png|center|312px]]
</span>

<span style="text-align: center; font-size: 75%;">f. Time = 3.72s</span>
|-
| <span style="text-align: center; font-size: 75%;">
[[Image:draft_García-Espinosa_900443996-image11.png|center|312px]]
</span>

<span style="text-align: center; font-size: 75%;">g. Time = 3.77s</span>
| <span style="text-align: center; font-size: 75%;">
[[Image:draft_García-Espinosa_900443996-image12.png|center|312px]]
</span>

<span style="text-align: center; font-size: 75%;">h. Time = 3.82s</span>
|}

<div class="center" style="width: auto; margin-left: auto; margin-right: auto;">
<span style="text-align: center; font-size: 75%;">'''Figure 5.''' Simulation of the spray generated by the 14 ft WAM-V hull at 25 knots using the PFEM/X-IVAS method.</span></div>

<br/>

{| style="width: 100%;border-collapse: collapse;" 
|-
|  style="vertical-align: top;"|
[[Image:draft_García-Espinosa_900443996-image13.jpeg|center|600px]]

|-
|  style="text-align: center;vertical-align: top;"|<span style="text-align: center; font-size: 75%;">'''Figure 6.''' Details of the spray generated by the 14 ft WAM-V hull at 15 knots; YZ plane; Time = 3.66s.</span>
|}

{| style="width: 100%;border-collapse: collapse;" 
|-
|  style="vertical-align: top;"|
[[Image:draft_García-Espinosa_900443996-image14-c.jpeg|center|600px]]

|-
|  style="text-align: center;vertical-align: top;"|<span style="text-align: center; font-size: 75%;">'''Figure 7.''' Details of the spray generated by the 14 ft WAM-V hull at 20 knots; YZ plane; Time = 3.59s</span>
|}

{| style="width: 100%;border-collapse: collapse;" 
|-
|  style="vertical-align: top;"|
[[Image:draft_García-Espinosa_900443996-image15.jpeg|center|600px]]

|-
|  style="text-align: center;vertical-align: top;"|<span style="text-align: center; font-size: 75%;">'''Figure 8.''' Details of the spray generated by the 14 ft WAM-V hull at 25 knots; YZ plane; Time = 3.52s</span>
|}


Next, a video of the spray generated by the 14 ft WAM-V hull at 25 knots is presented.


{| style="font-size:75%; color: #222222; border: 1px solid darkgray; background: #f3f3f3; table-layout: fixed; width:100%;"
|- 
| {{#evt:service=youtube|id=https://youtu.be/d7vzDiugF18|alignment=center}}
|- style="text-align: center;" 
| '''Video 1.''' Spray generated by the 14 ft WAM-V hull at 25 knots. 
|}


As discussed in the previous section, the blue isosurface drawn in the pictures corresponds to the piecewise planar surface where the projected identity of the particles on the nodes takes a value 0. It has to be emphasized that this projection process implies a relevant loss of resolution; the flow is actually solved with twenty particles per element (on average) but the blue surface is constructed only with the weighted values on the nodes of the linear elements. In the zones where the developed flow involves a complex mixture of air and water, the free surface cannot be identified. In that case, the calculated isosurface should be understood only as a reference; below that isosurface we will find water with a relatively small proportion of air bubbles, and above it, we will find an increasing presence of air. In particular, the isosurface shows big bulges below the hull that should not be understood as air pockets, but as volumes where the air-water mixing phenomenon is quite complex.


{| style="width: 100%;border-collapse: collapse;" 
|-
|  style="vertical-align: top;"|
[[Image:García-Espinosa_et_al_2016a_3013_Imagen1.jpg|center|600px]]

|-
|  style="text-align: center;vertical-align: top;"|<span id='_Ref451883427'></span>

<span style="text-align: center; font-size: 75%;">'''Figure 9.''' XY plane view of the spray generated by the 14 ft WAM-V hull at 25 knots and time = 3:57 s.</span>
|}


All the run cases took nearly 24 hours each to perform the computations using a workstation with an Intel®Core<sup>TM</sup>i7 - 3820 CPU and 32 GB RAM. The computations were performed using the four (available) cores thanks to the parallelization capabilities of the implemented solver. The memory requirements of these simulations were nearly 17 GB of RAM.

==CONCLUSIONS==

In the different simulations of the WAM-V, we have seen that the semi-Lagrangian Particle Finite Element method can a suitable tool for the analysis of incompressible flows subjected to challenging physical conditions, e.g. violent interface motions, spray generating conditions, etc.

Due to the Lagrangian treatment of the advective processes and the Lagrangian data storage strategy in the SL-PFEM, the interfaces are accurately tracked.

Further, the computational task associated to advective transport is mutually exclusive and hence scalable on parallel computers.

The paper has presented a proposal to represent the generation of water spray based on the discrepancies between the intrinsic information stored in the particles (which is either air or water) and the projected information into the mesh nodes. Nevertheless, the number of particles that compose the spray in the calculation is not a representation of the mass of water in the spray. However, this result can be an indicator of the intensity of the generated spray, which could be used to design and evaluate strategies to reduce the spray generation of the vessel.

Reproducing the complex physics of the different phenomena involved in the spray generation was out of the scope of this work.

==ACKNOWLEDGEMENTS ==

This study was partially supported by the WAM-V project funded under the Navy Grant N62909-12-1-7101 issued by Office of Naval Research Global, the SAFECON project (ref. 267521, FP7-IDEAS-ERC), the FORECAST project (ref. 664910, H2020-ERC-2014-PoC) and the X-SHEAKS project (ref. ENE2014-59194-C2-1-R). The United States Government has a royalty-free license throughout the world in all copyrightable material contained herein.

Permission to use the image shown in Figure 2 has been granted by Prof. Mehdi Ahmadian, VirginiaTech, USA. This image has appeared earlier in Andrew William Peterson’s Ph.D. thesis (2014), figure 3.12, page 55.

==REFERENCES==

Becker, P., Idelsohn, S.R., Oñate, E.: A unified monolithic approach for multi-fluid flows and fluid–structure interaction using the Particle Finite Element Method with fixed mesh. Computational Mechanics, Vol. 55, Issue 6, June 2015, pp. 1091-1104. DOI 10.1007/s00466-014-1107-0.

Becker, P. An enhanced Particle Finite Element Method with special emphasis on landslides and debris flows. Ph.D. thesis, Barcelona Tech (2015).

Celledoni, E., Kometa, B.K., Verdier, O.: High Order Semi-Lagrangian Methods for the Incompressible Navier–Stokes Equations. Journal of Scientific Computing, Vol. 66, Issue 1, Jan. 2016, pp. 91-115. DOI 10.1007/s10915-015-0015-6

Courant, R., Friedrichs, K., Lewy, H.: On the Partial Difference Equations of Mathematical Physics. IBM Journal of Research and Development, Vol. 11, No. 2, 1967, pp. 215–234. DOI 10.1147/rd.112.0215.

Courant, R., Isaacson, E., Rees, M.: On the solution of nonlinear hyperbolic differential equations by finite differences. Communications on Pure and Applied Mathematics Vol. 5, No. 3, 1952, pp. 243–255. DOI 10.1002/cpa.3160050303.

Dadvand, P., Rossi, R., Oñate, E.: An Object-Oriented Environment for Developing Finite Element Codes for Multi-disciplinary Applications. Archives of Computational Methods in Engineering, Vol. 17, No. 3, 2010, pp. 253–297. DOI 10.1007/s11831-010-9045-2.

Dupont, T.F., Liu, Y.: Back and forth error compensation and correction methods for removing errors induced by uneven gradients of the level set function. Journal of Computational Physics, Vol. 190, No. 1, 2003, pp. 311–324. DOI 10.1016/S0021-9991(03)00276-6.

Dupont, T.F., Liu, Y.: Back and forth error compensation and correction methods for semi-lagrangian schemes with application to level set interface computations. Mathematics of Computation. Vol. 76, No. 258, 2007, pp. 647–669. DOI 10.1090/S0025-5718-06-01898-9.

Evans, M.W., Harlow, F.H.: “The Particle-in-Cell Method for Hydrodynamic Calculations”. LA-2139 Nov. 1957, Los Alamos National Laboratory, Los Alamos, New Mexico.

García-Espinosa, J., Oñate, E., Serván-Camas, B., Nadukandi, P. and Becker, P.A. “Advanced Numerical Simulation and Performance Evaluation of WAM-V ® in Spray Generating Conditions. Final Report. Navy Grant N62909-12-1-7101.” July 2015. CIMNE, Barcelona, Spain.

Gelet, R.M., Nguyen, G., Rognon, P.: Modelling interaction of incompressible fluids and deformable particles with the Material Point Method. In: The 6th International Conference on Computational Methods (ICCM2015) (2015)

García-Espinosa, J., Onate, E.: An unstructured finite element solver for ship hydrodynamics problems. Journal of Applied Mechanics Vol. 70, No. 1, 2003, pp. 18-26.

Harlow, F.H.: Hydrodynamic Problems Involving Large Fluid Distortions. Journal of the ACM, Vol. 4, No 2, 1957, pp. 137–142. DOI 10.1145/ 320868.320871.

Idelsohn, S.R., Marti, J., Becker, P., Oñate, E.: Analysis of multifluid flows with large time steps using the particle finite element method. International Journal for Numerical Methods in Fluids, Vol. 75, No 9, 2014, pp. 621–644.

Idelsohn, S.R., Nigro, N., Limache, A., Oñate, E.: Large time-step explicit integration method for solving problems with dominant convection. Computer Methods in Applied Mechanics and Engineering No. 217-220, 2012, pp. 168–185. DOI 10.1016/j.cma.2011.12.008.

Idelsohn, S.R., Nigro, N.M., Gimenez, J.M., Rossi, R., Marti, J.M.: A fast and accurate method to solve the incompressible Navier-Stokes equations. Engineering Computations Vol. 30, No. 2, 2013, pp. 197–222. DOI 10.1108/02644401311304854.

Idelsohn, S.R., Oñate, E., Del Pin, F.: The particle finite element method: a powerful tool to solve incompressible flows with freesurfaces and breaking waves. International Journal for Numerical Methods in Engineering, Vol. 61, No. 7, 2004, pp. 964–989. DOI 10.1002/nme.1096.

Idelsohn, S.R., Oñate, E., Nigro, N., Becker, P., Gimenez, J.: Lagrangian versus Eulerian integration errors. Computer Methods in Applied Mechanics and Engineering Vol. 293, 2015, pp. 191–206. DOI 10.1016/j.cma.2015.04.003.

MacCormack, R.W.: The Effect of Viscosity in Hypervelocity Impact Cratering. Journal of Spacecraft and Rockets Vol. 40, No. 5, 2003, pp. 757–763. DOI 10.2514/2.6901.

Min, C., Gibou, F.: A second order accurate projection method for the incompressible Navier-Stokes equations on non-graded adaptive grids. Journal of Computational Physics 219(2), 912–929 (2006). DOI 10.1016/j.jcp.2006.07.

Nadukandi, P.: Numerically stable formulas for a particle-based explicit exponential integrator. Computational Mechanics, Vol. 55, No. 5, 2015, pp. 903–920. DOI 10.1007/s00466-015-1142-5.

Nadukandi, P., Serván-Camas, B., Becker, P.A. and García-Espinosa, J., Seakeeping with the semi-Lagrangian Particle Finite Element Method. Published online in Computational Particle Mechanics (2016). DOI 10.1007/s40571-016-0127-2

Nielson, G.M., Jung, I.H.: Tools for computing tangent curves for linearly varying vector fields over tetrahedral domains. IEEE Transactions on Visualization and Computer Graphics, Vol. 5, No. 4, 1999, pp. 360–372. DOI 10.1109/2945.817352.

Onate, E., García-Espinosa, J., Idelsohn, S. R. “Ships Hydrodynamics”. In Stein, de Borst and Hughes (eds.), Encyclopedia of computational mechanics. John Wiley & Sons, 2004. DOI: 10.1002/ 0470091355.ecm070

Peterson, A.W.: Simulation and Testing ofWave-Adaptive Modular Vessels. Ph.D. thesis, Virginia Polytechnic Institute and State University (2014).

Robert, A.: A stable numerical integration scheme for the primitive meteorological equations. Atmosphere-Ocean, Vol. 19, No. 1, 1981, pp. 35-46. DOI 10.1080/07055900.1981.9649098.

Sawyer, J.S.: A semi-Lagrangian method of solving the vorticity advection equation. Tellus, Vol. 15, No. 4, 1963, pp. 336–342. DOI 10.1111/j.2153-3490.1963.tb01396.x.

Selle, A., Fedkiw, R., Kim, B., Liu, Y., Rossignac, J.: An Unconditionally Stable MacCormack Method. Journal of Scientific Computing. Vol. 35, No. 2-3, 2008, pp. 350–371. DOI 10.1007/s10915-007-9166-4.

Sulsky, D., Zhou, S.J., Schreyer, H.L.: Application of a particle-in-cell method to solid mechanics. Computer Physics Communications. Vol. 87, No. 1-2, 1995, pp. 236–252. DOI 10.1016/0010-4655(94)00170-7.

WAM-V: The Wave Adaptive Modular Vessel. URL [http://www.wam-v.com http://www.wam-v.com]

Zhang, D.Z., Zou, Q.,VanderHeyden,W.B., Ma, X.: Material point method applied to multiphase flows. Journal of Computational Physics, Vol. 227, No. 6, 2008, pp. 3159–3173. DOI 10.1016/ j.jcp.2007.11.021.