A FEM fluid-structure interaction algorithm for analysis of the seal dynamics of a Surface-Effect Ship

<!-- metadata commented in wiki content
dvdffdfddf

<span id='_Toc337283010'></span><span id='_Toc337283180'></span><span id='_Toc337283553'></span><div id="_Toc337479082" class="center" style="width: auto; margin-left: auto; margin-right: auto;">
<big>'''A FEM fluid-structure interaction algorithm for<br/>analysis of the seal dynamics of a Surface-Effect Ship '''</big></div>

<div class="center" style="width: auto; margin-left: auto; margin-right: auto;">
'''Julio García-Espinosa<sup>1,2</sup>, Daniel Di Capua<sup>1,2</sup>, Borja Serván-Camas<sup>1,*</sup>,<br/>Pere-Andreu Ubach<sup>1</sup> and Eugenio Oñate<sup>1,2</sup>'''</div>
fdfdfdfdfdfd
<div class="center" style="width: auto; margin-left: auto; margin-right: auto;">
<sup>1</sup> Centre Internacional de Mètodes Numèrics en Enginyeria (CIMNE)</div>
dccdcdcdcdcdcd
<div class="center" style="width: auto; margin-left: auto; margin-right: auto;">
Gran Capitán s/n, 08034 Barcelona, Spain</div>

<div class="center" style="width: auto; margin-left: auto; margin-right: auto;">
e-mail : [mailto:julio@cimne.upc.edu julio@cimne.upc.edu], web page: [http://www.cimne.upc.edu http://www.cimne.upc.edu]</div>
dffffffffffffffffffffffffffffffffff
<div class="center" style="width: auto; margin-left: auto; margin-right: auto;">
<sup>2</sup> Universitat Politècnica de Catalunya - BarcelonaTech (UPC)</div>

<div class="center" style="width: auto; margin-left: auto; margin-right: auto;">
Campus Nàutica, Edif. NT3, C. Escar 6-8, 08039 Barcelona, Spain</div>

<div class="center" style="width: auto; margin-left: auto; margin-right: auto;">
web page: [http://www http://www].upc.edu</div>

<div class="center" style="width: auto; margin-left: auto; margin-right: auto;">
<sup>*</sup> Corresponding author. Tel.: +34 93 405 4627.</div>

<div class="center" style="width: auto; margin-left: auto; margin-right: auto;">
E-mail address: [mailto:bservan@cimne.upc.edu bservan@cimne.upc.edu] (B. Serván-Camas).</div>

<div class="center" style="width: auto; margin-left: auto; margin-right: auto;">
URL: [http://www.cimne.com/ http://www.cimne.com/] (B. Serván-Camas).</div>
-->
ffffffffffffffffffffffff
==ABSTRACT==
ffgggggggggggggggggg
This paper shows the recent work of the authors in the development of a time-domain FEM model for evaluation of the seal dynamics of a surface effect ship. The fluid solver developed for this purpose, uses a potential flow approach along with a stream-line integration 
of the free surface. The paper focuses on the free surface-structure algorithm that has been developed to allow the simulation of the complex and highly dynamic behavior of the seals in the interface between the air cushion, and the water.

The developed fluid-structure interaction solver is based, on one side, on an implicit iteration algorithm, communicating pressure forces and displacements of the seals at memory level and, on the other side, on an innovative wetting and drying scheme able to predict the water action on the seals. The stability of the iterative scheme is improved by means of relaxation, and the convergence is accelerated using Aitken’s method.

Several validations against experimental results have been carried out to demonstrate the developed algorithm.

'''Keywords:''' surface effect ship; seal dynamics; finite element method; fluid-structure interaction.

==1. Introduction==

The Transformable Craft (T-Craft) is a novel ship concept of the US Office of Naval Research, operative in multiple modes. T-Craft can deploy in an unloaded condition from the intermediate support base to the seabase, and then be used as a high speed connector to the shore, transporting wheeled and tracked vehicles through the surf zone and onto the beach.

T-Craft has been conceived as a Surface-Effect Ship (SES). A SES is a non-amphibious vehicle supported by an air cushion, with flexible seals at the bow and stern, and twin hulls, like a catamaran, at the sides. Due to the lack of air leakage at the craft sides, lift power can be reduced significantly compared with other type of Air-Cushion Vehicles (ACV). Also, it is possible to install conventional water propellers or waterjet propulsion, with rather smaller machinery space requirements compared to that for air propellers or fans used on ACVs. Furthermore, the SES can operate in modes of full displacement, partial air-cushion support, and full aircushion support.

Predicting the overall performance of a SES is of paramount importance to support the design phase, as the motion of the ship can be affected by the interaction between the air, the cushion, the ship structure, the seals, the sea waves and the sea bottom in the shallow water region. Different approaches with different types of complexity and accuracy have been taken to cope with this type of analyses.

<span id='_Ref325361953'></span>In the last decade, there have been extensive applications of Navier-Stokes models to naval hydrodynamics problems. For example, Oñate and García-Espinosa [<span id='cite-1'></span>[[#1|1]]] presented a stabilized FEM for fluid structure interaction with free surface. In [<span id='cite-2'></span>[[#2|2]]] Löhner et al. developed a FEM capable of tracking violent free surface flows interacting with objects. Also García-Espinosa et al. [<span id='cite-3'></span>[[#3|3]]] developed a new technique to track complex free surface shapes. More recently, in [<span id='cite-4'></span>[[#4|4]]], an application for the calculation of the flow about a SES in still water, using a commercial Volume of Fluid model, has been presented. While, in [<span id='cite-5'></span>[[#5|5]]], Mousaviraad et al. uses an URANS solver for evaluating the manoeuvring performance of a SES. While the outcome of the analyses is outstanding, the CPU-time reported in this paper, makes this model quite unaffordable for being used during design stages.

<span id='_Ref345086558'></span><span id='_Ref346118523'></span>Actually, it is a common consensus that solvers based on the Navier-Stokes equations are too expensive computationally speaking when it comes to simulate unsteady naval hydrodynamics problems. These sorts of problems can be more efficiently calculated using potential flow theory. This approach, jointly with the Stokes perturbation approximation, is widely used for analysis of seakeeping problems [<span id='cite-6'></span>[[#6|6]]]. In [<span id='cite-7'></span>[[#7|7]]], Connell et al., uses a boundary-element time-domain potential flow solver to calculate the multi-body seakeeping behaviour of a T-Craft SES and a LMSR in different scenarios. While, in [<span id='cite-8'></span>[[#8|8]]], the same computational solver is adapted to calculate the manoeuvre of a SES.

<span id='_Ref346124508'></span>Despite the complexity of the above referred SES computational models, none of them takes into account the seal dynamics, or the effect of free surface-seal interaction. However, it is well known the relevance of this interaction in the unsteady dynamics of a SES [<span id='cite-9'></span>[[#9|9]]][<span id='cite-10'></span>[[#10|10]]]. The complexity of this phenomenon makes impossible to develop a theoretical background, and prompts many design parameters to be traditionally decided by empirical formulas [<span id='cite-_Ref346124508'></span>[[#_Ref346124508|9]]]. Actually, only limited theoretical and computational models have been developed to analyze seal dynamics [<span id='cite-11'></span>[[#11|11]]][<span id='cite-12'></span>[[#12|12]]][<span id='cite-13'></span>[[#13|13]]].

This paper shows an extension of the work presented by Serván-Camas and García-Espinosa [<span id='cite-_Ref345086558'></span>[[#_Ref345086558|6]]] in the development of an efficient seakeeping solver. In particular, it is focused in the recent work of the authors in the development of a computational model for the analysis of the complex and highly dynamic behavior of the seals in the interface between the air cushion, and the water of a T-Craft [<span id='cite-14'></span>[[#14|14]]]. The fluid solver developed for this purpose, uses a potential flow approach along with a stream-line integration of the free surface. While this approximation is much simpler than using RANS computations, significant outcomes can be obtained as well, allowing to significantly reducing computational time by 2 or 3 orders of magnitude even when computing on a regular desktop or laptop.

The developed fluid-structure interaction solver is based, on one side, on an implicit iteration algorithm, using a TCP/IP sockets link, able to communicate pressure forces and displacements of the seals at memory level and, on the other side, on an innovative wetting and drying scheme able to predict the water action on the seals.

The outline of this paper is as follows. In section 2, we summarize the mathematical modeling and the details of the integration of the fluid equations. Section 3 describes de structural solver. Section 4 introduces the wetting and drying scheme and the fluid-structure interaction solver. In section 5 different validation cases are presented. Finally, we arrive at some conclusions in section 6.

<span id='_Toc337282423'></span><span id='_Toc337282762'></span><span id='_Toc337283014'></span><span id='_Toc337283184'></span><span id='_Toc337283557'></span><span id='_Toc337479086'></span>

==2. Hydrodynamic model==

===2.1 Governing equations===
xxxxxxx
We consider the first order diffraction-radiation problem of a ship moving with forward speed <math display="inline">\mathit{\boldsymbol{u}}</math>. We assume that the flow is incompressible and irrotational, so that the velocity field <math display="inline">\boldsymbol{v}</math>, is derived from a potential function <math display="inline">\phi</math>  as: <math display="inline">\boldsymbol{\, v=}\nabla \phi</math> . The following assumptions are made on the order of magnitude of the velocity components and free surface elevation <math display="inline">\eta</math> :

{| style="width: 100%;border-collapse: collapse;" 
|-
|  style="text-align: center;"|<math>\mathit{\boldsymbol{u}}=O\left( 1\right) ;\, {\phi }_{x}=O\left( 1\right) ;\, {\phi }_{y}=</math><math>O\left( 1\right) ;\, {\phi }_{z}=O\left( \epsilon \right)</math>

<math>{\phi }_{\alpha \beta }=O\left( \epsilon \right) ;\, \, {\eta }_{x}=O\left( \epsilon \right) ;\, \, {\eta }_{y}=</math><math>O\left( \epsilon \right)</math> 
|  style="text-align: right;"|(1)
|}


Based on the previous assumption, the governing equations for the first order diffraction-radiation wave problem on a frame of reference fixed on the ship are:

{| class="formulaSCP" style="width: 100%;border-collapse: collapse;width: 100%;text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;width: 100%;"
|-
| <math>{\nabla }^{2}\phi =0</math>

<math>{\partial }_{t}\phi -\mathit{\boldsymbol{u}}\cdot {\nabla }_{h}\phi +{{\frac{1}{2}}\nabla }_{h}\phi \cdot {\nabla }_{h}\phi +</math><math>\, \frac{P}{\rho }+g\eta =0</math>

<math>{\partial }_{t}\eta +\left( \mathit{\boldsymbol{-u}}+{\nabla }_{h}\phi \right) \cdot {\nabla }_{h}\eta -</math><math>{\partial }_{z}\phi =0</math>

<math>\nabla \phi \cdot {\mathit{\boldsymbol{n}}}_{b}=-\mathit{\boldsymbol{u}}\cdot {\mathit{\boldsymbol{n}}}_{b}</math>

<math>{\partial }_{z}\phi =0</math>
|}
| 
{| style="text-align: center; margin:auto;width: 100%;"
|-
| in &#x03a9;

in z = 0

in z = 0

in &#x0393;<sub>b</sub>

in ''z = -H''
|}
|  style="text-align: right;width: 5px;text-align: right;white-space: nowrap;"|(1)
|}


where <math display="inline">{\nabla }_{h}\,</math> is the gradient in the horizontal plane, <math display="inline">\Omega \,</math> is the fluid domain, <math display="inline">{\Gamma }_{b}\,</math> represents the wetted surface of the ship,  [[Image:Draft_1_338089811-image1.png|18px]] is the water depth, <math display="inline">{\mathit{\boldsymbol{n}}}_{b}</math> is the normal of the ship wetted surface pointing outwards the ship, and  [[Image:Draft_1_338089811-image2.png|12px]] is the free surface pressure.

===2.2  FEM solver===

The finite element method, FEM, is adopted in this work to solve the system of equations (11-12). The FEM statement of the problem is as follows:

Let <math display="inline">{Q}_{h}^{\ast }</math>  be the finite element space to interpolate functions, constructed in the usual manner. From this space, we can construct the subspace <math display="inline">{Q}_{h,\phi }</math>, that incorporates the Dirichlet conditions for the potential <math display="inline">\, \phi</math> . The space of test functions, denoted by <math display="inline">{\, Q}_{h}</math>, is constructed as <math display="inline">{\, Q}_{h,\phi }</math>, but with functions vanishing on the Dirichlet boundary. The weak form of the problem can be formulated as:

Find <math display="inline">{\, [{\phi }_{h}]\epsilon Q}_{h}</math>, by solving the discrete variational problem:

{| class="formulaSCP" style="width: 100%;border-collapse: collapse;width: 100%;text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;width: 100%;"
|-
| <math>\int_{\Omega }^{}\nabla {\upsilon }_{h}\cdot \nabla {\phi }_{h}d\Omega =</math>

<math>=\int_{{\Gamma }^{B}}^{}{\upsilon }_{h}\cdot {\hat{\phi }}_{n}^{B}d\Gamma +\int_{{\Gamma }^{{Z}_{0}}}^{}{\upsilon }_{h}\cdot {\hat{\phi }}_{n}^{{Z}_{0}}d\Gamma +\int_{{\Gamma }^{{Z}_{-H}}}^{}{\upsilon }_{h}\cdot {\hat{\phi }}_{n}^{{Z}_{-H}}d\Gamma \quad \, \forall {\upsilon }_{h}\in {Q}_{h}</math> 
|}
|  style="text-align: right;width: 5px;text-align: right;white-space: nowrap;"|<span id='_Ref382924947'></span><span id='_Ref382924951'></span>(2)
|}


Where <math display="inline">{\hat{\phi }}_{n}^{B},\, {\hat{\phi }}_{n}^{{Z}_{0}}\,</math> and <math display="inline">{\hat{\phi }}_{n}^{{Z}_{-H}}</math> are the potential normal gradients corresponding to the Neumann boundary conditions on body, radiation boundary, free surface and bottom, respectively. At this point, it is useful to introduce the associated matrix form

{| class="formulaSCP" style="width: 100%;border-collapse: collapse;width: 100%;text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;width: 100%;"
|-
| <math>\overline{\overline{\mathit{\boldsymbol{L}}}}\mathit{\boldsymbol{\phi }}={\mathit{\boldsymbol{b}}}^{B}+</math><math>{\mathit{\boldsymbol{b}}}^{\Gamma o}+{\mathit{\boldsymbol{b}}}^{{\Gamma }_{-H}}</math> 
|}
|  style="text-align: right;width: 5px;text-align: right;white-space: nowrap;"|<span id='_Ref382927175'></span>(3)
|}


Where <math display="inline">\overline{\overline{\mathit{\boldsymbol{L}}}}</math> is the standard laplacian matrix, and <math display="inline">\, {\mathit{\boldsymbol{b}}}^{B}</math>''', ''' <math display="inline">{\mathit{\boldsymbol{b}}}^{\Gamma o}\,</math> and''' ''' <math display="inline">{\mathit{\boldsymbol{b}}}^{{\Gamma }_{-H}}</math>''' '''are the vector resulting of integrating the corresponding boundary condition term. Further details on the FEM solver can be found in [<span id='cite-15'></span>[[#15|15]]].

<span id='_Toc337282425'></span><span id='_Toc337282764'></span><span id='_Toc337283016'></span><span id='_Toc337283186'></span><span id='_Toc337283559'></span><span id='_Toc337479088'></span>

===2.3  Free surface boundary condition===

Solving the free surface boundary condition efficiently is essential when dealing with free surface-structure interaction (FSSI) problems. The free surface conditions can be rewritten as:

{| class="formulaSCP" style="width: 100%;border-collapse: collapse;width: 100%;text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;width: 100%;"
|-
| <math>{\partial }_{t}\phi +\mathit{\boldsymbol{U}}\cdot {\nabla }_{h}\phi -{{\frac{1}{2}}\nabla }_{h}\phi \cdot {\nabla }_{h}\phi +</math><math>\, \frac{P}{\rho }+g\eta =0</math>

<math>{\partial }_{t}\eta +\mathit{\boldsymbol{U}}\cdot {\nabla }_{h}\eta -{\partial }_{z}\phi =</math><math>0</math>
|}
|  style="text-align: right;width: 5px;text-align: right;white-space: nowrap;"|(4)
|}


where <math display="inline">\mathit{\boldsymbol{U=-u+}}{\nabla }_{h}\phi</math>  is the base flow. In this work no linearization is assumed since it is considered that <math display="inline">{\nabla }_{h}\phi =</math><math>O(1)</math>. Therefore, the convective velocity <math display="inline">\mathit{\boldsymbol{U}}</math> will be different in each time step, and the numerical scheme will be adapted accordingly. Three main issues must be kept in mind:

The numerical schemes adopted for solving the kinematic-dynamic free surface boundary conditions are based on an implicit midpoint rule for the kinematic condition, and backward Euler for the dynamic condition. Then <math display="inline">{\phi }^{n+1}</math> can be imposed as a Dirichlet boundary condition in the system <span id='cite-_Ref382924951'></span>[[#_Ref382924951|(2)]]. The schemes read as follows:

{| class="formulaSCP" style="width: 100%;border-collapse: collapse;width: 100%;text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;width: 100%;"
|-
| <math>{\phi }^{n+1}={\phi }^{n}-\Delta t\, {\mathit{\boldsymbol{U}}}^{n+1}\cdot {\left( {\nabla }_{h}\phi \right) }^{n+1}-</math><math>\Delta t\cdot {\left( {\nabla }_{h}\phi \right) }^{n+1}\cdot {\left( {\nabla }_{h}\phi \right) }^{n+1}-</math><math>\Delta t\, \left( {P}^{n+1}/\rho +g{\eta }^{n+1}\right) =0</math>

<math>{\eta }^{n+1}={\eta }^{n}-\Delta t\, {\mathit{\boldsymbol{U}}}^{n+1/2}\cdot {\left( {\nabla }_{h}\eta \right) }^{n+1/2}+</math><math>\Delta t{\left( {\partial }_{z}\phi \right) }^{n+1/2}=0</math>
|}
|  style="text-align: right;width: 5px;text-align: right;white-space: nowrap;"|(5)
|}


It is well known that the convective terms appearing in the above equations can result in solutions corrupted by spurious node-to-node oscillations, if no special care is taken when evaluating them. In this work, the convective term is obtained by differentiating along streamlines:

{| class="formulaSCP" style="width: 100%;border-collapse: collapse;width: 100%;text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;width: 100%;"
|-
| <math>{\mathit{\boldsymbol{U}}}^{n+1}\cdot {\left( {\nabla }_{h}\phi \right) }^{n+1}=\left| {\mathit{\boldsymbol{U}}}^{n+1}\right| {\left( {\partial }_{L}\phi \right) }^{n+1}</math>

<math>{\mathit{\boldsymbol{U}}}^{n}\cdot {\left( {\nabla }_{h}\eta \right) }^{n}=\left| {\mathit{\boldsymbol{U}}}^{n}\right| {\left( {\partial }_{L}\eta \right) }^{n}</math>
|}
|  style="text-align: right;width: 5px;text-align: right;white-space: nowrap;"|(6)
|}


where <math display="inline">{\partial }_{L}</math> denotes the derivative along the streamline. The streamline derivatives are estimated using a three points upstream and one point downstream differential operator inspired by the work of [<span id='cite-16'></span>[[#16|16]]]. The scheme has been conceived to work with unstructured meshes. The use of unstructured meshes enhances geometry flexibility and speed ups the initial modelling time. Furthermore it helps to reduce the number of elements by refining only in those areas of interest. An automatic unstructured grid generator based on the advancing front method is used to generate triangular surface grids and tetrahedral volumetric meshes.

<div class="center" style="width: auto; margin-left: auto; margin-right: auto;">
 [[Image:Draft_1_338089811-image3.jpg|600px]] </div>

<span id='_Ref326165045'></span><span id='_Ref320697593'></span><div id="_GoBack" class="center" style="width: auto; margin-left: auto; margin-right: auto;">
Figure 1: Streamline discretization.</div>

<span id='cite-_Ref326165045'></span>[[#_Ref326165045|Figure 1]] shows the tracing of the streamline at node C. The left (L) and forward left (FL) points are the upstream points, while the right (R) point corresponds to the downstream point. The values of the velocity potential <math display="inline">\phi</math>  and free surface elevation <math display="inline">\eta</math>  at L, FL and R points are obtained by linear interpolation between the nodes of the edges where they lie on. The streamline differential operator reads as:

{| class="formulaSCP" style="width: 100%;border-collapse: collapse;width: 100%;text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;width: 100%;"
|-
| <math>{\partial }_{L}\phi \approx {\delta }_{L}\phi ={\alpha }_{R}{\phi }_{R}+{\alpha }_{C}{\phi }_{C}+</math><math>{\alpha }_{L}{\phi }_{L}+{\alpha }_{FL}{\phi }_{FL}</math>

<math>{\partial }_{L}\eta \approx {\delta }_{L}\eta ={\alpha }_{R}{\eta }_{R}+{\alpha }_{C}{\eta }_{C}+</math><math>{\alpha }_{L}{\eta }_{L}+{\alpha }_{FL}{\eta }_{FL}</math>
|}
|  style="text-align: right;width: 5px;text-align: right;white-space: nowrap;"|(7)
|}


where the stencils are obtained as a one dimensional third order finite difference scheme:

{| class="formulaSCP" style="width: 100%;border-collapse: collapse;width: 100%;text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;width: 100%;"
|-
| <math>\left[ \begin{matrix}\begin{matrix}1&\quad \, \, 1\\{\Delta x}_{R}&\quad \, \, 0\end{matrix}&\begin{matrix}1&1\\{-\Delta x}_{L}&\quad \, \, -{\Delta x}_{L}-{\Delta x}_{FL}\end{matrix}\\\begin{matrix}{\left( {\Delta x}_{R}\right) }^{2}&0\\{\left( {\Delta x}_{R}\right) }^{3}&0\end{matrix}&\quad \begin{matrix}{\left( {\Delta x}_{L}\right) }^{2}&\, {\left( {\Delta x}_{L}+{\Delta x}_{FL}\right) }^{2}\\{-\left( {\Delta x}_{L}\right) }^{3}&-{\left( {\Delta x}_{L}+{\Delta x}_{FL}\right) }^{3}\end{matrix}\end{matrix}\right] \left[ \begin{matrix}\begin{matrix}{\alpha }_{R}\\{\alpha }_{C}\end{matrix}\\\begin{matrix}{\alpha }_{L}\\{\alpha }_{FL}\end{matrix}\end{matrix}\right] =</math><math>\left[ \begin{matrix}\begin{matrix}0\\1\end{matrix}\\\begin{matrix}0\\0\end{matrix}\end{matrix}\right]</math> 
|}
|  style="text-align: right;width: 5px;text-align: right;white-space: nowrap;"|(8)
|}


==3. Structural solver==

For the structural model of the seals it has been considered membrane and shell finite elements.

The membrane finite element used in this work is a 3-noded membrane triangle. This element uses a total Lagrangian formulation and assumes a linear-elastic behaviour of the material [<span id='cite-17'></span>[[#17|17]]]. The kinematic description is formulated in large displacements and large deformations, making use of the right Cauchy-Green deformation tensor.

The shell finite element used in this paper is a corotational 3-noded triangle shell element [<span id='cite-18'></span>[[#18|18]]]. The core formulation of the bending part corresponds to a classical Discrete Kirchhoff plate Triangle (DKT). The core formulation of the membrane part corresponds to the Assumed Natural Deviatoric Strain (ANDES) approach [<span id='cite-19'></span>[[#19|19]]]. The ANDES formulation represents a linear strain triangle with 3 corner nodes and 3 DOFs per node (two in-plane displacements and a drilling rotation).

An energy conserving/decaying time-stepping algorithm has been used for the dynamic analysis of the problem [<span id='cite-20'></span>[[#20|20]]] [<span id='cite-21'></span>[[#21|21]]] [<span id='cite-22'></span>[[#22|22]]]. This approach ensures the stability in the non-linear geometric range.

In general, shell elements are the preferred option to model the seals behaviour, because they can consider the actual flexural rigidity of the seals, despite being small. However, since the thickness of the seals is very small, numerical instabilities can appear. This happens because the tangent matrix is poorly conditioned by the large ratio between membrane and flexural rigidity of the seals. In these cases, the only solution is to neglect the flexural rigidity by using membrane elements. In a dynamic analysis the membrane elements no longer have these problems since the tangent stiffness matrix is properly conditioned by the mass matrix.

<span id='_Toc337282427'></span><span id='_Toc337282766'></span><span id='_Toc337283018'></span><span id='_Toc337283188'></span><span id='_Toc337283561'></span><span id='_Toc337479090'></span>

==4. Coupled fluid-structure solver==

<span id='_Toc337282434'></span><span id='_Toc337282773'></span><span id='_Toc337283025'></span><span id='_Toc337283195'></span><span id='_Toc337283568'></span><span id='_Toc337479097'></span>

===4.1  Free surface-seal interaction (FSSI)===

Modelling the behaviour of the seals of aircushion vehicles is not a trivial issue, due to the complex interaction of the seals with free surface. In this work, a new algorithm for handling the free surface-seal interaction problem is formulated. It is based on finding an equivalent pressure field to be applied over the free surface such that the elevation of this one is limited by the location of the seal. That is to say, the seal act as an upper limit for the free surface elevation.

The free surface boundary conditions are applied in different ways depending on whether the free surface is in contact with the seal or not. If the free surface is not in contact, the boundary conditions are applied as if there is no seal, but if there is contact, the implementation will be different in order to ensure that the free surface does not penetrate the seal and the necessary pressure to fulfil this condition is calculated, and imposed as boundary condition on the free surface.

It will be said that the free surface node where the algorithm is to be applied is dry if the seal is not in contact with the free surface at that location, and wet if it is. <span id='cite-_Ref326165111'></span>[[#_Ref326165111|Figure 2]] illustrates the wet and dry concepts.

The main challenge for an algorithm like this is to be capable of capturing when a node goes from dry to wet and vice versa, as well as estimating the pressure field on the wet nodes. <span id='cite-_Ref326165139'></span>[[#_Ref326165139|Figure 3]] shows the wetting and drying iterative algorithm. For a dry node, the implementation of both, the kinetic and dynamic boundary condition is the same as for any other node not interacting with the seal. However, for wet nodes, the free surface boundary condition is imposed via imposing that the free surface elevation matches the seal elevation, <math display="inline">\, {H}_{seal}</math>, and ensuring that there is no flow across the seal. These two conditions are represented by the following equations:

{| class="formulaSCP" style="width: 100%;border-collapse: collapse;width: 100%;text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;width: 100%;"
|-
| <math>{\eta }_{wet}^{n+1}={H}_{seal}^{n+1}</math>

<math>{\phi }_{z}^{n+1}=2\left( \frac{{\eta }^{n+1}{-\eta }^{n}}{\Delta t}+\, \left| {\mathit{\boldsymbol{U}}}^{n+1/2}\right| {\left( {\delta }_{L}\eta \right) }^{n+1/2}\right) -</math><math>{\phi }_{z}^{n}</math>
|}
|  style="text-align: right;width: 5px;text-align: right;white-space: nowrap;"|(9)
|}


where <math display="inline">{\left( {\partial }_{z}\phi \right) }^{n+1}</math> is to be imposed as Neumann condition instead of imposing <math display="inline">{\phi }^{n+1}</math> as Dirichlet condition.

{| style="width: 100%;" 
|-
|  style="text-align: center;vertical-align: top;width: 100%;"|[[Image:Draft_1_338089811-image4.jpeg|600px]] 
|-
|  style="text-align: center;vertical-align: top;"|<span id='_Ref326165111'></span><span id='_Ref320876970'></span>Figure 2: Wet and dry seal regions for free surface boundary condition implementations.
|}


The change from being a dry node to become a wet node is identified via the kinematic BC though the condition  <math display="inline">{\eta }_{wet}^{n+1}>{H}_{seal}^{n+1}\,</math>  . On the other hand, the switch from being a wet node to become a dry one is carried out by comparing the dynamic pressure with the reference pressure. The dynamic pressure on wet nodes is obtained by the dynamic BC:

{| class="formulaSCP" style="width: 100%;border-collapse: collapse;width: 100%;text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;width: 100%;"
|-
| <math>{P}_{wet}^{n+1}=-\rho \left[ \frac{{\phi }^{n+1}{-\phi }^{n}}{\Delta t}+\, \left| {\boldsymbol{U}}^{n+1}\right| \cdot {\left( {\delta }_{L}\phi \right) }^{n+1}-\right. </math><math>\left. {\frac{1}{2}}{\left( {\nabla }_{h}\phi \right) }^{n+1}\cdot {\left( {\nabla }_{h}\phi \right) }^{n+1}+\right. </math><math>\left. {\frac{1}{2}}{\phi }_{z}^{n+1}{\phi }_{z}^{n+1}+\, g{\eta }^{n+1}\right]</math> 
|}
|  style="text-align: right;width: 5px;text-align: right;white-space: nowrap;"|(10)
|}


We found that when the seal is interacting with the free surface, the term <math display="inline">{\phi }_{z}</math> might become of the same order of magnitude than <math display="inline">\boldsymbol{U}</math>. Hence the term <math display="inline">{\frac{1}{2}}{\phi }_{z}^{2}</math>has been kept when calculating the pressure on a free surface node in contact with the seal.

<div class="center" style="width: auto; margin-left: auto; margin-right: auto;">
 [[Image:Draft_1_338089811-image5.jpg|600px]] </div>

<span id='_Ref326165139'></span><div id="_Ref320876960" class="center" style="width: auto; margin-left: auto; margin-right: auto;">
Figure 3: Wetting and drying algorithm.</div>

<span id='_Toc337282430'></span><span id='_Toc337282769'></span><span id='_Toc337283021'></span><span id='_Toc337283191'></span><span id='_Toc337283564'></span><span id='_Toc337479093'></span>

===4.2  Free surface-flexible seals coupling algorithm===

<span id='_Ref383156396'></span>The fluid-structure coupling is performed by an implicit interfield iteration algorithm [<span id='cite-23'></span>[[#23|23]]], based on a Block Gauss-Seidel method [<span id='cite-24'></span>[[#24|24]]]. For the ''i-''th iteration of the current time step ''n'', let <math display="inline">{\boldsymbol{\phi }}^{n,i}</math>, <math display="inline">{\boldsymbol{\eta }}^{n,i}</math>, <math display="inline">{\boldsymbol{p}}^{n,i}</math> be the vectors of the velocity potential, free surface elevation and pressure in the corresponding nodes of the fluid, <math display="inline">{\boldsymbol{x}}^{n,i}</math> the vector containing the three components of the displacement of the nodes of the seals structure, with <math display="inline">{\boldsymbol{x}}_{k}^{n,i}</math> the vector of the ''k''-th component of the displacement field, and <math display="inline">{\boldsymbol{H}}_{seal}^{n,i-1}</math> the vector of the elevation field of the seal nodes, then the fluid-structure interaction algorithm, can be written as:

{| class="formulaSCP" style="width: 100%;border-collapse: collapse;width: 100%;text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;width: 100%;"
|-
| <math>\begin{matrix}\left( {\boldsymbol{p}}^{n,i},\, {\boldsymbol{\phi }}^{n,i}\right) ={\boldsymbol{F}}^{p}\left( {\boldsymbol{\phi }}^{n,i-1},{\boldsymbol{\eta }}^{n,i-1},{\boldsymbol{H}}_{seal}^{n,i-1},{{\boldsymbol{\eta }}^{n-1,i-1}\boldsymbol{,\phi }}^{n-1}\right) \\{\boldsymbol{x}}^{n,i}={\boldsymbol{S}}^{q}\left( {\boldsymbol{x}}^{n,i-1},{\boldsymbol{p}}^{n,i},{\boldsymbol{x}}^{n-1},{\boldsymbol{x}}^{n-2}\right) \\{\boldsymbol{H}}_{seal}^{n,i-1}={\boldsymbol{H}}_{seal}^{0}+{\boldsymbol{x}}_{k}^{n,i}\end{matrix}</math>
|}
|  style="text-align: right;width: 5px;text-align: right;white-space: nowrap;"|(11)
|}


Where <math display="inline">\boldsymbol{F,\, S}</math>'' '' represent the fluid and solid solver, and ''p'' and ''q'' are the number of internal iterations done within the non-linear loop of each solver. This way, the pressure field computed in the hydrodynamic solver is sent to the structural solver to compute the seals deformation. The resulting displacements are used to compute the new seal elevation on the fluid interface, and then calculate the pressure field for the following iteration.

{| style="width: 100%;border-collapse: collapse;" 
|-
|  style="text-align: center;vertical-align: top;width: 100%;"|[[Image:Draft_1_338089811-image6.jpeg|600px]] 
|-
|  style="text-align: center;vertical-align: top;"|<span id='_Ref381981831'></span><span id='_Ref381981819'></span>Figure 4: Fluid-structure interaction algorithm.
|}


Since the hydrodynamics and structural solvers used in this work are independent, the strategy developed to communicate both is based on the interchange of information at memory level by means of TCP-IP sockets. This way, only minor adaptations of the solvers were required. The required functionality has been implemented in C++, in form of library. The communication protocol developed for the library is quite simple. Actually, it just requires sending and receiving functions for vectorial fields, including an implicit waiting instruction in the receiving function for synchronization purposes. Furthermore, this library has been extended to be able to interpolate data between the boundary meshes of the fluid and solid interfaces (see <span id='cite-_Ref381981831'></span>[[#_Ref381981831|Figure 4]]). The interpolating algorithm is based on an octtree search to find the interpolating element of every node of the other mesh.

Let <math display="inline">{\boldsymbol{I}}^{FS}{\boldsymbol{,\, I}}^{SF}</math> be the interpolator operators from the fluid boundary mesh to the solid interface, and vice-versa, then the fluid-structure interaction algorithm reads

{| class="formulaSCP" style="width: 100%;border-collapse: collapse;width: 100%;text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;width: 100%;"
|-
| <math>\begin{matrix}\left( {\boldsymbol{p}}^{n,i},\, {\boldsymbol{\phi }}^{n,i}\right) ={\boldsymbol{F}}^{p}\left( {\boldsymbol{\phi }}^{n,i-1},{\boldsymbol{\eta }}^{n,i-1},{\boldsymbol{H}}_{seal}^{n,i-1},{{\boldsymbol{\eta }}^{n-1,i-1}\boldsymbol{,\phi }}^{n-1}\right) \\{\check{\boldsymbol{p}}}^{n,i}={\boldsymbol{I}}^{FS}\left( {\boldsymbol{p}}^{n,i}\right) \\{\boldsymbol{x}}^{n,i}={\boldsymbol{S}}^{q}\left( {\boldsymbol{x}}^{n,i-1},{\check{\boldsymbol{p}}}^{n,i},{\boldsymbol{x}}^{n-1},{\boldsymbol{x}}^{n-2}\right) \\{\check{\boldsymbol{x}}}^{n,i}={\boldsymbol{I}}^{SF}\left( {\boldsymbol{x}}^{n,i}\right) \\{\boldsymbol{H}}_{seal}^{n,i-1}={\boldsymbol{H}}_{seal}^{0}+{\check{\boldsymbol{x}}}_{k}^{n,i}\end{matrix}</math>


|}
|  style="text-align: right;width: 5px;text-align: right;white-space: nowrap;"|(12)
|}


Even though the iterations within each time step lead to the monolithic solution of the coupled problem, sometimes relaxation techniques are combined with this method to ensure the convergence of the coupling scheme. In this work, a relaxation method based on Aitken’s acceleration has been adopted [<span id='cite-25'></span>[[#25|25]]]. Let <math display="inline">{\omega }_{i}</math> be the optimal relaxation parameter, obtained via Aitken’s method, the resulting algorithm can be written as:

{| class="formulaSCP" style="width: 100%;border-collapse: collapse;width: 100%;text-align: center;" 
|-
| 
{| style="text-align: center; margin:auto;width: 100%;"
|-
| <math>\begin{matrix}\left( {\boldsymbol{p}}^{n,i},\, {\boldsymbol{\phi }}^{n,i}\right) ={\boldsymbol{F}}^{p}\left( {\boldsymbol{\phi }}^{n,i-1},{\boldsymbol{\eta }}^{n,i-1},{\boldsymbol{H}}_{seal}^{n,i-1},{{\boldsymbol{\eta }}^{n-1,i-1}\boldsymbol{,\phi }}^{n-1}\right) \\{\hat{\boldsymbol{p}}}^{n,i}=\left( 1-{\omega }_{i}\right) {\hat{\boldsymbol{p}}}^{n,i-1}+{\omega }_{i}{\boldsymbol{p}}^{n,i}\\{\check{\boldsymbol{p}}}^{n,i}={\boldsymbol{I}}^{FS}\left( {\hat{\boldsymbol{p}}}^{n,i}\right) \\{\boldsymbol{x}}^{n,i}={\boldsymbol{S}}^{q}\left( {\boldsymbol{x}}^{n,i-1},{\check{\boldsymbol{p}}}^{n,i},{\boldsymbol{x}}^{n-1},{\boldsymbol{x}}^{n-2}\right) \\{\check{\boldsymbol{x}}}^{n,i}={\boldsymbol{I}}^{SF}\left( {\boldsymbol{x}}^{n,i}\right) \\{\boldsymbol{H}}_{seal}^{n,i-1}={\boldsymbol{H}}_{seal}^{0}+{\check{\boldsymbol{x}}}_{k}^{n,i}\end{matrix}</math>


|}
|  style="text-align: right;width: 5px;text-align: right;white-space: nowrap;"|(13)
|}


As stated above, the free surface-seals interaction algorithm is based on finding an equivalent pressure field to be applied over the free surface such that the elevation of this one is limited by the location of the seal. For this purpose, the elevation of the seal over the reference free surface is required. The algorithm <span id='cite-_Ref382927175'></span>[[#_Ref382927175|(3)]] assumes that the variation of this elevation is only due to the vertical displacements of the structure. This is a good enough approximation to predict the forces of the water on the seals and the free surface elevation in the seals region.

<span id='_Toc337282452'></span><span id='_Toc337282791'></span><span id='_Toc337283043'></span><span id='_Toc337283213'></span><span id='_Toc337283586'></span><span id='_Toc337479115'></span>

==5. Validation of FSSI algorithm with flexible bow seals ==

===5.1  Validation against experiments at University of Michigan===

<span id='_Ref382560783'></span><span id='_Ref383159193'></span>Zalek and Doctors [<span id='cite-26'></span>[[#26|26]]] carried out at University of Michigan an extensive experimental study of free surface and flexible seal interaction. Some of this data has been already used for validation purpose against numerical results based on a smoothed particle hydrodynamics approach [<span id='cite-27'></span>[[#27|27]]]. <span id='cite-_Ref382987144'></span>[[#_Ref382987144|Figure 5]] sketches the basic arrangements of the experiments, such as seal length, location of stiffeners and cables and undeflected seal angle. <span id='cite-_Ref382560620'></span>[[#_Ref382560620|Table 1]] provides the seal properties used for the experiments [<span id='cite-_Ref383159193'></span>[[#_Ref383159193|27]]].

<div class="center" style="width: auto; margin-left: auto; margin-right: auto;">
 [[Image:Draft_1_338089811-image7.jpeg|372px]] </div>

<div id="_Ref382987144" class="center" style="width: auto; margin-left: auto; margin-right: auto;">
Figure 5: Experiments setup at University of Michigan</div>

{| style="width: 39%;margin: 1em auto 0.1em auto;border-collapse: collapse;" 
|-
|  style="border: 1pt solid black;text-align: center;"|Thickness
|  style="border: 1pt solid black;text-align: center;"|3.175 mm
|-
|  style="border: 1pt solid black;text-align: center;"|Young Modulus
|  style="border: 1pt solid black;text-align: center;"|12.74 MPa
|-
|  style="border: 1pt solid black;text-align: center;"|Poisson coefficient
|  style="border: 1pt solid black;text-align: center;"|0.25
|-
|  style="border: 1pt solid black;text-align: center;"|Density
|  style="border: 1pt solid black;text-align: center;"|1107 kg/m<sup>3</sup>
|}


<div id="_Ref382560620" class="center" style="width: auto; margin-left: auto; margin-right: auto;">
Table 1: Seal properties</div>

In this work, the experimental cases shown in <span id='cite-_Ref382560623'></span>[[#_Ref382560623|Table 2]] were analysed with the proposed fluid-structure interaction model. The test runs were simulated modelling the left half part of seal and tank, using the actual width and depth of the channel of the University of Michigan. The seals were modelled with shells elements, with the properties provided in <span id='cite-_Ref382560620'></span>[[#_Ref382560620|Table 1]].

{| style="width: 78%;margin: 1em auto 0.1em auto;border-collapse: collapse;" 
|-
|  rowspan='2' style="border: 1pt solid black;text-align: center;"|Test run [<span id='cite-_Ref382560783'></span>[[#_Ref382560783|26]]]
|  colspan='2'  style="border: 1pt solid black;text-align: center;"|Seal Immersion
|  colspan='2'  style="border: 1pt solid black;text-align: center;"|Velocity
|  colspan='2'  style="border: 1pt solid black;text-align: center;"|Internal Pressure
|-
|  style="border: 1pt solid black;text-align: center;"|(in)
|  style="border: 1pt solid black;text-align: center;"|(m)
|  style="border: 1pt solid black;text-align: center;"|(ft/s)
|  style="border: 1pt solid black;text-align: center;"|(m/s)
|  style="border: 1pt solid black;text-align: center;"|(in. H<sub>2</sub>O )
|  style="border: 1pt solid black;text-align: center;"|(Pa)
|-
|  style="border: 1pt solid black;text-align: center;"|1
|  style="border: 1pt solid black;text-align: center;"|9
|  style="border: 1pt solid black;text-align: center;"|0.2286
|  style="border: 1pt solid black;text-align: center;"|6
|  style="border: 1pt solid black;text-align: center;"|1.8288
|  style="border: 1pt solid black;text-align: center;"|1.42
|  style="border: 1pt solid black;text-align: center;"|353.3
|-
|  style="border: 1pt solid black;text-align: center;"|10
|  style="border: 1pt solid black;text-align: center;"|9
|  style="border: 1pt solid black;text-align: center;"|0.2286
|  style="border: 1pt solid black;text-align: center;"|8
|  style="border: 1pt solid black;text-align: center;"|2.4384
|  style="border: 1pt solid black;text-align: center;"|1.01
|  style="border: 1pt solid black;text-align: center;"|251.3
|-
|  style="border: 1pt solid black;text-align: center;"|11
|  style="border: 1pt solid black;text-align: center;"|9
|  style="border: 1pt solid black;text-align: center;"|0.2286
|  style="border: 1pt solid black;text-align: center;"|8
|  style="border: 1pt solid black;text-align: center;"|2.4384
|  style="border: 1pt solid black;text-align: center;"|1.91
|  style="border: 1pt solid black;text-align: center;"|475.3
|-
|  style="border: 1pt solid black;text-align: center;"|16
|  style="border: 1pt solid black;text-align: center;"|9
|  style="border: 1pt solid black;text-align: center;"|0.2286
|  style="border: 1pt solid black;text-align: center;"|8
|  style="border: 1pt solid black;text-align: center;"|2.4384
|  style="border: 1pt solid black;text-align: center;"|4.24
|  style="border: 1pt solid black;text-align: center;"|1055
|-
|  style="border: 1pt solid black;text-align: center;"|32
|  style="border: 1pt solid black;text-align: center;"|7.5
|  style="border: 1pt solid black;text-align: center;"|0.1905
|  style="border: 1pt solid black;text-align: center;"|6
|  style="border: 1pt solid black;text-align: center;"|1.8288
|  style="border: 1pt solid black;text-align: center;"|1.28
|  style="border: 1pt solid black;text-align: center;"|318.5
|-
|  style="border: 1pt solid black;text-align: center;"|45
|  style="border: 1pt solid black;text-align: center;"|7.5
|  style="border: 1pt solid black;text-align: center;"|0.1905
|  style="border: 1pt solid black;text-align: center;"|9
|  style="border: 1pt solid black;text-align: center;"|2.7432
|  style="border: 1pt solid black;text-align: center;"|4.01
|  style="border: 1pt solid black;text-align: center;"|997.8
|}


<div id="_Ref382560623" class="center" style="width: auto; margin-left: auto; margin-right: auto;">
Table 2: Test run particulars</div>

<span id='cite-_Ref382988606'></span>[[#_Ref382988606|Table 3]] provides the main characteristics of the computational mesh used for the simulation, as well as the length of the numerical channel used. For test run 45, a larger channel was necessary to avoid upstream and downstream interactions with the inlet and outlet boundaries respectively.

{| style="width: 100%;border-collapse: collapse;" 
|-
|  rowspan='2' style="border: 1pt solid black;text-align: center;"|Test run
|  colspan='2'  style="border: 1pt solid black;text-align: center;"|Seal
|  colspan='2'  style="border: 1pt solid black;text-align: center;"|Water volume
|  rowspan='2' style="border: 1pt solid black;text-align: center;"|Time step (s)
|  rowspan='2' style="border: 1pt solid black;text-align: center;"|Domain length (m)
|-
|  style="border: 1pt solid black;text-align: center;"|Number of Nodes
|  style="border: 1pt solid black;text-align: center;"|Number of Elements
|  style="border: 1pt solid black;text-align: center;"|Number of Nodes
|  style="border: 1pt solid black;text-align: center;"|Number of Tetrahedra
|-
|  style="border: 1pt solid black;text-align: center;"|1
|  style="border: 1pt solid black;text-align: center;"|234
|  style="border: 1pt solid black;text-align: center;"|424
|  style="border: 1pt solid black;text-align: center;"|58136
|  style="border: 1pt solid black;text-align: center;"|297188
|  style="border: 1pt solid black;text-align: center;"|0.002
|  style="border: 1pt solid black;text-align: center;"|5.3
|-
|  style="border: 1pt solid black;text-align: center;"|10
|  style="border: 1pt solid black;text-align: center;"|234
|  style="border: 1pt solid black;text-align: center;"|424
|  style="border: 1pt solid black;text-align: center;"|58315
|  style="border: 1pt solid black;text-align: center;"|297888
|  style="border: 1pt solid black;text-align: center;"|0.001
|  style="border: 1pt solid black;text-align: center;"|5.3
|-
|  style="border: 1pt solid black;text-align: center;"|11
|  style="border: 1pt solid black;text-align: center;"|234
|  style="border: 1pt solid black;text-align: center;"|424
|  style="border: 1pt solid black;text-align: center;"|58315
|  style="border: 1pt solid black;text-align: center;"|297888
|  style="border: 1pt solid black;text-align: center;"|0.001
|  style="border: 1pt solid black;text-align: center;"|5.3
|-
|  style="border: 1pt solid black;text-align: center;"|16
|  style="border: 1pt solid black;text-align: center;"|234
|  style="border: 1pt solid black;text-align: center;"|424
|  style="border: 1pt solid black;text-align: center;"|58240
|  style="border: 1pt solid black;text-align: center;"|297643
|  style="border: 1pt solid black;text-align: center;"|0.002
|  style="border: 1pt solid black;text-align: center;"|5.3
|-
|  style="border: 1pt solid black;text-align: center;"|32
|  style="border: 1pt solid black;text-align: center;"|234
|  style="border: 1pt solid black;text-align: center;"|424
|  style="border: 1pt solid black;text-align: center;"|56112
|  style="border: 1pt solid black;text-align: center;"|287493
|  style="border: 1pt solid black;text-align: center;"|0.002
|  style="border: 1pt solid black;text-align: center;"|5.3
|-
|  style="border: 1pt solid black;text-align: center;"|45
|  style="border: 1pt solid black;text-align: center;"|234
|  style="border: 1pt solid black;text-align: center;"|424
|  style="border: 1pt solid black;text-align: center;"|80325
|  style="border: 1pt solid black;text-align: center;"|419569
|  style="border: 1pt solid black;text-align: center;"|0.001
|  style="border: 1pt solid black;text-align: center;"|15.3
|}


<div id="_Ref382988606" class="center" style="width: auto; margin-left: auto; margin-right: auto;">
Table 3: Computational domain and numerical particulars</div>

In the experiments, three stiffeners (rods) were fixed to the seal to reduce the transversal deformation of the seals. Since we had no data regarding the properties of such stiffeners, we simulated them by introducing very stiff beams at the same locations, but weightless, assuming that the ones used in the experiments were light enough to avoid significant perturbation in the seal behaviour. Furthermore several cables were attached to the stiffeners rods to avoid the forward displacement of the seal due to the internal pressure. Since there is no available information about the length and exact disposition of these cables, they were simulated in the computational model by a boundary condition, constraining the longitudinal displacement of the points of the seal where the cables were attached.

<span id='cite-_Ref382812431'></span>[[#_Ref382812431|Figure 6]] shows the quasi-steady deformed front seal and free surface elevation respectively for test run number 10. <span id='cite-_Ref382812522'></span>[[#_Ref382812522|Figure 7]] compares de average seal deformation obtained numerically and experimentally. Although the computational results do not match exactly the experimental ones, they agree well in terms of the main trends under different conditions of velocity, seal immersion and internal pressure. Anyhow, the differences are probably within the uncertainty range of the data, since the unsteady and three-dimensional effects are relevant in most of the studied cases.

<div class="center" style="width: auto; margin-left: auto; margin-right: auto;">
 [[Image:Draft_1_338089811-image8.jpeg|522px]] </div>

<span id='_Ref382812431'></span><div id="_Ref382812425" class="center" style="width: auto; margin-left: auto; margin-right: auto;">
Figure 6: Seal deformation and structural mesh for test run 10.</div>

<div class="center" style="width: auto; margin-left: auto; margin-right: auto;">
 [[Image:Draft_1_338089811-image9.jpeg|522px]] </div>

<div id="_Ref382812522" class="center" style="width: auto; margin-left: auto; margin-right: auto;">
Figure 7: Seal deformation: experimental results (squares) versus numerical (solid line).</div>

<span id='cite-_Ref382812527'></span>[[#_Ref382812527|Figure 8]] provides the top view of the hydrodynamic pressure distribution acting on the seal due to the free surface seal interaction. It can be observed that while in some cases the distribution is almost uniform in the transverse direction it is not in others, resulting in relevant three dimensional effects, as mentioned before.

<div class="center" style="width: auto; margin-left: auto; margin-right: auto;">
 [[Image:Draft_1_338089811-image10.jpeg|342px]] </div>

<div id="_Ref382812527" class="center" style="width: auto; margin-left: auto; margin-right: auto;">
Figure 8: Top view of pressure distribution over seal</div>

===5.2  Validation against experiments at the U.S. Navy´s Large Cavitation Channel (LCC)===

Wiggins et al. [<span id='cite-_Ref382813140'></span>[[#_Ref382813140|28]]] developed a large scale testing platform for finger bow seals at the U.S. Navy´s Large Cavitation Channel. The purpose of those test were to provide experimental results to be used for validation purposes of analytical/numerical approaches to the problem. <span id='cite-_Ref383012663'></span>[[#_Ref383012663|Figure 9]] sketches the model setup at the LC facilities

[[Image:Draft_1_338089811-image11.jpeg|600px]]

<div id="_Ref383012663" class="center" style="width: auto; margin-left: auto; margin-right: auto;">
Figure 9: LCC model setup sketch</div>

The seal material used is the same one that the one used in the University of Michigan experiments and described in the previous section. In this case, a membrane model was used for the computational analysis of the structure, since preliminary tests with a shell model showed an unstable behavior, probably due to the poor condition of the tangent matrix.

With regards to the test condition, we used a free stream flow of 6 m/s, seal immersion of 15 cm, and internal pressure of 1.5 kPa [<span id='cite-_Ref382813140'></span>[[#_Ref382813140|28]]]. Other data needed for the setup of the case under study are provided in <span id='cite-_Ref382814840'></span>[[#_Ref382814840|Table 4]].

{| style="width: 53%;margin: 1em auto 0.1em auto;border-collapse: collapse;" 
|-
|  style="border: 1pt solid black;text-align: center;"|Breadth of cushion (m)
|  style="border: 1pt solid black;text-align: center;"|1.52
|-
|  style="border: 1pt solid black;text-align: center;"|Sidewall depth (m)
|  style="border: 1pt solid black;text-align: center;"|1.68
|-
|  style="border: 1pt solid black;text-align: center;"|Sidewall draft at zero immersion (m)
|  style="border: 1pt solid black;text-align: center;"|0.42
|-
|  style="border: 1pt solid black;text-align: center;"|Channel width (m)
|  style="border: 1pt solid black;text-align: center;"|3.05
|-
|  style="border: 1pt solid black;text-align: center;"|Number of Fingers 5
|  style="border: 1pt solid black;text-align: center;"|5
|-
|  style="border: 1pt solid black;text-align: center;"|Finger height (m)
|  style="border: 1pt solid black;text-align: center;"|1.00
|-
|  style="border: 1pt solid black;text-align: center;"|Finger width (m)
|  style="border: 1pt solid black;text-align: center;"|0.31
|-
|  style="border: 1pt solid black;text-align: center;"|Bow seal angle (º)
|  style="border: 1pt solid black;text-align: center;"|50.00
|}


<span id='_Ref382814840'></span><div id="_Ref382814836" class="center" style="width: auto; margin-left: auto; margin-right: auto;">
Table 4: Tests particulars [<span id='cite-_Ref382813140'></span>[[#_Ref382813140|28]]]</div>

<span id='cite-_Ref383007174'></span>[[#_Ref383007174|Table 5]] provides the main characteristics of the computational mesh used for the simulation, as well as the length of the numerical channel used. <span id='cite-_Ref383007645'></span>[[#_Ref383007645|Figure 10]] shows the computational mesh used for solving the fluid problem.

{| style="width: 57%;margin: 1em auto 0.1em auto;border-collapse: collapse;" 
|-
|  rowspan='2' style="border: 1pt solid black;text-align: center;"|Fingers
|  style="border: 1pt solid black;text-align: center;"|Number of Nodes
|  style="border: 1pt solid black;text-align: center;"|7023
|-
|  style="border: 1pt solid black;text-align: center;"|Number of Elements
|  style="border: 1pt solid black;text-align: center;"|13436
|-
|  rowspan='2' style="border: 1pt solid black;text-align: center;"|Water domain
|  style="border: 1pt solid black;text-align: center;"|Number of Nodes
|  style="border: 1pt solid black;text-align: center;"|120942
|-
|  style="border: 1pt solid black;text-align: center;"|Number of Tetrahedra
|  style="border: 1pt solid black;text-align: center;"|659687
|-
|  colspan='2'  style="border: 1pt solid black;text-align: center;"|Time step (s)
|  style="border: 1pt solid black;text-align: center;"|0.0005
|-
|  colspan='2'  style="border: 1pt solid black;text-align: center;"|Domain length (m)
|  style="border: 1pt solid black;text-align: center;"|21
|}


<div id="_Ref383007174" class="center" style="width: auto; margin-left: auto; margin-right: auto;">
Table 5: Computational domain and numerical particulars</div>

<div class="center" style="width: auto; margin-left: auto; margin-right: auto;">
 [[Image:Draft_1_338089811-image12.jpeg|510px]] </div>

<div id="_Ref383007645" class="center" style="width: auto; margin-left: auto; margin-right: auto;">
Figure 10: Computational mesh for potential flow solver</div>

<span id='cite-_Ref337479284'></span>[[#_Ref337479284|Figure 11]] (a) shows the steady seal deformation obtained in the experiments under the abovementioned test conditions. <span id='cite-_Ref337479284'></span>[[#_Ref337479284|Figure 11]] (b), (c) and (d) show seal deformation obtained by the numerical approach proposed in this work under the same test conditions. Free surface elevation at the seal and aircushion areas are shown in <span id='cite-_Ref382813833'></span>[[#_Ref382813833|Figure 12]] (a). A top view of the pressure distribution over finger seals is shown in <span id='cite-_Ref382813833'></span>[[#_Ref382813833|Figure 12]] (b).

<div class="center" style="width: auto; margin-left: auto; margin-right: auto;">
 [[Image:Draft_1_338089811-image13.jpeg|402px]] </div>

<span id='_Ref337479284'></span><div id="_Ref382813140" class="center" style="width: auto; margin-left: auto; margin-right: auto;">
Figure 11: Experimental results (a) versus numerical results (b), (c) and (d).  Picture obtained in the experimental tests reported in [<span id='cite-28'></span>[[#28|28]]]  (a); Front view of deformed fingers (b); deformed finger detail (c); top view (d).</div>

<div class="center" style="width: auto; margin-left: auto; margin-right: auto;">
<big>''' [[Image:Draft_1_338089811-image14.jpeg|600px]] '''</big></div>

<div id="_Ref382813833" class="center" style="width: auto; margin-left: auto; margin-right: auto;">
Figure 12: Numerical results: free surface elevation (a) and top view of pressure distribution on fingers (b)</div>

As can be seen, the bow seals deformation obtained in this analysis agrees with the qualitative experimental information available in reference [<span id='cite-_Ref382813140'></span>[[#_Ref382813140|28]]]. In fact, the deformation mode obtained in the computational model is quite similar to that found experimentally.

==6. Conclusions==

This paper reports the advances in the development of a fluid-structure interaction algorithm for evaluation of a SES seal dynamics.

The fluid solver developed for this purpose, uses a potential flow approach along with a stream-line integration of the free surface. The work has been focused on the FSSI algorithm that has been developed to allow simulating complex and highly dynamic behavior of the seals in the interface between the air cushion, and the water.

The developed fluid-structure interaction solver is based, on one side, on an implicit iteration algorithm, communicating pressure forces and displacements of the seals at memory level and, on the other side, on an innovative wetting and drying scheme able to predict the water action on the seals. The seal dynamics are solved using a Lagrangian FEM formulation for a 3-noded membrane/shell triangle.

Special care has to be taken to ensure the convergence of the fluid-structure interaction algorithm due to the highly flexible and non-linear nature of the seal. We found that relaxation method based on Aitkens acceleration methods worked well under all the case studies analysed.

Different validation and demonstration cases have confirmed the validity of the approach to study the complex dynamic behaviour of the seals.

==7. Acknowledgements==

This work relates to Department of the Navy Grant N62909-10-1-7053 issued by Office of Naval Research Global. The United States Government has a royalty-free license throughout the world in all copyrightable material contained herein.

<span id='_Toc337282459'></span><span id='_Toc337282798'></span><span id='_Toc337283050'></span><span id='_Toc337283220'></span><span id='_Toc337283593'></span><span id='_Toc337479124'></span>

==8. References==

<div id="1"></div>
[[#cite-1|[1]]] E. Oñate, J. García-Espinosa, A finite element method for fluid-structure interaction with surface waves using a finite calculus formulation, Comp. Methods Appl. Mech. and Eng. 2001; 191: 635-660.

<div id="2"></div>
[[#cite-2|[2]]] R. Löhner, C. Yang, E. Oñate, On the simulation of flows with violent free surface motion and moving objects using unstructured meshes, Comp. Methods Appl. Mech. Engng. 2007; 53: 1315-1338.

<div id="3"></div>
[[#cite-3|[3]]] J. García-Espinosa, A. Valls, E. Oñate, ODDLS: A new unstructured mesh finite element method for the analysis of free surface flow problems, Int. J. Numer.  Meth.  Fluids 2008; 76 (9): 1297-1327.

<div id="4"></div>
[[#cite-4|[4]]] D. J. Donnelly, W. L. Neu, Numerical Simulation of Flow About a Surface-Effect Ship, 11th International Conference on Fast Sea Transportation FAST 2011, Honolulu, Hawai. September 2001.

<div id="5"></div>
[[#cite-5|[5]]] S. M. Mousaviraad, S. Bhushan, F. Stern, CFD Prediction of Free-Running SES/ACV Deep and Shallow Water Maneuvering in Calm Water and Waves,  MARSIM 2012. Singapore, April 23-27, 2012.

<div id="6"></div>
[[#cite-6|[6]]] B. Serván-Camas, J. García-Espinosa, Accelerated 3D multi-body seakeeping simulations using unstructured finite elements, J. Comp. Phys. 252 (2013) 382–403.

<div id="7"></div>
[[#cite-7|[7]]] B. S. H. Connell,  W. M. Milewski, B. Goldman, D. C. Kring, Single and Multi-Body Surface Effect Ship Simulation for T-Craft Design Evaluation, 11th International Conference on Fast Sea Transportation FAST 2011, Honolulu, Hawai. September 2001.

<div id="8"></div>
[[#cite-8|[8]]] D. C. Kring, M. K. Parish, W. M. Milewski, B. S. H. Connell, Simulation of Maneuvering in Waves for a High-Speed Surface Effect Ship, 11th International Conference on Fast Sea Transportation FAST 2011, Honolulu, Hawai. September 2001.

<div id="9"></div>
[[#cite-9|[9]]] L. Yun, A. Bliault, Theory &Design of Air Cushion Craft, Elsevier 2005.

<div id="10"></div>
[[#cite-10|[10]]] N. Hirata, O. M. Faltinsen, Computation of Cobblestone effect with unsteady viscous flow under a stern seal bag of a SES, Journal of Fluids and Structures, 2000: 14, 1053–1069.

<div id="11"></div>
[[#cite-11|[11]]] L. J. Doctors, Nonlinear motion of an air-cushion vehicle over waves, Journal of Hydronautics 1975: 9 (2), 44–57.

<div id="12"></div>
[[#cite-12|[12]]] P. A. Sullivan, P. A. Charest, T. Ma, Heave stiffness of an air cushion vehicle bag and finger skirt, Journal of Ship Research 1994: 38 (4), 302–307.

<div id="13"></div>
[[#cite-13|[13]]] Q. Yang, V. Jones, L. McCue, Investigation of Skirt Dynamics of Air Cushion Vehicles under Non-linear Wave Impact Using a SPH-FEM Model, 11th International Conference on Fast Sea Transportation FAST 2011, Honolulu, Hawaii. September 2001.

<div id="14"></div>
[[#cite-14|[14]]] Evaluating Manoeuvering and Seakeeping Performance of a Surface Effect Ship, Technical Report, CIMNE IT-630, October 2012.

<div id="15"></div>
[[#cite-15|[15]]] O. C. Zienkiewicz, R. L. Taylor, J. Z. Zhu, The Finite Element Method: Its Basis and Fundamentals, Butterworth–Heinemann, ISBN 0750663200, 2005.

<div id="16"></div>
[[#cite-16|[16]]] S. R. Idelsohn, E. Oñate, C. Sacco, Finite element solution of free-surface ship-wave problems,  Int. J. Numer. Meth. Engng. 45, 50-528 (1999)

<div id="17"></div>
[[#cite-17|[17]]] R. L. Taylor, E. Oñate, P.-A. Ubach, Finite element analysis of membrane structures, Textile composites and inflatable structures. Springer (2005) 47-68.

<div id="18"></div>
[[#cite-18|[18]]] C. Felippa, B. Haugen, A unified formulation of small-strain corotational finite elements: I. Theory, Computer Methods in Applied Mechanics and Engineering 194 (2005) 2285–2335

<div id="19"></div>
[[#cite-19|[19]]] C. A. Felippa, C. Militello, Membrane triangles with corner drilling freedoms: II. The ANDES element, Finite Elem Anal Des, (1992) 12:189–201.

<div id="20"></div>
[[#cite-20|[20]]] M. A. Crisfield, Non-linear Finite Element Analysis of Solid and Structures, Advanced Topics, vol. 2, Wiley, 1997.

<div id="21"></div>
[[#cite-21|[21]]] F. S. Almeida, A. M. Awruch, Corotational nonlinear dynamic analysis of laminated composite shells, Finite Elements in Analysis and Design 47 (2011) 1131–1145.

<div id="22"></div>
[[#cite-22|[22]]] A. Ibrahimbegovic, S. Mamouri, Energy conserving/decaying implicit time-stepping scheme for nonlinear dynamics of three-dimensional beams undergoing finite rotations, Computer Methods in Applied Mechanics and Engineering 191 (2002) 4241–4258.

<div id="23"></div>
[[#cite-23|[23]]] C. A. Felippa, K. C. Park, C. Farhat, “Partitioned analysis of coupled mechanical systems,” Comput. Methods Appl. Mech. Eng., vol. 190, pp. 3247–3270, 2001.

<div id="24"></div>
[[#cite-24|[24]]] J. G. Valdés, “Nonlinear analysis of orthotropic membrane and shell structures including fluid-structure interaction,” PhD thesis, Barcelona Tech., 2007.

<div id="25"></div>
[[#cite-25|[25]]] B. M. Irons, R. C. Tuck, “A Version of the Aitken Accelerator for Computer Iteration,” Int. J. Numer. Methods Eng., vol. 1, pp. 275–277, 1969.

<div id="26"></div>
[[#cite-26|[26]]] S. F. Zalek , L. J. Doctors, Experimental Study of the Resistance of Surface-Effect-Ship Seals, Proceeding of the 28th Symposium on Naval Hydrodynamics, 2010.

<div id="27"></div>
[[#cite-27|[27]]] R. J. Paredes, Smoothed particle hydrodynamics applied to fluid structure interaction problems involving hydroelastic response, Ph.D. thesis, Stevens Institute of Technology, United States of America, 2013.

<div id="28"></div>
[[#cite-28|[28]]] A. D. Wiggins, S. F. Zalek, M. Perlin, S. L. Ceccio, Development of large scale surface effect ship bow seal testing plataform. 11th International Conference on Fast Sea Transportation FAST 2011, Honolulu, Hawai. September 2001.