Nothing Special   »   [go: up one dir, main page]

Thesis

Download as pdf or txt
Download as pdf or txt
You are on page 1of 137

ECOLE

CENTRALE DE MARSEILLE

Ecole
Doctorale :
Sciences pour lIngenieur : Mecanique, Physique, Micro et Nanoelectronique
Unite de recherche :

Institut de Recherche sur les Phenom`enes Hors Equilibre

`
THESE
DE DOCTORAT
pour obtenir le grade de

DOCTEUR de lECOLE
CENTRALE de MARSEILLE

Discipline : Mecanique des Fluides

Numerical and Experimental Study of the Wave


Response of Floating Support with Partially Filled Tank
Par
Yan SU

Directeur de th`
ese : Bernard MOLIN
Co-directeur de th`
ese : Olivier KIMMOUN

Soutenance prevue le 30 Septembre 2014


devant le jury compose de :
M. Pierre FERRANT
M. Stephane ABADIE
M. Vincent REY
M. Bernard MOLIN
M. Olivier KIMMOUN
M. Julien TOUBOUL

Professeur, Ecole Centrale Nantes, rapporteur


Professeur, Universite de Pau, rapporteur
Professeur, SEATECH
HdR, Ecole Centrale Marseille
Matre de Conference, Ecole Centrale Marseille
Matre de Conference, SEATECH

Contents
1 Introduction
1.1 The practical problem . . . . . . . . . . . . . .
1.2 Previous researches . . . . . . . . . . . . . . . .
1.3 The developments of Boussinesq-type equations
1.4 Outline of the thesis . . . . . . . . . . . . . . .

.
.
.
.

.
.
.
.

.
.
.
.

.
.
.
.

.
.
.
.

.
.
.
.

.
.
.
.

.
.
.
.

.
.
.
.

.
.
.
.

.
.
.
.

.
.
.
.

.
.
.
.

.
.
.
.

.
.
.
.

2 Linear theories for the sloshing in a two-dimensional rectangular tank


2.1 Natural sloshing modes in a two-dimensional rectangular tank . . . . . . .
2.2 Elementary forced motion problems . . . . . . . . . . . . . . . . . . . . . .
2.3 The numerical experiments . . . . . . . . . . . . . . . . . . . . . . . . . . .
2.4 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

.
.
.
.
.
.
.
.

.
.
.
.
.
.
.
.

.
.
.
.

1
1
2
6
7

.
.
.
.

9
10
11
12
14

3 Extended Boussinesq-type models for the sloshing in two-dimensional rectangular tank


3.1 The extended Boussinesq-type models . . . . . . . . . . . . . . . . . . . . . . . .
3.2 The theoretical models of sloshing in a moving flat bottom rectangular tank . .
3.3 Introduction of an energy dissipation term in the dynamic free surface equation
3.4 The numerical models of sloshing in a moving flat bottom rectangular tank . . .
3.5 Integral equation method . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
3.5.1 The theoretical models of the fluid in a moving rectangular tank . . . . .
3.5.2 The discretization of the integral equation . . . . . . . . . . . . . . . . .
3.5.3 The tangential derivative of the velocity potential . . . . . . . . . . . . .
3.6 Validation: numerical experiments for a flat bottom tank with intermediate
water depth . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
3.6.1 Comparisons with linear theory . . . . . . . . . . . . . . . . . . . . . . .
3.6.2 Comparisons with integral equation method . . . . . . . . . . . . . . . .
3.7 Validation: numerical experiments for a flat bottom tank with shallow water
and finite water depth conditions . . . . . . . . . . . . . . . . . . . . . . . . . .
3.7.1 Shallow water case . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
3.7.2 Finite water depth case . . . . . . . . . . . . . . . . . . . . . . . . . . . .
3.8 The numerical models of the sloshing in a moving inclined bottom rectangular
tank . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
3.9 Comparisons between flat bottom tank and inclined bottom tank . . . . . . . .
3.10 Experimental campaigns . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
3.11 Comparisons between numerical solutions and experimental results . . . . . . .
3.11.1 Flat bottom tank . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
3.11.2 Inclined bottom tank . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
i

19
20
21
23
24
26
26
28
29
29
29
30
43
43
47
52
53
56
58
58
60

3.12 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 70
4 The transformation of the equations of motions between the frequency domain and the time domain
4.1 The boundary value problem for a rectangular floating barge in the linear frequency domain theory . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
4.2 Eigen-function expansions and Matching adjacent regions . . . . . . . . . . . . .
4.3 The results of hydrodynamic coefficients and loads . . . . . . . . . . . . . . . . .
4.4 The equations of motions in the time domain . . . . . . . . . . . . . . . . . . . .
4.5 Comparisons between frequency domain and time domain . . . . . . . . . . . . .

71
72
73
76
79
81

5 The comparisons of coupled motions between numerical results and experimental results
85
5.1 The description of the experiments . . . . . . . . . . . . . . . . . . . . . . . . . 86
5.2 The frequency domain computations . . . . . . . . . . . . . . . . . . . . . . . . 88
5.3 The time domain computations . . . . . . . . . . . . . . . . . . . . . . . . . . . 90
5.4 Validations of the coupled motions based on linearized free surface conditions . . 91
5.5 Comparisons between Boussinesq model and linear theory . . . . . . . . . . . . 101
5.6 Comparisons between numerical results and experimental results . . . . . . . . . 110
5.7 The stronger nonlinear case . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 120
5.8 Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 123
6 Concluding remarks and perspectives

ii

125

Chapter 1
Introduction
1.1

The practical problem

The worldwide utilization of natural gas has led producing countries to look for transportation solutions. Marine transportation becomes economically feasible for long distances or when
pipe laying becomes unpractical. In order to ship natural gas in sufficient quantities to make
a complete energy supply project viable, it is liquefied at -163 C thereby reducing its volume
by a factor of 600. The first liquefied natural gas (LNG) carrier Methane Pioneer left the
Calcasieu River on the Louisiana Gulf coast on 25 January 1959. With a capacity of 5,000 m3
distributed over 5 tanks, it sailed to the UK where the cargo was delivered. Nowadays, the
largest LNG carriers in the world, have capacities of 266,000 m3 and larger vessels are under
development.
Generally, LNG carriers sail in either completely filled or empty tanks condition. But for the
gas loading and unloading in offshore locations, the filling ratios in the tanks vary continuously
from 0 up to 100%. Within a range of tank filling levels, the wave induced pitching and rolling
motions of the ship, and the liquid free surface effect, can cause the liquid to move within the
tank. It is possible for considerable liquid motion to take place, creating high impact pressure
on the tank walls. This sloshing phenomenon can be found in industrial areas widely.

Figure 1.1: LNG tanker with Membrane containment systems.

Sloshing is a resonant phenomenon with a strong nonlinear behavior. Including ship applications, sloshing must be considered for almost any moving vehicle or structure containing a
liquid with a free surface, such as harbor resonance, land transportations, onshore tanks and
space applications. Resonant free surface flows in tanks in aircrafts, missiles and rockets have
been the focus of many researches (Abramson [1963], Bauer et al. [1976], etc.). In this case,
sloshing will have a strong influence on the dynamic stability. Large scale sloshing in a lake
with steep sides may be the result of a landslide or earthquake. Shallow water sloshing in a
container can be used to dampen out wind induced motions of tall buildings. A tuned liquid
damper, TLD, is a well known concept for suppressing motions of structures (Fujino et al.
[1992]).
A variety of ship tank shapes exist, comprising rectangular, spherical as well as prismatic
tanks. The liquid may be oil, liquefied gas, or water. Sloshing may appear, for instance, as
a consequence of a collision between two ships, grounding or the collision of a ship with ice.
However, wave-induced ship motions are the primary excitation mechanism. Now the coupled
motions of the ship with tanks are affected by the liquid in the tanks and external fluid which
should be simulated simultaneously.
Because sloshing is a typical resonance phenomenon, it is not necessarily the most extreme
ship motion or external wave loads that cause the most severe sloshing. Therefore, external
wave-induced loads can in many practical cases be described by linear theory. For the violent
liquid motion in the tank, the nonlinear theories must be adopted due to the complicated free
surface elevations. The lateral and angular tank motions cause the largest liquid response
in the coupled motions, but vertical tank motions cannot excite sloshing according to linear
theory and so it is secondary.

1.2

Previous researches

The sloshing phenomena have been studied by numerous researchers through a wide spectrum of different analytical, numerical and experimental approaches. Due to its strong nonlinearity with large fluid motions, wave breaking, spray and mixing of liquid and gas, sloshing is
difficult to predict theoretically.
In general analytical solutions to sloshing problems are based on potential flow theory. Some
analytical solutions can be given for tanks of simple geometrical shapes. A good summary given
by Abramson [1966] presents linear and nonlinear analytical solutions of sloshing in tanks of a
variety of different geometries undergoing harmonic oscillations. Liquid sloshing dynamics with
elastic containers can be found in Bauer et al. [1976]. Recently a fairly comprehensive discussion
of the almost every aspect of liquid sloshing dynamics was given by Ibrahim [2005]. Taking more
attentions in marine hydrodynamics, Faltinsen et al. [2000] use the work of Moiseyev [1958] and
derive a nonlinear analytical theory for sloshing in a rectangular tank with finite water depth.
Moiseyev suggests a general nonlinear method based on potential flow for determination of free
and forced oscillations of the liquid in generally shaped tanks. More systematic discussions of
this analytical theory can be found in Faltinsen and Timokha [2009].
For the numerical methods, the field equations governing the fluid flow may be the complete
Navier-Stokes equations, the Euler equation or the Laplace equation when potential flow is
assumed. Either an Eulerian or Lagrangian description, or a mixture of these, is used. In
the Eulerian description, the grid is fixed relative to the reference frame and the fluid moves
through the grid, as opposed to the Lagrangian description where the coordinate system and
2

the grid points move with the fluid. The type of grid is dependent on the method used to
create an equivalent discrete form of the continuous field equations.
One part of numerical methods for simulating sloshing phenomenon are based on the volume
tracking method. The basic features of volume tracking methods were described by Rider and
Kothe [1998].
The finite difference method, FDM, adopts a structured grid and flow variables are computed for fixed discrete points within this grid using an Eulerian approach. The partial derivatives of the governing equations are reformed by algebraic difference quotients which give a
system of algebraic equations (Mitchell and Griffiths [1980]). Chen and Nokes [2005] use timeindependent finite-difference method to analyze two-dimensional sloshing motions in a rectangular tank based on the primitive 2D Navier-Stokes equations. Both the fully non-linear free
surface condition and fluid viscosity are included. This method has been extended to simulate
the sloshing motions in a three dimensional tank excited by coupled surge and sway motions
in Wu and Chen [2009]. Internal structures baffles are considered in the sloshing motions in
Wu et al. [2012].
The finite volume method, FVM, divides the spatial domain into a finite number of discrete
contiguous control volumes. The solution of the governing equations is given in two steps. First,
an integration over each control volume is performed and secondly, the resulting cell boundary
values are approximated (Versteeg and Malalasekera [2007]). Peric et al. [2009] demonstrated
the application of a procedure to predict internal sloshing loads on partially filled tank walls of
LNG ships subject to realistic wave excitations. A moving grid approach and a finite volume
solution method are designed to allow for arbitrary ship motions.
Another widely used free surface tracking method is the Volume of Fluid, VOF, method. A
thorough description of this method is given by Hirt and Nichols [1981]. The pilot studies are
taken by Van Daalen et al. [2000] which use the VOF method to track the free surface in antiroll tank. Jes
us et al. [2013] compared the multimodal method with an open source package
OpenFOAM and a commercial CFD solver STAR-CCM+ for two phases analysis of sloshing
in a rectangular tank. The discretization of the flow governing equations in the OpenFOAM
package is based on the FVM and the VOF method is used to govern the air and water free
surface elevation in the commercial solver STAR-CCM+.
Finite element methods, FEM, adopt a different discretization process than the finite difference method. A typical work out of the method involves (1) dividing the domain of the
problem into a collection of sub-domains, with each sub-domain represented by a set of element equations to the original problem, followed by (2) systematically recombining all sets of
element equations into a global system of equations for the final calculation. An introduction
of this method is given by Reddy [2005]. Wu et al. [1998] analyzed the sloshing waves in a
three dimensional rectangular tank using a finite element method based on the fully non-linear
wave potential theory. A two dimensional nonlinear random sloshing problem is analyzed by
fully nonlinear wave velocity potential theory based on the finite element method in Wang and
Khoo [2005].
The constrained interpolation profile (CIP) method, which was developed by Yabe and
Aoki [1991] and Yabe et al. [1991] for solving hyperbolic equations, has attracted a great
deal of attention in marine engineering. Hu and Kashiwagi [2004] presented a finite difference
method based CIP method for the numerical simulations of violent free surface flow. The
key points of the CIP method are a compact upwind scheme with subcell resolution for the
advection calculation and a pressure-based algorithm that can treat liquid, gas and solid phases,
irrespective of whether the flow is compressible or incompressible by solving one set of governing
3

equations. The method for hydrodynamic problems using the pressure-based algorithm is called
the CIP combined and unified procedure (CCUP) method. In Kishev et al. [2006], the sloshing
problem is considered as a multiphase problem that includes water and air flows. The CCUP
scheme was adopted for the flow solver and both the CIP scheme and the CIP conservative
semi-Lagrangian with cubic interpolation polynomial scheme were used for interface capturing.
Three dimensional sloshing studies based on the CIP method can be found in Hu et al. [2010].
Boundary element methods, BEM, are based on a potential flow assumption, i.e. the effect of viscosity is neglected and the fluid is assumed incompressible and the flow irrotational.
Applying Greens second identity, the singularities (infinite fluid sources and normal dipoles)
representing the velocity potential are distributed all over the boundary of the fluid domain.
With the solution of the linear system, singularity densities can be given. A careful treatment
of the contact point is needed, and post-breaking flow is difficult to handle. The details of the
programming process can be found in Brebbia and Dominguez [1992]. Numerous improved
versions of this method are given. Gedikli and Erg
uven [2003] presents a variational boundary
element method based on Hamiltons principle that produces symmetric matrices. The boundary integral equations can be regularized by moving the singular points outside the domain. A
modal approach is developed for the second-order analysis of sloshing problem using boundary element method in Firouz-Abadi et al. [2010]. A developed successive boundary element
method is given by Ebrahimian et al. [2013] to determine the symmetric and antisymmetric
sloshing natural frequencies and mode shapes for multi baffled axis-symmetric containers with
arbitrary geometries.
When dealing with violent free-surface flows accompanied by large deformations and fragmentations, the particle methods, including smoothed particle hydrodynamics (SPH) and moving particle semi-implicit methods (MPS), appear to be appropriate candidates due to their
Lagrangian and mesh-free features. The SPH method was firstly proposed by Gingold and
Monaghan [1977] and Lucy [1977] in astrophysics. A good review for this method can be found
in Monaghan [2005]. Rafiee et al. [2009] applied this method to violent sloshing flows in a
2D rectangular tank under harmonic sway excitation motion. The pressure field in forced roll
motions is investigated experimentally and numerically in Delorme et al. [2009]. Comparative
studies on the accuracy and stability of SPH schemes in energetic free surface flows have been
presented in Rafiee et al. [2012]. In this paper, the accuracy, stability and efficiency of improved
techniques based on the SPH method are compared with the experimental results.
Spectral methods were firstly proposed by Cooley and Tukey [1956] in FFT. A review
of applications of spectral methods in general computational fluid dynamics can be found in
Fornberg [1995]. Ferrant and Le Touze [2001] use a pseudo-spectral method based on fully
nonlinear potential theory to study sloshing. The velocity potential is expanded in series of the
natural modes of the tank geometry. The results given agree well with finite element analysis
results by Wu et al. [1998] for a very shallow water case.
Some experimental data can be found in Abramson [1966] and Ibrahim [2005] with rectangular, spherical and prismatic tanks. Kuo et al. [2009] reviewed the relevant Exxon Mobils
sloshing assessment methods. The most significant development is the introduction of a probabilistic based framework that facilitates modeling of the high variability of sloshing impact
pressures. Bunnik and Huijsmans [2009] described a study of model test experiments on a large
scale 2D section (scale 1:10) of an LNG tank in various loading conditions without depressurisation. Using high speed video observations the wave front formed by the bore of the LNG in
resonance was related to measured impacts on the tank hull. Hirdaris et al. [2010] discussed
some of the recent advances in model test tank investigations and reviewed the relevance of
4

allowing unrestricted filling levels in LNG membrane tank ships. The assessment presented in
this work allowed for a reduction in the barred filling range from 80% to 70% of the internal
height of the tank. Xu et al. [2012] carried out sloshing model tests with 1:55 scale membranetype LNG tank oscillating in regular harmonic motions. The purpose of this study has been
to investigate the characteristics of free surface motions, impact pressures and structural responses. During the sloshel project Brosset et al. [2009], the membrane containment systems
have been tested at scale 1 in a long wave tank. Space-focusing waves have been generated in
order to obtain single impact on the membrane (see figure 1.2 right).

Figure 1.2: The experiments are carried out in the Maritime Research Institute Netherlands. These two figures
are the large scale (1:10) tests (left) and full scale tests (right) of membrane containment systems respectively.

The differences in the methods of simulating the sloshing motions in the tanks and the
seakeeping of the ships result in different ways of coupling the computations. One method of
coupling computations in the frequency domain was proposed by Molin et al. [2002]. The forces
and moments produced by the liquid in the tanks are calculated by a linear modal method and
taken into the equations of motions in the frequency domain. The responses of the motions
of the fluid in the tanks and the motions of the floating body can be obtained after solving
the equations of motions. Kim [2002] presented a computational study on sloshing flows, the
anti-rolling tank problem, coupled with ship motions. The three dimensional sloshing flow was
simulated by the finite difference method and the ship motion based on a time-domain panel
method. At each time step, the instantaneous displacements, velocities and accelerations of ship
motions were applied to the excitations of liquid motions, and corresponding sloshing induced
forces and moments were added into the wave-induced excitations. With time stepping, the
time histories of the responses can be given. Another approach for the coupling computations
in the frequency domain was described by Newman [2005]. The interior wetted tank surfaces
are considered as an extension of the exterior ship hull to form one global boundary surface
with two independent fluid domains. This method has been applied to the panel code WAMIT
(Wave Analysis at Massachusetts Institute of Technology). The principal advantage of this
approach is that the exterior panel code can be extended to the internal tanks with relatively
few modifications.
Bunnik and Veldman [2010] proposed two coupling computational methods in both frequency domain and time domain, which are compared to experimental results. In the frequency
domain, it is a linear diffraction method in which the effect of the liquid in the tank is modeled
as a solid inertia of the fluid mass, an added mass and damping of the sloshing fluid. In the
time domain, the sloshing motions in the tank are modeled with the CFD code ComFLOW.
5

The method for coupling computations is similar to the method mentioned in Kim [2002]. Both
methods are applied to sloshing model tests, published by Molin et al. [2008].
Clauss et al. [2010] analyzed a LNG carrier with four membrane tanks numerically with
the radiation-diffraction panel code WAMIT. The numerical method has been described in
Newman [2005]. Tank filling heights as well as the density of the internal fluid can be chosen
arbitrarily. The comparisons of the roll motion RAOs for different filling levels of LNG and
ship drafts are given. And later the frequency shift between the theoretical natural frequency
of the tank and the respective motion peak for a vessel with four tanks mounted to the hull
is investigated in Clauss et al. [2011]. Clauss et al. [2012] presented the relationships between
the frequency peaks and the roll motion, the first transverse sloshing mode of the tank and the
asymmetry of the tank.

1.3

The developments of Boussinesq-type equations

In this thesis, Boussinesq-type models have been used for the sloshing problem, and the
historical developments of the Boussinesq-type model should be given.
The main issue of Boussinesq-type equations is to reduce a three-dimensional problem to a
two-dimensional one by introducing a polynomial approximation of the vertical distribution of
the flow field into the integral conservation laws, while accounting for non-hydrostatic effects
due to the vertical acceleration of water.
The idea was firstly proposed by Boussinesq [1872], and the problem discussed in that paper
was restricted to a horizontal bottom. The vertical velocity was assumed to vary linearly from
zero at the bottom to a maximum at the free surface.
The classical Boussinesq equations are given by Peregrine [1967] in two horizontal dimensions for the case of uneven bottom. However for many practical applications, the range is
limited to kh < 0.75 (where k is wavenumber and h water depth). After that, a number of
enhanced and higher-order Boussinesq equations have been formulated to improve both dispersive and nonlinear properties. As a result the best forms have extended the linear range of
applicability up to kh 6 (Madsen and Schaffer [1998]; Gobbi et al. [2000]).
A breakthrough in treating nonlinearity was made by Agnon et al. [1999]. They presented
a new procedure by which it is possible to achieve the same accuracy in nonlinear properties
as in dispersive properties. But the problem is that an accurate vertical distribution of the
velocity field cannot be provided. Madsen et al. [2002] and Madsen and Agnon [2003] combine
a time-stepping of the exact surface boundary conditions with an approximate series expansion
solution to the Laplace equation in the interior domain. As a result, dispersive and nonlinear
wave characteristics become very accurate up to wavenumbers as high as kh = 40, while the
vertical variation of the velocity field becomes applicable for kh up to 12. All this work has
used the original velocity formulation which is in terms of the surface elevation, a vertical
velocity variable and a two-component horizontal velocity variable.
In order to match with other potential flow solvers and improve computational efficiency,
Jamois et al. [2006] re-casted Boussinesq-type equations in a potential formulation in terms of a
horizontal velocity potential variable which replaced velocity variable. The shoaling properties
of the potential formulation is not as well as that of velocity of formulation. A improved version
was given by Bingham et al. [2009]. In our works, the Boussinesq-type method for the internal
fluid in the tank is based on the theories in Bingham et al. [2009].
6

1.4

Outline of the thesis

The works reported in the thesis focus on the sloshing motions of the fluid in tanks with
forced motions and the coupling motions of a floating body with partially filled tanks. All the
researches are carried out in the framework of potential theories, however the viscous damping
term, is also considered in the numerical models. Generally, the thesis is divided into two parts:
Forced motions and Coupled motions. The tank with forced motions is studied from three
different perspectives: linear theories, Boussinesq-type models and integral equations. The
coupled motions are analyzed in both frequency domain and time domain. Some experiments
are done for comparisons. A coupled analysis of a floating support and partially filled tanks was
given by Molin et al. [2002] in the frequency domain based on linear theories. In this paper,
limitations of linearized free surface conditions appear in predicting the actual free surface
elevations in the tanks. So our discussions will start from this paper.
The chapters 2 and 3 present the computational models of the sloshing motions in a two
dimensional rectangular tank under forced motion.
The linearized theoretical models of the sloshing motions are presented in chapter 2. All
the solutions are in the framework of linearized boundary conditions at the free surface. First,
the computational models for the liquid in a still rectangular tank are given by superposing
natural sloshing modes in the frequency domain. And then the natural sloshing modes are
related to the forced motions of the tank in a pendulum equation. In the derivations, the
velocity potential is separated into two parts which are related by expanding on the same
basis, the natural sloshing modes. The damping term has been considered to take into account
the friction at the walls. Some computational results are presented about the sloshing of liquid
in a rectangular tank with forced harmonic motions.
The extended Boussinesq-type models for the sloshing motions in the tank are given in
chapter 3. Fully nonlinear free surface boundary conditions are adopted for simulating the
internal liquid motions. The theoretical derivations of the Boussinesq-type models are introduced and the numerical discretization processes are given. Both the flat bottom and the
inclined bottom rectangular tanks in forced motions are discussed. The Boussinesq model in
flat bottom tank with different filling levels is validated. The results are compared with linear
theory, integral equation method and results in literature. The nonlinear free surface boundary
conditions are taken for the integral equation method and the mixed Euler-Lagrange method is
used for the iterative process. The comparisons between flat bottom tank and inclined bottom
tank are given and both of them are compared with experimental results.
The following chapters are about the models of the coupled computations of the sloshing
in the tank and the sea-keeping of a floating body.
Chapter 4 presents the linearized boundary value problem of a rectangular barge and the
transformation of the equations of motions from the frequency domain to the time domain. The
hydrodynamic characteristics of the barge are calculated by a semi-analytical method in the
frequency domain. The external fluid domain is divided into a succession of rectangular subdomains where eigen-function expansions are used to express the velocity potential. Matching
of the velocity potentials and horizontal velocities at the successive boundaries gives the solution
to the problem. Based on the hydrodynamic characteristics of the barge in the frequency
domain, the infinite frequency added masses, retardation function and exciting forces in the
equations of motions in the time domain can be obtained. The motions of the barge in the
frequency domain are compared to the solutions in the time domain.
Chapter 5 discusses the coupled motions of a floating body, a rectangular tank fixed on
7

a rectangular barge. The numerical solutions are compared with experimental results. The
experiments were carried out in the wave tank BGO-FIRST, at la Seyne sur mer. A rectangular barge with two partially filled tanks is submitted to irregular beam waves with varying
significant wave heights and peak periods. Both frequency domain and time domain computations are carried out by solving the equations of motions of the floating body. In the frequency
domain, the linear theories for the liquid in the tank and eigen-function expansions method
for the external fluid domain are used. In the time domain computations, the Boussinesq-type
approach is adopted for the sloshing in the tank and the coefficients of the equations of motions are based on the frequency domain computations of the barge. Some comparisons of the
motions of the floating body and the free surface elevations in the tank are given between these
two computational methods and experimental results.
Chapter 6 presents the concluding remarks and perspectives.

Chapter 2
Linear theories for the sloshing in a
two-dimensional rectangular tank
In this chapter, the assumption of linearity is employed in the boundary conditions on the
free surface, and the problem is discussed in the frequency domain within potential theory.
A modal approach is followed to represent the motions of the liquid in the tank which are
combined with the motions of the tank by a pendulum equation. As a result, the transfer
functions of every natural sloshing mode are obtained. The free surface elevations of the liquid
in the tank can be given by superposing the natural sloshing modes.

2.1

Natural sloshing modes in a two-dimensional rectangular tank

The following equations in the frequency domain are employed for governing the liquid in
a tank with linearized free surface condition.
= 0 in the domain
~n = 0 on the wetted surface of tank walls SC
gz 2 = 0 on the mean free surface SF
where

z, t) = (y, z) sin(t + )
(y,

(2.1)
(2.2)

~n is the outward normal vector, g is the gravity acceleration and is the angular frequency.
z

O
O

y
h

l
Figure 2.1: The sketch of 2D rectangular tank

A rectangular tank, of length l and water depth h, is considered (figure 2.1). The origin of
the coordinate system Oyz is located on the middle of the free surface at rest. Following Molin
et al. [2002], the eigen-modes of the liquid in the tank are given by:
n = An0 g cosh n (z + h) cos n (y + l ) sin(n t + n )

n
cosh n h
2

(2.3)

where the natural frequencies n are


n2 = gn tanh n h with n =

n
l

(2.4)

and the free surface elevation is


l
n (y, t) = An0 cos n (y + ) cos(n t + n )
2

(2.5)

For a tank of general shape the eigen-modes can be obtained with a numerical software. It is
easy to show that the potentials n (y, 0) form an orthogonal and complete basis at the free
surface.
As a result the free surface elevation can be expanded as:
X
(y, t) =
An (t)n (y)
(2.6)
n

where
n (y) =

n
n (y, 0)
g
10

(2.7)

And the velocity potential within the tank is then given by


XZ t

(y, z, t) =
( An ( )d )n n (y, z)
n

2.2

(2.8)

Elementary forced motion problems

The tank under forced motion is considered in this section. The elementary problems to be
solved are
j = 0
j ~n = Nj
gjz 2 j = 0

fluid domain
on the wetted surface SC
on the mean free surface SF

(2.9)

where the angular frequency is given and Nj stands for one component of the generalized
normal vector. For the 3 degrees of freedom, Nj writes:
(N2 , N3 , N4 ) = (n2 , n3 , Y n3 Zn2 )

(2.10)

The potential j is divided into two parts:


j = j + j

(2.11)

where j is the infinite frequency radiation potential which satisfies:


j = 0
j ~n = Nj
j = 0

fluid domain
on the wetted surface SC
on the mean free surface SF

(2.12)

and the other term j satisfies


j = 0
fluid domain
j ~n = 0
on the wetted surface SC
2
gjz j = gjz on the mean free surface SF

(2.13)

After expanding j with n , the amplitude An of the sloshing motion n under motion amplitude
Xj is obtained.
2gjn
2
An /Xj =
=

Dnj
(2.14)
n (n2 2 )
n2 2
which can be reformed into a pendulum equation
X
j
An + n2 An =
Dnj X

(2.15)

The coefficients Dnj are defined by


Dnj

R
g SC Nj n dS
R
=
n SF 2n dS
11

(2.16)

which remains valid when Xj (t) is not harmonic since the Dnj do not involve the frequency .
Introducing a linear damping term B1n , the equation of motion can be expressed:
X
j
An + B1n A n + n2 An =
Dnj X
(2.17)
j

Following Molin et al. (2002), the damping coefficient due to laminar friction at the walls is
given by
r
1
n l 2n h
B1n = 2n
(b + l + b
)
(2.18)
2n lb
sinh 2n h
where b is the width in the x-direction and is the molecular viscosity.
The total velocity potential in the tank is given by
XZ t
X

(y, z, t) =
( An ( )d )n n (y, z) +
X j (t)j
n

(2.19)

~
And the generalized force tensor F~ = (F~ , C):
X
~
F~ =
An (t)f~n Ma ()X

(2.20)

where Ma () is the infinite frequency added mass matrix and


Z
fnj = n
n Nj dS

(2.21)

SC

It may be noticed that the fnj and Dnj coefficients are obtained through the same integral over
SC and that they are simply related:
Dnj =

g
f
R R nj
2
n
2n dS
SF

(2.22)

Here we give the expressions of Ma () = [mij ], along y-direction, around the x-direction and
their coupling terms.
1
8l2 X tanh n h
m22 = lh 3

n odd
n3

1
lh2 l3 8l3 X 1
2
1
m42 = m24 =
+ 4
(
+ n h tanh n h)
4

2
6
n odd n cosh n h
1
lh3 l3 h 8l4 X 1
4n h
m44 =
+
+ 5
(4
tanh

2n h2 tanh n h)
n
5

3
12
n odd n
cosh(n h)

2.3

(2.23)

The numerical experiments

As mentioned above, the motion of the tank Xj (t) is not necessary harmonic, which means
the time domain computations should be applied. Here, some results of a 2D partially filled
rectangular tank with forced harmonic motions are presented. Based on the motion equations,
the motions of the liquid in the tank can be solved in time domain.
12

The temporal derivatives of An and A n can be introduced as follows:


d
An = A n
dt
d
2 + Dn4 X
4 w 2 An B1n A n
An = An = Dn2 X
n
dt

(2.24)

Initially, the modal amplitudes An and A n are set to zeros, which means the right-hand side
of the two equations are known. The coefficients in these two equations can be derived analytically. With the temporal derivatives on the left-hand side, An and A n can be updated by an
iterative process where the classical 4th Runge-Kutta process is adopted. With the solution
of the amplitude of each sloshing mode, the free surface elevation and the forces and moments
can be calculated by the superposition of them.
For the frequency domain, the equation of motion should be transformed as follows:
2 + Dn4 X
4
An + B1n A n + n2 An = Dn2 X

(2.25)

Here only the sway and roll motions are considered. Due to the harmonic motions, the amplitudes should be satisfy the following equation:
[( 2 n2 ) + iB1n ]An = 2 Dn2 X2 + 2Dn4 X4

(2.26)

where is the angular frequency of the harmonic motion. Here An , X2 and X4 are the
amplitudes of each mode and the motions. The amplitudes of the free surface elevations can
be obtained by equation (2.6).
A 2D rectangular tank with length 1.08 m and filling level 0.15 m is considered. Due to
the computations of damping coefficient B1n , the width of the tank in the x-direction is taken
equal to 0.1 m. The motions of the tank are harmonic sway motions with amplitude 1 mm.
This is the basic case. The range of periods is from 0.6 rad/s to 2.0 rad/s. Based on the former
dispersion relation, the natural angular frequencies of the first and third sloshing modes are
3.423 rad/s and 8.6 rad/s which means the periods are 1.83 s and 0.73 s respectively.
In figures 2.2 and 2.3, the modal amplitudes An are calculated by the equation (2.26) using
the basic case. The range of excitation periods is from 0.6 second to 2.0 second. The first modal
amplitude A1 and the third modal amplitude A3 divided by the amplitude of the sway motion
1 mm give the values of the ordinates. The different values of the viscosity are compared in
these two figures. Obviously, with the increase of the molecular viscosity , the peak values of
the modal amplitudes around the first natural period and the third natural period decrease.
The time varying modal amplitudes An (t) can be solved by the equation (2.24). The free
surface elevations are calculated by equation (2.6). In our numerical computations, the number
of modes is taken equal to nine. Figure 2.4 and 2.5 present the time histories of the free surface
elevations at the position of the left wall during 200 seconds. The excitation periods are 0.7
second (figure 2.4) and 1.6 second (figure 2.5) respectively. As is well known, the beating wave
can be obtained by a linear superposition of two sinusoidal waves with different frequencies.
In figure 2.4, the beating period can be calculated by 2/(| 1 |) where is the excitation
frequency and 1 is the first natural frequency. The case of excitation period 1.6 second is
shown in figure 2.5.
The comparisons of the first mode and the third mode between the frequency domain
and the time domain are given in figures 2.6 and 2.7. The first modal amplitudes A1 in the
frequency domain are solved by the equation 2.26 with the excitation periods from 0.6 second to
13

2.0 second. The viscosity is taken equal to 1106 m2 s1 . When the excitation period is close
to the first natural period, the resonance phenomenon becomes significant. In the time domain
computations, the values of A1 (t) are calculated by the equation 2.24 and the amplitudes of
A1 (t) after the transients have died out are adopted. In figure 2.6, the first modal amplitudes
divided by the amplitude of sway motion are compared between the frequency domain and the
time domain computations. The comparison of the third sloshing modes A3 in the frequency
domain and the time domain is given in figure 2.7. The results in both of the figures coincide
well.

2.4

Conclusion

The modal approach is an efficient method for simulating the sloshing motions. The motion
is simply considered as a linearized superposition of the natural sloshing modes. Due to the
linearized assumption and the symmetric shape of the tank, the resonance in the even natural
sloshing modes cannot be predicted by this modal approach. The derivation process is carried
out in the frequency domain, but the excitation motions are not limited in harmonic motions.
From the time domain computations, the beating phenomenon can be found in the beginning
and the transients gradually disappear due to the damping term in the equations of motions.

14

Frequency domain analysis


50

Amplitude of mode / Amplitude of motion

45

1.E06
5.E06
1.E05

40
35
30
25
20
15
10
5
0
0.6

0.8

1.2
1.4
Period (s)

1.6

1.8

Figure 2.2: The amplitudes of the first sloshing mode divided by the amplitude of sway motion give the values
of ordinates. Three different values of the viscosity are adopted in the comparisons. [basic case]

Frequency domain analysis


70
1.E06
5.E06
1.E05

Amplitude of mode / Amplitude of motion

60

50

40

30

20

10

0
0.6

0.8

1.2
1.4
Period (s)

1.6

1.8

Figure 2.3: The amplitudes of the third sloshing mode divided by the amplitude of sway motion give the values
of ordinates. Three different values of the viscosity are adopted in the comparisons. [basic case]

15

The point on the wall

x 10

Free surface elevations (m)

2
0

4
6

20

40

60

80

100
120
Time (s)

140

160

180

200

Figure 2.4: The time history of the free surface elevations are given for the point on the left wall. The period
of the sway motion of the tank is 0.7 second. [basic case]
3

The point on the wall

x 10

Free surface elevations (m)

1
0

1
2

20

40

60

80

100
120
Time (s)

140

160

180

200

Figure 2.5: The time history of the free surface elevations are given for the point on the left wall. The period
of the sway motion of the tank is 1.6 second. [basic case]

16

Mode 1
50

Amplitude of mode / Amplitude of motion

45
40
35
30
25
20
15
10
5
0
0.6

0.8

1.2
1.4
Period (s)

1.6

1.8

Figure 2.6: The amplitude of sloshing mode divided by the amplitude of motion gives the value of the ordinate.
The results (red) in the time domain are compared to the solutions (blue) of the frequency domain computations.
In the time domain computations, the amplitudes of A1 (t) after the transients have died out are divided by
the amplitude of motion. [basic case]
Mode 3
70

Amplitude of mode / Amplitude of motion

60

50

40

30

20

10

0
0.5

1.5

Period (s)

Figure 2.7: The results in the time domain are compared to the solutions of the frequency domain. The third
sloshing mode is adopted here. [basic case]

17

18

Chapter 3
Extended Boussinesq-type models for
the sloshing in two-dimensional
rectangular tank
In this chapter, sloshing in a flat bottom rectangular tank and in an inclined bottom tank
will be considered. Both of them are simulated based on the Boussinesq-type approach, which
means non-linear free surface boundary conditions are adopted. Different filling levels in the flat
bottom tank are discussed. The results for the flat bottom rectangular tank will be compared
with the former linear theory, integral equation method and results from literature. Due to
asymmetry of the inclined bottom tank, a significant difference can be found in the second
sloshing mode between flat bottom tank and inclined bottom tank. The numerical results
calculated by the Boussinesq model are then compared with the experimental results.

19

3.1

The extended Boussinesq-type models

Firstly, the governing equations for a general three-dimensional nonlinear wave problem
are given. A Cartesian coordinate system is adopted here, with the origin located on the
still water plane and the z-axis pointing vertically upwards. The boundaries of fluid domain
are given by the bottom at z = h(x) and the free surface at z = (x, t) with x = [x, y].
Following Zakharov [1968], the free surface boundary conditions are written in terms of the
= (x, , t) and the vertical velocity w = (z )z= defined directly on the
velocity potential
free surface
w(1
t +
+ ) = 0
t + g + 1 ()
2 1 w 2 (1 + ) = 0

2
2

(3.1)
(3.2)

where = [/x, /y] is the horizontal gradient operator and g is the gravitational acceleration. The Laplace equation in the fluid domain and the kinematic boundary condition on the
bottom are given as follows:
2 + zz = 0
(3.3)
w + h = 0, z = h(x)

(3.4)

+ (z z(x))w + 1 (z z(x))2
(2) + 1 (z z(x))3
(3) + . . .
(x, z, t) =
2
6

(3.5)

A Boussinesq-type Taylor series expansion of the velocity potential about an arbitrary level
z = z(x) will be used. The Taylor series expansion of the solution (x, z, t) about an arbitrary
vertical position z = z(x) is given by

where
(x, z, t)
(1)

w = =

z
z=
z (x)

=
(0) = (x, z(x), t),

n (x, z, t)
(n)

=
,

z=
z (x)
z n

n = 2, 3, . . . , .

(3.6)

Introducing = z z(x) for brevity, we substitute the potential into the Laplace equation
(n) is
and set the coefficient of each power of equal to zero. Finally, a recursion relation for
obtained:
(n) =

1
(n2) + 2
(n1) + 2 z
(n1) ), n = 2, 3, . . . , .
(2
z
1 +
z2

(3.7)

Following Bingham et al. [2009], the recursion relation can be simplified based on assuming
small bottom slope that is z(x) where 1. Collecting terms at each order of and keeping
terms up to , we get:
(n) = 2
(n2) + 2
(n1) , n = 2, 3, . . . , .

(3.8)

(n) can be expressed in terms of


and w.
Successively solving for n = 2, 3, . . . , , each
In
order to improve the accuracy of the truncated approximation, an enhancement operator with
, w is adopted.
new expansion variables
= Lp (
, w = Lw (

z )
z )w .
20

(3.9)

where the operators Lp = L0 +


z L1 , Lw = L0 +
z L2 are given:
L0 (
z ) =
L1 (
z ) =

N
+1
X
n=1

2N
X
n=0

a2n1 z2n1 2n2 ,

2n z2n 2n
L2 (
z ) =

(3.10)
N
+1
X
n=1

b2n1 z2n1 2n2

(3.11)

Here, the coefficients ai , bi can be obtained by optimizing the linearized shoaling behavior.
Applying the enhancement operators, the velocity potential is given finally as follows
+ (J02 +
(x, z, t) = (J01 +
z J11p )
z J12p )w

(3.12)

The same procedure can be applied to the horizontal and vertical velocities
+ (J02 +
u(x, z, t) = (J01 +
z J11u )
z J12u )w
+ (J01 +
w(x, z, t) = (J02 2
z J12w )
z J11w )w
.

(3.13)
(3.14)

Setting N = 1, the operators are given as follows:


2 z2 2
3 z2 2
J01 = 1 + ( + ) J02 = + ( +
)
2
10
6
10
1 a1
3
J11p = a1 z + [ 2 z( + ) + z3 a3 ]2
3
5
2
4
2 2

z
1
b1
J12p = ( 2 + b1
z ) + [ +
3 z( + ) +
z 3 b3 ]2
6
10
15
6
2
2
2
J12w = [ +
z ( + a1 )]
5
z2
1 b1
2 3
+
2 z( + ) + z3 b3 ]2
J11w = (2 + zb1 ) + [
3
5
5
2
3
2
1

z
3
a1
1
J11u = [ + z( + a1 )]2 + [ +
2 z( + ) + z3 ( + a3 )]4
5
2
10
10
2
30
2
2
3
z
1
J12u = 1 + [

+
z ( + b1 )]2
(3.15)
2
10
5
where a1 = 0.465737, a3 = 0.0653131, b1 = 0.138568 and b3 = 0.0138568. These coefficients are optimized for a range of kh from 0 to 6.

3.2

The theoretical models of sloshing in a moving flat


bottom rectangular tank

A two-dimensional rectangular tank partially filled is considered here with length l and
water depth h, and a tank-fixed coordinate system Oyz lies on the free surface at rest in the
middle of the tank (figure 2.1). Small-amplitude motions in sway, heave and roll of the tank are

considered, with the axis of rotation in the middle of the bottom O . The velocity ~v = (0, vy , vz )
on the boundary of the tank is given:
~ = (0, 0, h) + (0, y, z) = (0, y, h + z)
~v0 = (0, v02 , v03 ), w
~ = (w1 , 0, 0), ~r = O~ O + OM
~v = ~v0 + w
~ ~r = (0, v02 w1 (z + h), w1 y + v03 )
21

(3.16)
(3.17)

where v02 is the translational velocity, v03 is the vertical velocity and w1 is the angular velocity.
Firstly, the governing equations in the tank-fixed coordinate system should be reformed.
The temporal derivatives in the governing equations defined in the inertial coordinate system
and in the non-inertial coordinate system are related by the following expressions:

~v ;
t

~v
t

(3.18)

The left-sided physical quantities are based on the inertial coordinate system and the rightsided are based on the non-inertial coordinate system. And ~v is the velocity of non-inertial
system relative to the inertial system.
The governing equations of the liquid in the tank are given by:
l
y = ; z = vz z = h;
2
t vy y z + y y = 0 z = ;
1
t vy y vz z + g + (2y + 2z ) = 0 z =
2
yy + zz = 0;

y = vy

(3.19)
(3.20)
(3.21)

Following Zakharov [1968], the free surface conditions are written in terms of the velocity
= (y, , t) and the vertical velocity w = (z )z= defined on the free surface. Based
potential
on the chain rule, the partial derivatives follow:
y = (y )z= + (z )z= y = v + w
()
y
t = (t )z= + (z )z= t = (t )z= + w(
()
w
+ (vy v)y )

(3.22)
(3.23)

Substitution in the boundary conditions on the free surface gives the evolution equations for
and :

y y + (1 + 2 )w;
t = vy y
y
2 + 1 w 2 (1 + 2 )
t = vy
y + vz w g 1

y
2 y 2

(3.24)
(3.25)

and the elevation can be advanced in time


In the numerical computations, the potential
t and t if the vertical velocity w is known.
from the derivatives
The potential is separated into two parts: = + , with a particular solution
satisfying the Laplace equation and the no-flow condition on the walls given by:
= v02 y +v03 z +w1

l
2 sinh(n z)
n cos(n (y + ))[
(z +h)];
2

n cosh(n h)
odd

n =

n
l

n =

4l
n2 2

(3.26)
Here includes two parts: one is directly proportional to the tank velocity, and the other is
due to the roll motion. The potential is the solution of the following equations.
yy + zz = 0;

l
y = 0 y = ; z = 0 z = h (3.27)
2
t vy y z z + y (y + y ) = 0

1
t + t vy y vy y vz z vz z + g + [(y + y )2 + (z + z )2 ] = 0;
2

z=
(3.28)

22

where
l
2 sinh(n z)
n n sin(n (y + ))[
(z + h)]
2 n cosh(n h)
n odd
X
l 2 cosh(n z)
1]
z = v03 + w1
n cos(n (y + ))[
2
cosh(n h)

y = v02 w1

(3.29)
(3.30)

n odd

Then the vertical velocity on the free surface is given by:


w
= (z )z= = (z )z= + (z )z=

(3.31)

The hydrodynamic loads are given by the integrations of the pressure on the wetted surface as
follows:
Z
Z
~
~
F = p ~n ds; M = p ~r ~n ds;
(3.32)
s

1
p = (t vy y vz z + (2y + 2z ) + gz)
2

3.3

(3.33)

Introduction of an energy dissipation term in the


dynamic free surface equation

According to the definitions in linear theory, a linear damping term of the first sloshing
mode can be obtained by comparing the amplitude. Firstly, the equation of motion (equation
2.17) without external excitation for the first sloshing mode can be written as follows:
A1 + B11 A 1 + 12 A1 = 0

(3.34)

with B11 obtained as equation (2.18) to account for friction on the tank wall.
The solution can be given as follows,
A1 = A0 eB11 t/(21 ) cos(1 t + )

(3.35)

The linearized dynamic and kinematic free surface conditions are taken as:
t + g = 0
t z = 0

(3.36)
(3.37)

A dissipative term is introduced in the first equation.


t + g + = 0
t z = 0

(3.38)
(3.39)

tt + t + gz = 0

(3.40)

Eliminating the variable , we obtain:

The expression of velocity potential is taken as:


= {

l cosh 1 (z + h) it
igA(t)
cos(1 y + )
e
}

2
cosh 1 h
23

(3.41)

and the temporal derivatives of velocity potential are given:

iAg
) cos(1 y +

2
i Ag ) cos(1 y +
tt = {(iAg 2Ag

t = {(Ag

l cosh 1 (z + h) it
)
e
}
2
cosh 1 h
l cosh 1 (z + h) it
)
e
}
2
cosh 1 h

(3.42)
(3.43)

The coefficient of term tt + t + gz should be equal to zero. After eliminating the term
O(2 ), the amplitude A satisfies,
Ag = 0
2Ag
(3.44)

hence

A(t) = A0 et/2

(3.45)

Compared with equation (3.35), the coefficient can be obtained:


=

B11
1

(3.46)

The linear damping term of the Boussinesq model is added in right-hand side of equation
(3.25):
B11

(3.47)
1

3.4

The numerical models of sloshing in a moving flat


bottom rectangular tank

Based on the above sections and setting N = 1, the velocity potential and the vertical
velocity for constant water depth in the tank are given:
(y, z, t) = J01 + J02 w
w(y, z, t) =

(3.48)

= J02 2 + J01 w
z

(3.49)

where

2 z2 2
3 z2 2
+ ) ; J02 = + ( +
) ;
2
10
6
10
The potential on the free surface is given by:
J01 = 1 + (

F S
FS
= J01
+ J02
w

= z

= z z

(3.50)

(3.51)

And the vertical velocity on the bottom


b
b
w = 0
w = J02
2 + J01

(3.52)

The notation F S expresses that the operators have been evaluated with z = and the notation
a linear
b expresses that the operators have been evaluated with z = h. With the knowns ,
system of unknowns w is obtained.
   
 FS
FS
J01
J02

(3.53)
=
b
b
J02
2 J01
0
w
24

The domain is discretized and a finite difference scheme is adopted for the spatial derivatives
in the numerical computations (f : arbitrary function).
fi
1
=
(fi2 8fi1 + 8fi+1 fi+2 )
y
12(dy)
2 fi
1
=
(fi2 + 16fi1 30fi + 16fi+1 fi+2 )
2
y
12(dy)2
1
3 fi
=
(fi2 + 2fi1 2fi+1 + fi+2 )
3
y
2(dy)3
4 fi
1
=
(fi2 4fi1 + 6fi 4fi+1 + fi+2 )
4
y
(dy)4

(3.54)

Here the assumption y = 0, y = l/2 is used on the walls.


1
1
=
(1 80 + 82 3 ) = 0
y
12dy
n
1
=
(n2 8n1 + 8n+1 n+2 ) = 0
y
12dy
1 = 3 , 0 = 2 ;

n2 = n+2 , n1 = n+1

(3.55)
(3.56)

And the same process for the potential on the walls


1
1
=
(1 80 + 82 3 ) = 0
y
12dy
1
n
=
(n2 8n1 + 8n+1 n+2 ) = 0
y
12dy
1 = 3 , 0 = 2 ; n2 = n+2 , n1 = n+1
With a flat bottom, the unknowns w have the following relationship

(3.57)
(3.58)

1 = 3 ; 0 = 2 ; n2 = n+2 ; n1 = n+1

(3.59)

w1
= w3 ; w
0 = w2; wn2
= wn+2
; wn1
= wn+1

(3.60)

In order to obtain a band matrix, the linear system should be rebuilt with the order of unknowns

like [1 , w1, 2 , w2 , . . . , n1, wn1


, n , wn ]. For n = 5, the band width is equal to 11. The open
source codes LAPACK is used for solving the linear system.
With the solutions of the linear system, the derivative z on the free surface can be obtained.
The vertical velocity on the free surface can be separated into two parts:
w = (z )z= + (z )z=

(3.61)

The second part is the derivative of the particular solution on the free surface.
Based on the dynamic and kinematic free surface boundary conditions, the temporal deriva t can be given, which will be used for updating in time the values of and .

tives t ,
With the free surface elevations and the solutions of linear system, the pressure distributed
on the wetting surface can be calculated. For the spatial derivative y , the first part due to
the potential is given by:

= J01 + J02 w
(3.62)
y
25

where , w are the solutions of the above linear system. And the derivative y can be
calculated by adding the second part y .
The temporal derivative in the computations of the pressure on the wetted surface is calculated by the backward difference:
t =

(t) (t t)
t

(3.63)

The classical forth Runge-Kutta method is adopted for the time-stepping processes and the
Savitzky-Golay filter is used for smoothing the data in each time step.
For the beginning of the process, the ramp function is taken as follows:

 
t
1
1 cos
t 6 nT
ramp(t) =
2
nT

1
t > nT

where T is characteristic period of the motion and n is number of period where the ramp takes
place.

3.5
3.5.1

Integral equation method


The theoretical models of the fluid in a moving rectangular
tank
z

~n
y

~s

l
Figure 3.1: Sketch of sloshing in a rectangular tank.

A two dimensional rectangular tank (figure 3.1) with partially filled liquid is considered
here with the length l and filling level h. The coordinate system is fixed on the still free surface
in the middle of the tank and the normal vector at the boundary of the fluid points to the
outside.
For the liquid in the tank, the boundary of the fluid domain is divided into two parts: the
free surface and the wetted surface. At an instant, we assume that the velocity potential
on the free surface is known and the governing equations of the fluid in the tank are given as
follows:
2 = 0 in fluid domain
n = f1 on the wetted surface
= f2 on the free surface
26

(3.64)
(3.65)
(3.66)

where f1 is the normal velocity and f2 is the velocity potential.


According to the Greens theorem, the domain problem can be transformed into the boundary value problem based on the following boundary integral equation:
Z Z
Z
2
2
( G G )d = (Gn Gn )dS = 0
(3.67)

where is the fluid domain and S is the boundary of the fluid domain. G is the Greens function
which is selected with the basic solution for the Laplace equation. In the two dimensional
problem, the function takes the form as follows:
G = lnr(P, Q)

(3.68)

Here r(P, Q) represents the distance between the source point Q distributed along the contour
S and the field point P belonging to the contour S. Gn and n are the normal derivatives of
the Greens function and of the velocity potential.
The Greens function is singular when the source point Q coincides with the field point P .
In this case, a small circle of radius is used to bypass the point Q and the integration will
be performed over the new boundary. Taking the limit as 0, the equation (3.67) can be
reformed as follows:
Z
(P ) = (Gn Gn )dS
(3.69)
S

where is the interior angle at the field point P lying on the boundary.
If the linear segments are used for approximating the boundary of the fluid domain, the
integral equation can be discretized as follows:
XZ
XZ
i i
Gn dS
(3.70)
Gn dS =
+
Sj

Sj

The situations in the two parts of the boundary are different. On the free surface, the velocity
potential is known but the normal derivative of the velocity potential is unknown. Otherwise,
on the wetted surface, the normal derivative of the velocity potential is known but the velocity
potential is unknown. Fortunately, the quantities of unknowns in the integral equation equal
to the quantities of the equations in the system.
Finally, a linear system with unknowns on the wetted surface and n on the free surface
can be built.
N
N
X
X
Aij j =
Bij nj
(3.71)
j=1

j=1

where Aij and Bij are the coefficients.


Rebuilding the linear system with unknowns on the left-hand side, the normal derivatives
on the free surface can be given with the solutions of the linear system. Based on the normal
derivatives of the velocity potential n and tangential derivatives of the velocity potential
s , the spatial derivatives of the velocity potential y and z are calculated by the following
equations:
y = s nz + n ny
z = s ny + n nz
27

(3.72)
(3.73)

where (ny ,nz ) are the components of the normal vector ~n.
The kinematic and dynamic free surface conditions for the tank-fixed coordinate system
are the same as mentioned in the section 3.2
y vy )y (1 + 2 )w = 0
t + (
y
1
1
t vy
y vz w + g +
2 w 2 (1 + y2 ) = 0

2 y 2

(3.74)

where vy and vz are the velocity components. The vertical velocity w = (z )z= and the
= z= are defined on the free surface. With the initial conditions on the
velocity potential
free surface and the temporal derivative of the free surface elevation and the velocity potential,
the time histories of the free surface elevation and the velocity potential can be calculated with
the Mixed Eulerian Lagrangian method.

3.5.2

The discretization of the integral equation

The boundary of the liquid will be approximated by a series of linear elements which means
the unknowns vary linearly over the elements. In the global coordinate system (y-z), P is a
field point and Q is a source point. Here, the local coordinate system (, ) is adopted. ~n and
~s denote the normal and tangential (clockwise) unit vectors of this element respectively. The
origin is fixed on the field point P and -axis is in the direction of ~s.
= P~Q ~s

= P~Q ~n

(3.75)

The boundary integrals can be approximated as the sum of individual integrations on each
element in a local coordinate system,

Z
X Z j+1 (lnr)
X  j+1 I2 I4
(lnr)
I4 j I2

dS =
dS =
j +
j+1

n
n

j+1
j
j+1
j
S
j
j
j
Z

(lnr)dS =
n

XZ

j+1

(lnr)dS =
n
j

j+1I1 I3
j+1 j

I3 j I1
j+1 j

(3.76)
!

j+1

(3.77)

where
I1
I2
I3
I4

 j+1

1

2
2
1
=
lnrd =
,
ln( + ) 2 + 2 tan
2
j
j
 j+1

Z j+1

1
,
(ln r)d = tan
=
n
j
j
Z j+1
1
=
ln rd = [( 2 + 2 )(ln( 2 + 2 ) 1)]j+1
j ,
4
j
Z j+1

1
=
(ln r)d = [ln( 2 + 2 )]j+1
j .
n
2
j
Z

j+1

(3.78)
(3.79)
(3.80)
(3.81)

When the boundary of the fluid domain is discretized into linear elements, the points at
corners require special attentions because the conditions on both sides may be not the same.
28

Since the potential is unique at any point of the boundary, no problem will be introduced into
the point at corner. But the normal derivatives of the points at corner are not the same value
from the left and right limits. this problem can be arranged by treating the connect point as
two points (Brebbia and Dominguez [1992]).

3.5.3

The tangential derivative of the velocity potential

High order numerical differentiations are usually not as good as expected, and most of the
time, they are not necessary. In the lower curvature region, three-point differentiation is usually
enough, while in the region where high curvature occurs, higher-order differentiation always
loses its accuracy and is not accurate as lower-order differentiation. Here, the three-point
Lagrangian interpolation is employed.
 

= Lj1 j1 + Lj j + Lj+1 j+1


(3.82)
s j
where
Lj1 =

sj
(sj1 + sj )sj1

Lj =
Lj+1 =

sj sj1
sj1 sj

(3.83)
(3.84)

sj1
(sj1 + sj )sj

(3.85)

and
sj1 =

(xj xj1

)2

+ (yj yj1

)2 ,

sj =

(xj+1 xj )2 + (yj+1 yj )2

(3.86)

The connections of the free surface and the walls should be considered with care. The tangential
derivatives are calculated by the forward difference and backward difference.
(s22 + 2s1 s2 )1 + (s1 + s2 )2 2 s21 3
s1 s2 (s1 + s2 )
2
2
sn1 n2 (sn2 + sn1 ) n1 + (s2n2 + 2sn2 sn1 )n
=
sn2 sn1 (sn2 + sn1 )
s,1 =

s,n

3.6
3.6.1

(3.87)
(3.88)

Validation: numerical experiments for a flat bottom


tank with intermediate water depth
Comparisons with linear theory

A rectangular tank is adopted here with length l=1.08 m and filling level h=0.15 m. Forced
harmonic sway and roll motions are used. The results calculated by the Boussinesq model are
compared with the linear theory (nine modes). The viscosity is not considered here.
Figure 3.2 presents the time histories of the free surface elevations on the left wall, one
third and two thirds of length away from left wall. The horizontal harmonic motion is used
with excitation period 0.7 s and excitation amplitude 1 mm. The excitation period is close to
29

the third natural sloshing mode 3 which includes one and a half wavelength in the tank. The
beating period 2/| 3 | is around 17 seconds. The results shown in the figure coincide well.
The excitation period 1.6 second of forced motion with amplitude 1 mm is shown in figure 3.3,
which is close to the first natural sloshing mode.
From the figures, the Boussinesq model gives an overall good agreement. It should be noted
that acceleration of forced motion is used for linearized modal method, but Boussinesq model
adopts the velocity for external excitations. The phase lag of half of the excitation period
should be treated in the comparisons.
The larger excitation amplitude 2 mm is used for the next comparisons. Figure 3.4 presents
the comparison between Boussinesq model and linear theory with forced period T =0.7 second.
The discrepancies between the two models can be found on the positions of troughs. The
nonlinear free surface conditions used in Boussinesq model give more flat trough than linear
theories. The excitation period 1.6 second is shown in figure 3.5. The first sampling point
is on the left wall, and more flat trough can also be found. Due to steeper crest and flatter
trough, the crests of free surface elevation on the second and third sampling points calculated by
Boussinesq model are lower than that with linear theories. The free surface heights (elevation +
filling level) in the tank at instant 7.15 s and 7.95 s are shown in figures 3.6 and 3.7 respectively.
The length of the tank has been normalized. From the figures, the nonlinear characteristic of
the free surface can be found clearly.
The forced harmonic roll motions are used and the axis of rotation is located on the origin
of coordinate system. Figure 3.8 presents the free surface elevations at 0 cm, 36 cm and 72
cm away from the left wall. The excitation period is 0.7 second and motion amplitude is 0.5
degree. The results calculated by the Boussinesq model give an overall good agreement with
linear theory. The flatter trough can be found in the Boussinesq model. The higher excitation
amplitude 1.0 degree is compared in figure 3.9 which the nonlinearity is more significant. The
excitation period 1.6 second which is close to the first sloshing mode is shown in figure 3.10.
The motion amplitude is 0.5 degree. The free surface elevations on the left wall calculated
by Boussinesq model present steeper crests than linear theory. At one third and two thirds
of length away from left wall, the amplitude of the crest is lower than linear theory due to
nonlinear free surface profile.

3.6.2

Comparisons with integral equation method

The rectangular flat bottom tank with filling level h=0.15 m and length l=1.08 m is also
considered. The initial free surface profile is set to (y) = A cos(2y/l) with the amplitude
A=0.004 m. The integral equation method is firstly validated by fixed tank. Figure 3.11
presents the time histories of the free surface elevations at the wall during 50 seconds. The
energy attenuation is not present and the signal is keeping sinusoidal waves because no damping
terms exists.
The tank with forced harmonic horizontal motions is calculated with the free surface at rest
initially. According to the dispersion relation, the natural periods of the first and third sloshing
modes are close to 1.8356 second and 0.7306 second. In figure 3.12 and 3.13, the excitation
periods 0.73 and 1.83 second are taken with excitation amplitude a=1 mm. The viscosity
is not considered here. The Boussinesq model and integral equation method with nonlinear
free surface conditions are compared with linear modal approach (nine modes). From the
comparisons, the nonlinear free surface elevations are steeper at crests and flatter at troughs.

30

0.015
linear
Boussinesq

free surface elevation (m)

0.01

0.005

0.005

0.01

0.015

10

15

20

25
time (s)

30

35

40

45

50

0.015
linear
Boussinesq

free surface elevation (m)

0.01

0.005

0.005

0.01

0.015

10

15

20

25
time (s)

30

35

40

45

50

0.015
linear
Boussinesq

free surface elevation (m)

0.01

0.005

0.005

0.01

0.015

10

15

20

25
time (s)

30

35

40

45

50

Figure 3.2: The time histories of the free surface elevations at 0 cm, 36 cm and 72 cm away from the left wall.
The excitation period is 0.7 s and the excitation amplitude is 1 mm.

31

x 10

linear
Boussinesq

free surface elevation (m)

10

15

20

25
time (s)

30

35

40

45

50

x 10

linear
Boussinesq

free surface elevation (m)

10

15

20

25
time (s)

30

35

40

45

50

x 10

linear
Boussinesq

free surface elevation (m)

10

15

20

25
time (s)

30

35

40

45

50

Figure 3.3: The time histories of the free surface elevations at 0 cm, 36 cm and 72 cm away from the left wall.
The excitation period is 1.6 s and the excitation amplitude is 1 mm.

32

0.03
linear
Boussinesq

free surface elevation (m)

0.02

0.01

0.01

0.02

0.03

10

15

20

25

30

35

time (s)

0.03
linear
Boussinesq

free surface elevation (m)

0.02

0.01

0.01

0.02

0.03

10

15

20

25

30

35

time (s)

0.03
linear
Boussinesq

free surface elevation (m)

0.02

0.01

0.01

0.02

0.03

10

15

20

25

30

35

time (s)

Figure 3.4: The time histories of the free surface elevations at 0 cm, 36 cm and 72 cm away from the left wall.
The excitation period is 0.7 s and the excitation amplitude is 2 mm.

33

0.01
linear
Boussinesq

0.008

free surface elevation (m)

0.006
0.004
0.002
0
0.002
0.004
0.006
0.008
0.01

10

15

20

25
time (s)

30

35

40

45

50

0.01
linear
Boussinesq

0.008

free surface elevation (m)

0.006
0.004
0.002
0
0.002
0.004
0.006
0.008
0.01

10

15

20

25
time (s)

30

35

40

45

50

0.01
linear
Boussinesq

0.008

free surface elevation (m)

0.006
0.004
0.002
0
0.002
0.004
0.006
0.008
0.01

10

15

20

25
time (s)

30

35

40

45

50

Figure 3.5: The time histories of the free surface elevations at 0 cm, 36 cm and 72 cm away from the left wall.
The excitation period is 1.6 s and the excitation amplitude is 2 mm.

34

0.16
linear
Boussinesq

0.158

free surface height (m)

0.156
0.154
0.152
0.15
0.148
0.146
0.144
0.142
0.14

0.1

0.2

0.3
0.4
0.5
0.6
0.7
nondimensional length of tank

0.8

0.9

Figure 3.6: The free surface height in the tank at instant 7.15 s.

0.16
linear
Boussinesq

0.158

free surface height (m)

0.156
0.154
0.152
0.15
0.148
0.146
0.144
0.142
0.14

0.1

0.2

0.3
0.4
0.5
0.6
0.7
nondimensional length of tank

0.8

0.9

Figure 3.7: The free surface height in the tank at instant 7.95 s.

35

0.02
linear
Boussinesq

0.015

free surface elevation (m)

0.01

0.005
0

0.005

0.01
0.015

0.02

10

15
time (s)

20

25

30

0.015
linear
Boussinesq

free surface elevation (m)

0.01

0.005

0.005

0.01

0.015

10

15
time (s)

20

25

30

0.015
linear
Boussinesq

free surface elevation (m)

0.01

0.005

0.005

0.01

0.015

10

15
time (s)

20

25

30

Figure 3.8: The time histories of the free surface elevations at 0 cm, 36 cm and 72 cm away from the left wall.
The excitation period is 0.7 s and the excitation amplitude is 0.5 degree.

36

0.05
linear
Boussinesq

0.04

free surface elevation (m)

0.03
0.02
0.01
0
0.01
0.02
0.03
0.04

10

15
time (s)

20

25

30

0.03
linear
Boussinesq

free surface elevation (m)

0.02

0.01

0.01

0.02

0.03

10

15
time (s)

20

25

30

0.03
linear
Boussinesq

free surface elevation (m)

0.02

0.01

0.01

0.02

0.03

10

15
time (s)

20

25

30

Figure 3.9: The time histories of the free surface elevations at 0 cm, 36 cm and 72 cm away from the left wall.
The excitation period is 0.7 s and the excitation amplitude is 1.0 degree.

37

0.05
linear
Boussinesq

0.04

free surface elevation (m)

0.03
0.02
0.01
0
0.01
0.02
0.03
0.04

10

15
time (s)

20

25

0.03
linear
Boussinesq

free surface elevation (m)

0.02

0.01

0.01

0.02

0.03

10

15
time (s)

20

25

0.03
linear
Boussinesq

free surface elevation (m)

0.02

0.01

0.01

0.02

0.03

10

15
time (s)

20

25

Figure 3.10: The time histories of the free surface elevations at 0 cm, 36 cm and 72 cm away from the left wall.
The excitation period is 1.6 s and the excitation amplitude is 0.5 degree.

38

x 10

Free surface elevation (m)

3
2
1
0
1
2
3
4
5

10

15

20

25
Time (s)

30

35

40

45

50

Figure 3.11: The time histories of the free surface elevations on the position of the wall.

The filling level h=0.1 m and the same length of tank is considered here. The excitation
periods T = 2.1 second and T =2.4 second are selected for the comparisons at the left wall and
36 cm away from the left wall. The motion amplitude is 1 mm. The excitation period 2.1 s is
close to the first natural period. On the left wall (figure 3.14), steeper crests and flatter troughs
are given by nonlinear free surface conditions. On the position of 36 cm away from the left
wall (figure 3.15), the amplitude of crests are lower than linear results. The excitation period
2.4 second is shown in figures 3.16 and 3.17. The free surface elevations are more complex.
Harmonics are clearly present in the Boussinesq model and in the integral equation method.
However differences can be noticed between the two methods and these differences are larger
after the first period of the envelope. This phenomenon can not be predicted by linear theory.
T = 0.73 s
0.04
linear
Boussinesq
BEM

free surface elevation (m)

0.03

0.02

0.01

0.01

0.02

0.03

10
time (s)

12

14

16

18

Figure 3.12: Comparison of time histories of the free surface elevations at the left wall.

39

20

T = 1.83 s
0.04
linear
Boussinesq
BEM

free surface elevation (m)

0.03

0.02

0.01

0.01

0.02

0.03

10

15
time (s)

20

25

30

Figure 3.13: Comparison of time histories of the free surface elevations at the left wall.
3

10

T = 2.1 s

x 10

linear
Boussinesq
BEM

free surface elevation (m)

6
4
2
0
2
4
6
8

10

15

20

25

time (s)

Figure 3.14: Comparison of time histories of the free surface elevations at the left wall. The period is 2.1 second
and the filling level is 0.1 m.

40

T = 2.1 s

x 10

linear
Boussinesq
BEM

free surface elevation (m)

1
0

2
3

10

15

20

25

time (s)

Figure 3.15: Comparison of time histories of the free surface elevations at 36 cm away from the left wall. The
period is 2.1 second and the filling level is 0.1 m.
3

T = 2.4 s

x 10

linear
Boussinesq
BEM

free surface elevation (m)

4
0

10

15

20

25
time (s)

30

35

40

45

50

Figure 3.16: Comparison of time histories of the free surface elevations at the left wall. The period is 2.4 second
and the filling level is 0.1 m.

41

x 10

linear
Boussinesq
BEM

free surface elevation (m)

10

15

20

25
time (s)

30

35

40

45

50

Figure 3.17: Comparison of time histories of the free surface elevations at 36 cm away from the left wall. The
period is 2.4 second and the filling level is 0.1 m.

42

3.7

3.7.1

Validation: numerical experiments for a flat bottom


tank with shallow water and finite water depth conditions
Shallow water case

According to different filling levels in the tank, the regimes of the motions of the liquid can
be divided roughly into shallow water theory (h/l < 0.1), intermediate water depth regime
(0.1 h/l 0.2 0.25) and finite water depth condition (0.2 0.25 < h/l < 1.0) (Faltinsen
and Timokha [2009]). Here l is the length of the tank and h is the filling level. The shallow
water sloshing motions in a rectangular tank are considered and compared with Antuono et al.
[2012]. In this paper, Antuono et al proposed to solve a simple Boussinesq model in terms of
depth-averaged fluid velocity. The results are compared with the experimental data and a SPH
model. Their comparisons proved that the proposed shallow-water model is accurate, fast and
robust. A linear viscous term is added into the momentum equation of Boussinesq model
of the paper.
r
X
einx
n
2h l
4n2 2
=
n Un
; n =
(1 + ) +
(3.89)
2i
2Re
b
h
Re
n

where n is the damping coefficient of sloshing


mode proportional to the modal amplitude Un of
depth-averaged fluid velocity. And Re = (l gh)/ is the Reynolds number. is the kinematic
viscosity which equals to 0.0094 cm2 s1 and b is the length along x-axis. The first term in
the right-hand side takes into account the dissipation due to the boundary layers while the
second one describes the dissipation inside the fluid bulk. In our Boussinesq model, the linear
damping term proposed in section 3.3 is also adopted in the dynamic free surface condition.
Firstly, a rectangular tank, length l = 1.175 m and filling level h = 0.06 m, is considered.
The width of the tank is b = 12 cm. The first resonant frequency based on the linear dispersion
relation is r = 2.0425 s1 . The horizontal harmonic forced motions are taken for the tank
with the amplitude 0.0039 m. The kinematic viscosity is 0.0094 cm2 s1 . The simulations time
was t/T = 200 in the paper to go beyond the transient regime and ensure the establishment of
a steady state. Here T is the excitation period. From figure 3.18 to 3.21, different frequencies
of the forcing excitation are compared respectively between the Boussinesq-type equations
in terms of velocity and the Boussinesq model in terms of velocity potential. The small discrepancy of the phases have been removed for the comparisons of amplitudes through shifting
curves. The non-dimensional time t/T is used for the values of the abscissa. From the figures,
the results given by the Boussinesq model in terms of the velocity potential present an overall
good agreement with the Boussinesq model in terms of the velocity.
The second tank with length l = 0.6095 m and filling level h = 0.06 m is used for our
computations. The ratio h/l = 0.098 is close to the intermediate water depth regime. The
width of the tank b is 0.23 m. The horizontal forcing amplitude is 0.00196 m and the first
resonant frequency predicted through the linear theory wr is 3.893 rad/s. The viscosity term
taken in the numerical model is the same as before. The evolution of the free surface at
the left wall are compared with different excitation frequency /r = 0.97 (figure 3.22) and
/r = 1.04 (figure 3.23). In these two cases, a two-wave regime and a single-bump-wave
regime are observed during one excitation period. When the excitation frequency is close to
the first natural frequency, a nonlinear subharmonic free surface elevation motion can be found
in figure 3.24. The results give an overall good agreement with Antuonos model.
43

0.7
Boussinesq
Antuono

free surface elevation / filling level

0.6
0.5
0.4
0.3
0.2
0.1
0
0.1
0.2

89.5

90

90.5

91
t/T

91.5

92

92.5

Figure 3.18: The evolution of the free surface at the wall with /r = 0.98; T is the period of the forcing
excitation. The length of tank is 1.175 m and filling level is 0.06 m.

0.7
Boussinesq
Antuono

free surface elevation / filling level

0.6
0.5
0.4
0.3
0.2
0.1
0
0.1
0.2

89.5

90

90.5

91
t/T

91.5

92

92.5

Figure 3.19: The evolution of the free surface at the wall with /r = 1.02. T is the period of the forcing
excitation. The length of tank is 1.175 m and filling level is 0.06 m.

44

0.7
Boussinesq
Antuono

free surface elevation / filling level

0.6

0.5
0.4

0.3
0.2
0.1
0

0.1

89.5

90

90.5

91
t/T

91.5

92

92.5

Figure 3.20: The evolution of the free surface at the wall with /r = 1.08; T is the period of the forcing
excitation. The length of tank is 1.175 m and filling level is 0.06 m.

0.7
Boussinesq
Antuono

free surface elevation / filling level

0.6
0.5
0.4
0.3
0.2
0.1
0
0.1
0.2
88

88.5

89

89.5
t/T

90

90.5

91

Figure 3.21: The evolution of the free surface at the wall with /r = 1.10. T is the period of the forcing
excitation. The length of tank is 1.175 m and filling level is 0.06 m.

45

0.7
Boussinesq
Antuono

0.6

free surface elevation / filling level

0.5
0.4
0.3
0.2
0.1
0
0.1
0.2

89.5

90

90.5

91
t/T

91.5

92

92.5

Figure 3.22: The evolution of the free surface at the wall with /r = 0.97. The length of tank is 0.6095 m
and filling level is 0.06 m.

0.7
Boussinesq
Antuono

0.6

free surface elevation / filling level

0.5
0.4
0.3
0.2
0.1
0
0.1
0.2

89.5

90

90.5

91
t/T

91.5

92

92.5

Figure 3.23: The evolution of the free surface at the wall with /r = 1.04. The length of tank is 0.6095 m
and filling level is 0.06 m.

46

0.6
Boussinesq
Antuono

free surface elevation / filling level

0.5
0.4

0.3
0.2

0.1

0
0.1

0.2
80

82

84

86

88

90
t/T

92

94

96

98

100

Figure 3.24: The evolution of the free surface at the wall with /r = 0.99. The length of tank is 0.6095 m
and filling level is 0.06 m.

3.7.2

Finite water depth case

In Wu et al. [1998], a finite element method based on fully non-linear wave potential theory
is used for the simulation of the sloshing motions in a rectangular tank. The viscous term is
not used in this paper. A 2D rectangular tank with l/h = 2 and b/h = 0.2 is adopted. Here l
is the length of the tank, h is the filling level of the tank and b is the width along x-axis. For
our computations, the length of the tank is set to 1.08 m and the filling level is 0.54 m.
The free surface profiles are compared in figure 3.25 and figure 3.26 with different computational models. The non-dimensional amplitude a/h of thepharmonic horizontal motion equals
to 0.0186. The non-dimensional time is defined as = t g/h. The excitation period equals
to 0.999 times the first natural period. The free surface height is measured from the bottom
of the tank. The linear theories given in chapter 2 are used here for comparisons. From the
figure, we can find numerical solutions based on nonlinear free surface conditions coincide well
and the linear solution still gives good results in this case.
In Frandsen [2004], liquid motions in 2-D tanks excited by periodic loadings are considered.
A fully nonlinear numerical model is based on inviscid flow equations and solutions are obtained
using finite differences. The sloshing motion is discussed based on potential theory. A 2D
rectangular tank with aspect ratio h/l = 0.5 is considered. Here the length of the tank is set to
2.0 m and the filling level equals to 1.0 m. The non-dimensional excitation amplitude equals
to a 2 /g where a is excitation amplitude and is excitation frequency. Two kind of external
forcing excitation frequencies /r = 0.7 and /r = 1.3 are taken in our comparisons. The
first resonance frequency r is calculated by the linear dispersion relation. According to the
linear dispersion equation, the first natural frequency is 3.7594 s1 .
47

linear
Boussinesq
BEM
Wu

free surface height / filling level

1.3
1.2
1.1
1
0.9
0.8
0.7
0.6

0.1

0.2

0.3

0.4
0.5
0.6
0.7
Nondimensional tank length

0.8

0.9

Figure 3.25: The comparison of the free surface elevations at the non-dimensional time = 13.0667.

1.4
linear
Boussinesq
BEM
Wu

free surface height / filling level

1.3

1.2

1.1

0.9

0.8

0.7

0.1

0.2

0.3

0.4
0.5
0.6
0.7
Nondimensional tank length

0.8

0.9

Figure 3.26: The comparison of the free surface elevations at the non-dimensional time = 15.725.

48

4
Boussinesq
Frandsen
3

elevation / amplitude

10

20

30

40
50
60
nondimensional time

70

80

90

100

Figure 3.27: The free surface elevations on the left wall. /r = 0.7 and = 0.0036.

4
Boussinesq
Frandsen
3

elevation / amplitude

10

20

30

40
50
60
nondimensional time

70

80

90

100

Figure 3.28: The free surface elevations on the left wall. /r = 0.7 and = 0.036.

49

Boussinesq
Frandsen

elevation / amplitude

10

20

30

40
50
60
nondimensional time

70

80

90

100

Figure 3.29: The free surface elevations on the left wall. /r = 1.3 and = 0.0036.

Boussinesq
Frandsen

elevation / amplitude

10

20

30

40
50
60
nondimensional time

70

80

90

100

Figure 3.30: The free surface elevations on the left wall. /r = 1.3 and = 0.072.

50

In figure 3.27 and 3.28, the frequency /r = 0.7 is used in the horizontal forcing motions
with different amplitudes = 0.0036 and = 0.036. The non-dimensional time is defined
as t r . The free surface elevations on the left wall divided by the amplitude of the forcing
motion give the values of the ordinates. In figure 3.29, the forcing frequency ratio is /r = 1.3.
The same excitation frequency and more stronger motions = 0.072 is shown in figure 3.30.
From the figures, the free surface elevations present an overall good agreement between the
Boussinesq model and the Frandsens model.
Another discussion in the range of finite depth condition can be found in Faltinsen et al.
[2000]. In this paper, a multidimensional modal analysis of nonlinear sloshing in a rectangular
tank is given based on the potential theories. Some analytical and experimental results of
the sloshing motions in a two dimensional rectangular tank are presented. The left figure in
figure 3.31 adopts a rectangular tank with length l = 1.73 m and filling level h = 0.6 m. The
harmonic horizontal motions are used with an amplitude equals to 0.025 m. The period of
the forcing excitation is 1.3 s. The right figure takes the same length of the tank but different
filling level h = 0.5m. The amplitude of the forced motion is 0.029 m and the period of the
forcing excitation is 1.4 s. The free surface elevations at 0.05 m away from the left wall are
compared between Boussinesq model and multidimensional modal in both figures.
The amplitude calculated by the Boussinesq model are much lower than the results in
Faltinsen et al. [2000]. And much stronger difference can be found with larger amplitudes of
motions. The phenomenon maybe is due to the boundary condition on the wall used in our
numerical model. The slope of tangent of the free surface profile close to the wall is not equal
to zero any more and the limitation should be revised for large amplitude of sloshing motions.
The problem could be solved by the moving boundary algorithm. The free surface profile will
be extended linearly to the outside of the tank and the slope of the external linear free surface
is determined by the internal free surface close to the wall. The details can be found in the
chapter Concluding remarks and perspectives.
0.5

0.2
Boussinesq
Faltinsen

0.15

Boussinesq
Faltinsen
0.4

free surface elevation (m)

free surface elevation (m)

0.1

0.05

0.05

0.3

0.2

0.1

0.1
0.1

0.15

0.2

10

0.2

15

time (s)

10
time (s)

Figure 3.31: The free surface elevations on the position of 0.05 m away from the left wall.

51

15

3.8

The numerical models of the sloshing in a moving


inclined bottom rectangular tank

A 2D rectangular tank of length l with an inclined bottom of bathymetry h(y) is adopted.


A tank-fixed coordinate system Oyz lies on the free surface at rest in the middle of the tank
and only the sway motion vy is considered here.
The governing equations of the liquid in the tank are given:
l
y = ; z + hy y = ~v ~n z = h(y);
2
2 + 1 w 2 (1 + 2 )
y y + (1 + 2 )w;
t = vy
y + vz w g 1
t = vy y

y
y
2 y 2
yy + zz = 0;

y = vy

(3.90)
(3.91)

w are the velocity potential and


where ~n is the normal vector of the inclined bottom and ,
vertical velocity defined directly on the free surface. The process of solving the problem is the
same as the case of flat bottom tank. Here is the computational process of the linear system.
Firstly, the velocity potential and the vertical derivative of the velocity potential are given as
follows:
+ (J02 +
(y, z, t) = (J01 +
z J11p )
z J12p )w
(3.92)
+ (J01 +
w(y, z, t) = (J02 2
z J12w )
z J11w )w

(3.93)

on the free surface is given by (the notation F S expresses that the operators
The potential
have been evaluated with z = ):
FS
FS
FS
FS
= (J01
+ (J02

+
z J11p
)
+
z J12p
)w

Based on the assumption h = O() ( =

= z

(3.94)

) on the bottom (Bingham et al. [2009]):


y

+ (J02 +
= (J01 +
z J11p )
z J12p )w
+ (J02 + (
= (J01 + (
z J11p ))
z J12p ))w
J01 = J01 +
z

J01 ;
z

J02 = J02 +
z

J02
z

Retaining the terms at O():


b
b
b
b
b
+ (J b +
z + hy y = (J02
2
z J12w
+ h J01
)
z J11w
+ h J02
)w
01
(3.95)

and
Together with velocity potential on the free surface, the linear system with unknowns
w is given as follows. The notation b expresses that the operators have been evaluated with
z = h(y).
   

FS
FS
FS
FS

J01
+
z J11p

J02
+
z J12p

(3.96)
=
b
b
b
b
2
b
b
0
z J12w + h J01 J01 +
z J11w + h J02 w
J02
The finite difference scheme is adopted for the spatial derivatives in the numerical computations.
52

Based on the assumption of slowly varying bathymetry, the following equations are used
for the calculations of the derivatives of the water depth.
h1 = h3 , h0 = h2 , (hy)1 = (hy)3, (hy)0 = (hy)2

(3.97)

with the same the conditions on the walls as the flat bottom, we can get
1 = 3 , 0 = 2

w1
= w3 , w0 = w2

(3.98)

These conditions are just an approximate treatment for the small slope of the varying bottom.

3.9

Comparisons between flat bottom tank and inclined


bottom tank

A rectangular tank with slowly varying bottom is considered here. The water depth varies
linearly with a bottom slope equals to 0.074 (see Figure 3.32 right). In order to compare with
the flat bottom tank, the length of the tank is set to 1.08 m. The water depth in the left side
is 0.19 m and the right side is 0.11 m. So the mass of the liquid is the same as in the flat
case with the filling level h=0.15 m. Forced harmonic sway motions are considered for the
comparisons. The periods are from 0.6 second to 2.0 seconds and the amplitude is 1 mm.
The Boussinesq models are used for the simulation of the sloshing motions in the tanks. This
means the nonlinear free surface conditions will be adopted. The viscosity is not considered
here. According to the linear theories in the chapter 2, there is a peak if the tank is not
symmetric. Based on the nonlinear free surface conditions, the symmetric and asymmetric
rectangular tanks will be compared for the first, second and third sloshing modes.

0.19

0.15

Figure 3.32: Both rectangular tanks are forced by harmonic sway motions with amplitude of 1 mm. The hight
of the inclined bottom is 0.08 m with slope 0.074.

In figure 3.33, the maximum free surface elevations during 40 seconds (consistent with
experimental data) divided by the amplitude of motion are calculated. Due to the slope
of the bottom, the amplitudes of the peaks on the first sloshing mode (1.82 s) and the third
sloshing mode (0.73 s) are a little higher and natural frequencies change a little. The significant
differences can be found in the second sloshing mode which means one wave length in the
tank. For the flat bottom, the motions of the liquid is anti-symmetrical and the free surface
elevations are small. But for the inclined bottom, the asymmetrical shape of the tank gives
larger amplitude of free surface elevations and so there is a small peak in the second sloshing
mode.
Figure 3.34 gives the time histories of the free surface elevations at the wall. The red line
is from the solutions of the inclined bottom tank and the blue line is the elevation in the flat
53

0.11

bottom tank. Significantly, there is a resonance in the inclined bottom when the period is 1.0
second.
Figure 3.35 and 3.36 are the results for the inclined bottom tank. The free surface elevations
on the left side and the right side of the tank are compared during 40 seconds. Due to the
asymmetric shape of the tank, the free surface elevations on the right side of the tank are
higher than the left side of the tank. Due to the anti-phase of the two sides, one curve has
been shifted by half of the period for the comparisons of the amplitudes.
70

Maximum elevations / amplitude of motion

flat
slope
60
X: 0.74
Y: 60.36

50

X: 1.83
Y: 43.33

X: 0.73
Y: 48.31

40
X: 1.82
Y: 38.77

30

20
X: 1
Y: 8.545

10

0.8

1.2
1.4
Period (s)

1.6

1.8

Figure 3.33: The figure gives a comparison between the flat bottom tank and the inclined bottom tank. The
maximum free surface elevations during 40 seconds divided by the amplitude of motion gives the values of
ordinates.

x 10

flat
inclined

Free surface elevations (m)

3
2
1
0
1
2
3
4
5

10

15

20
Time (s)

25

30

35

40

Figure 3.34: The period of the forced harmonic motion is 1.0 second. This figure gives a comparison of free
surface elevations on the right wall between flat bottom tank and inclined bottom tank.

54

The period is 0.75 second


0.02
right
left

Free surface elevations (m)

0.015

0.01

0.005

0.005

0.01

0.015

10

15

20
Time (s)

25

30

35

40

Figure 3.35: The comparisons of the free surface elevations between the left side and the right side.

The period is 1.8 second


0.025
right
left

0.02

Free surface elevations (m)

0.015
0.01
0.005

0
0.005
0.01

0.015

10

15

20
Time (s)

25

30

35

40

Figure 3.36: The comparisons of the free surface elevations between the left side and the right side.

55

3.10

Experimental campaigns

Forced motions are simulated by a dynamic hexapod which is designed by the company
SYMETRIE. The hexapod, MISTRAL (figure 3.37), gives the motion and simulation of high
payloads following 6 dof (Rx,Ry,Rz,Tx,Ty,Tz). The robust mechanical structure of the MISTRAL hexapod allows payloads up to 1 ton. The dynamic forces and moments can be measured
by the sensors embedded in hexapod directly. All the input and output signals of movements,
forces and moments can be operated by controlling computer. The coordinates system of the
hexapod is given in figure 3.38. And the range of the motions of the hexapod can be found in
the table 3.1.

Figure 3.37: The MISTRAL hexapod from http://www.symetrie.fr.

Figure 3.38: The Oxyz coordinates system. The translational motions Tx, Ty and Tz are along the x,y and z
axis. The rotational motions Rx, Ry and Rz are around former three axis respectively.

56

Because the numerical computations are two dimensional, in our experiments, a tank with
internal length 1.08 m and width 0.1 m is adopted. The material of the tank wall is made of
glass which can be considered as smooth. The glass walls are connected by a metal frame and
a smooth PVC plate is used for the sloping bottom. In figure 3.39, the tank is fixed on the
hexapod with longer side along the y-axis of the hexapod. The liquid in the tank is dyed with
fluorescein and a laser light is put on the top of the tank. So the free surface of the liquid in
the tank can be en-lighted due to the laser sheet. A high-speed camera, fixed on the hexapod,
is set in the front of the tank. An opaque sheet is used to stop the outside lights in order to
enhance the automatic detection of the free surface by the camera.
Axis
Linear travel Tx
Linear travel Ty
Linear travel Tz
Angular travel Rx
Angular travel Ry
Angular travel Rz

Travel range
460 mm
460 mm
400 mm
30
30
40

Speeds
1000 mm/s
1000 mm/s
700 mm/s
50 /s
50 /s
70 /s

Accelerations
10000 mm/s2
10000 mm/s2
8000 mm/s2
500/s2
500/s2
700/s2

Table 3.1: The range of the motions of the hexapod

Laser sheet

High lighted free surface

The viewfield of the camera

Figure 3.39: The left figure is the experimental set-up for the inclined bottom. The right figure is a sketch of
the setup.

57

From the setup shown in figure 3.39, the red dashed line corresponds to the lighted free
surface by the laser sheet. Due to the limited experimental condition, the laser sheet is not
wide enough to catch the whole free surface in the tank. The two ended sections of the free
surface are not considered. Through the high-speed camera, the free surface can be captured
clearly. Figure 3.40 gives some images of the free surface at different instants recorded by the
camera.

Figure 3.40: The free surfaces in the tank at different instants.

The monochrome images include the informations of the free surface elevations. The function im2bw can be used to convert the monochrome images into the binary images based on
the software Matlab. In the binary images, the pixel points of the free surface can be distinguished easily from the black background. Based on the initial free surface at rest, the free
surface elevations of each instant can be obtained.

3.11

Comparisons between numerical solutions and experimental results

3.11.1

Flat bottom tank

The rectangular flat bottom tank, length l=1.08 m and filling level h=0.15 m, is used.
The forced harmonic horizontal motions are taken with amplitude 2 mm. Both the numerical
results based on the Boussinesq model and the integral equation method will be compared with
the experimental results.
Due to the fact that the laser sheet is not wide enough to catch the whole free surface in
the tank, part of the free surface is used, which means 6.37 cm counted from the left wall and
12.05 cm counted from the right wall are not included. Figure 3.41 gives the time histories of
the free surface elevations at 6.37 cm away from the left wall during 30 seconds. The period
of the forced motions is 1.8356 second which means the first sloshing mode of the tank. Both
the Boussinesq model and the integral equation method are compared to the experimental
results. There is no damping terms in the numerical models so the values of the amplitudes
are somewhat higher than the experimental results.

58

T = 1.8356 s
0.06
Boussinesq
BEM
Experiment

0.05

free surface elevation (m)

0.04
0.03
0.02
0.01
0
0.01
0.02
0.03

10

15
time (s)

20

25

30

Figure 3.41: The time histories of the free surface elevations at 6.37 cm away from the left wall. The bottom
of tank is horizontal.

0.02

0.025
Boussinesq
BEM
Experiment

0.015

Boussinesq
BEM
Experiment

0.02

free surface elevation (m)

free surface elevation (m)

0.015
0.01

0.005

0.005

0.01

0.005
0
0.005

0.01

0.01

0.015

0.015

0.1

0.2

0.3

0.4
0.5
0.6
0.7
length of tank (m)

0.8

0.9

0.03

0.1

0.2

0.3

0.4
0.5
0.6
0.7
length of tank (m)

0.8

0.9

0.04
Boussinesq
BEM
Experiment

0.025

Boussinesq
BEM
Experiment

0.03

free surface elevation (m)

free surface elevation (m)

0.02
0.015
0.01
0.005
0

0.02

0.01

0.005
0.01
0.01
0.015

0.1

0.2

0.3

0.4
0.5
0.6
0.7
length of tank (m)

0.8

0.9

0.02

0.1

0.2

0.3

0.4
0.5
0.6
0.7
length of tank (m)

0.8

0.9

Figure 3.42: The free surface profiles in the tank corresponding to different instants: 15.11 s (upper left); 16.94
s (upper right); 18.76 s (down left); 20.58 s (down right).

59

T = 0.7306 s
0.12
Boussinesq
BEM
Experiment

0.1

free surface elevation (m)

0.08
0.06
0.04
0.02
0
0.02
0.04
0.06

10

15
time (s)

20

25

Figure 3.43: The time histories of the free surface elevations at 36 cm away from the left wall. The bottom of
the tank is horizontal.

The free surface profiles in the tank are given in figure 3.42 for different instants: 15.11 s
(upper left), 16.94 s (upper right), 18.76 s (down left) and 20.58 s (down right) corresponding
to the eighth, ninth, tenth and eleventh peaks in the figure 3.41. The period of the forced
motion is close to the first sloshing mode and half of the wavelength can be found in the tank.
From the figures, we can find that both numerical models give almost the same results. The
amplitude of the BEM model is higher for crest amplitude. The comparison of free surface
elevations with excitation period 0.7306 second is given in figure 3.43. The same conclusions
can be drawn as the previous case.

3.11.2

Inclined bottom tank

The comparisons of inclined bottom tank between the experimental results and the Boussinesq model are given. The filling level is 19 cm and the bathymetry varies linearly with a
bottom slope of 0.074. The excitation motions are harmonic sway motions with amplitude 1
mm. Due to the limitations in the experiments, the results are compared for the same window
as the experimental data. In figure 3.44, the maximum free surface elevations during 40 seconds
divided by the amplitude of motion give the values of ordinates. The linear damping term is
not considered here. The numerical results present an overall agreement with experimental
results. The viscosity effect is not very important during short time. Figure 3.45 gives the
time histories of the free surface at 6.37 cm away from the left wall during 40 seconds. The
period is 1.0 second which is close to the second sloshing mode. The excitation periods 0.73
second and 1.83 second are shown in figures 3.46 and 3.47. These are close to the first and
third sloshing modes respectively. We can see clearly that the numerical code overestimate the
wave elevations.
60

45
numerical
experimental

Maximum elevations / amplitude of motion

40
35
30
25
20
15
10
5
0
0.6

0.8

1.2
1.4
Period (s)

1.6

1.8

Figure 3.44: The amplitude of the motion is 1 mm. The free surface is from 6.37 cm away from left wall to
12.05 cm away from right wall.

x 10

Boussinesq
Experiment

Free surface elevations (m)

10

15

20
Time (s)

25

30

35

40

Figure 3.45: The comparison between the numerical and the experimental results at 6.37 cm away from the
left wall. The period is 1 second.

61

0.02
Boussinesq
Experiment
0.015

Free surface elevations (m)

0.01

0.005
0

0.005
0.01
0.015

0.02

10

15

20
Time (s)

25

30

35

40

Figure 3.46: comparison between the numerical and the experimental results at 6.37 cm away from the left
wall. The period is 0.73 second.

0.025
Boussinesq
Experiment

0.02

Free surface elevations (m)

0.015
0.01
0.005
0
0.005
0.01
0.015
0.02

10

15

20
Time (s)

25

30

35

40

Figure 3.47: The comparison between the numerical and the experimental results at 6.37 cm away from the
left wall. The period is 1.83 second.

62

Next, the harmonic horizontal motions are taken with excitation amplitude 2 mm. In figure
3.48, the period of the forced motions is 1.0 s corresponding to the second sloshing mode. The
figure presents the comparisons of the free surface elevations in the middle of the tank. The
period of the forced motion in figure 3.49 is 1.836 s corresponding to the first sloshing mode.
The comparisons of the free surface elevations at 6.37 cm away from the left wall are presented.
The subharmonic motions have been simulated in the troughs of the free surface elevations but
the amplitudes are much larger due to no viscous assumption. The free surface profiles in the
tank are given in figure 3.49 for different instants: 11.448 s (upper left), 13.284 s (upper right),
15.124 s (down left) and 16.956 s (down right) corresponding to the sixth, seventh, eighth and
ninth peaks in figure 3.50.
In figure 3.51 and 3.52, the period of the forced motions is 0.73 s corresponding to the
third sloshing mode. The free surface elevations at 6.37 cm and 36 cm away from the left wall
are shown in figures. The nonlinear effects, the steeper wave crest and more flat trough, can
also be found in these figures. The numerical results present an overall good agreement with
experimantal results.
T = 1.0 s
0.015
Boussinesq
Experiment

free surface elevation (m)

0.01

0.005

0.005

0.01

10

15
time (s)

20

25

30

Figure 3.48: The time histories of the free surface elevations in the middle of the tank. The bottom is inclined.
The excitation amplitude is 2 mm.

The linear damping term (section 3.3) is tested. The flat bottom tank with filling level 0.15
m is used. The numerical results without damping have been compared in figures 3.41 and 3.43.
The excitation amplitude is 2 mm. The first sloshing mode and the third sloshing mode are
shown in figure 3.53 and 3.54. The value of viscosity is set to 1105 . From the comparisons,
the numerical results present an overall good agreement with experimental results.

63

T = 1.836 s
0.05
Boussinesq
Experiment

0.04

free surface elevation (m)

0.03

0.02

0.01
0
0.01

0.02
0.03

10

15
time (s)

20

25

Figure 3.49: The time histories of the free surface elevations at 6.37 cm away from the left wall. The bottom
is inclined. The excitation amplitude is 2 mm.

0.02

0.02
Boussinesq
Experiment

0.015

0.015

0.01

0.01

free surface elevation (m)

free surface elevation (m)

Boussinesq
Experiment

0.005

0.005

0.01

0.015

0.005

0.005

0.01

0.1

0.2

0.3

0.4
0.5
0.6
0.7
length of tank (m)

0.8

0.9

0.015

0.025

0.1

0.2

0.3

0.4
0.5
0.6
0.7
length of tank (m)

0.8

0.9

0.03
Boussinesq
Experiment

Boussinesq
Experiment

0.025

0.02

0.02

free surface elevation (m)

free surface elevation (m)

0.015
0.01

0.005
0

0.015
0.01
0.005
0
0.005

0.005
0.01
0.01
0.015

0.015

0.1

0.2

0.3

0.4
0.5
0.6
0.7
length of tank (m)

0.8

0.9

0.02

0.1

0.2

0.3

0.4
0.5
0.6
0.7
length of tank (m)

0.8

0.9

Figure 3.50: The free surface profiles in the tank corresponding to different instants: 11.448 s (upper left);
13.284 s (upper right); 15.124 s (down left); 16.956 s (down right). The excitation amplitude is 2 mm.

64

T = 0.73 s
0.06
Boussinesq
Experiment

0.05

free surface elevation (m)

0.04
0.03
0.02
0.01
0
0.01
0.02
0.03
0.04

10

15
time (s)

20

25

30

Figure 3.51: The time histories of the free surface elevations at 6.37 cm away from the left wall. The bottom
is inclined. The excitation amplitude is 2 mm.

T = 0.73 s
0.08
Boussinesq
Experiment

free surface elevation (m)

0.06

0.04

0.02

0.02

0.04

10

15
time (s)

20

25

30

Figure 3.52: The time histories of the free surface elevations at 36 cm away from the left wall. The excitation
amplitude is 2 mm.

65

T = 1.8356 s
0.05
Boussinesq
Experiment

free surface elevation (m)

0.04

0.03

0.02

0.01

0.01

0.02

10

15
time (s)

20

25

30

Figure 3.53: The point is located 6.37 cm away from the left wall. The period of forced motion is 1.8356 second
and the amplitude is 2 mm. The viscosity is taken equal to 1105 .

T = 0.7306 s
0.08
Boussinesq
Experiment

free surface elevation (m)

0.06

0.04

0.02

0.02

0.04

10

15

20

25

time (s)

Figure 3.54: The point is located 36 cm away from the left wall. The excitation period is 0.7306 second and
motion amplitude is 2 mm. The viscosity is taken equal to 1105 .

66

The hexapod is equiped with 6 mono-axial force sensors mounted between the actuators
and the ball joints fixed on the top plate. To obtain the overall forces, a calculation including
the precise position of the hexapod is needed. Therefore in order to ensure the correctness
of measured results, a 6 axis balance is placed between the tank and hexapod to control the
efforts. As shown in figure 3.55, the tank in the left is placed on hexapod directly and the
efforts are measured by hexapod. In the right-hand, the tank is fixed on the balance which is
fixed on hexapod. The efforts can be measured by hexapod and balance simultaneously.

tank
tank
balance
hexapod

hexapod

Figure 3.55: The different ways to measure the efforts produced by the liquid in the tank.

Figure 3.56 shows the comparisons between the experimental results and numerical results
for a flat bottom tank. For the experimental results, the efforts measured by hexapod without
balance are shown by only hexapod. And hexapod with balance represents efforts measured by
hexapod when the balance is fixed. The results given by balance are shown by balance. The
harmonic horizontal motions are adopted with excitation amplitude 1 mm. The length of tank
is 1.08 m and filling level is 0.15 m. The width of tank is 0.1 m which provides two dimensional
sloshing motions. The excitation periods are from 0.6 second to 2.0 second. For each period,
the maximum horizontal forces during 100 seconds give the values of ordinates. As shown in
figure, the three different measurements give an overall agreement. The experimental results
are compared with numerical solutions based on Boussinesq model and linear theories. It
should be noted that numerical results multiplied by the width of tank are compared with
experimental results. The linear damping terms are considered in both numerical models.
Although the linear damping term of Boussinesq model still have a problem, the forces given
by Boussinesq model present relative good agreement.
The comparisons between numerical results and experimental results with excitation amplitude 2 mm are shown in figure 3.57. As shown in figure, the efforts given by linear theory are
much higher than experimental results. In order to decrease the peak values of linear theory,
higher value of viscosity = 3 106 is used. In both figures, part of the measured efforts are
negative before the third sloshing mode. And it can also be found in the following two figures.
No explanation is given on this behavior and further researches are needed.
The comparisons of efforts with inclined bottom tank are compared in figures 3.58 and
3.59. The filling level is 19 cm. Due to the asymmetric shape of tank, resonance can be found
around the second sloshing mode. The results based on Boussinesq model give an overall
agreement with experimental results. The forces produced by inclined bottom tank are a little
smaller than that produced by flat bottom tank with the same mass of liquid. That means the
horizontal acceleration of liquid is smaller due to the inclined bottom in the tank. However
the forces are overestimated by the model for the first sloshing mode. This is probably due to
the wrong calculation of the damping term.
67

bottom: flat; motion: sway; amplitude: 1mm


14
only hexapod
hexapod with balance
balance
Boussinesq
linear

12

Effort (N)

10
8
6
4
2
0
2
0.6

0.8

1.2
1.4
Period (s)

1.6

1.8

Figure 3.56: The comparisons of horizontal efforts with the case of flat bottom tank. The excitation amplitude
is 0.001 m.

bottom: flat; motion: sway; amplitude: 2mm


25
only hexapod
hexapod with balance
balance
Boussinesq
linear1.E06
linear3.E06

20

Effort (N)

15
10
5
0
5
0.6

0.8

1.2
1.4
Period (s)

1.6

1.8

Figure 3.57: The comparisons of horizontal efforts with the case of flat bottom tank. The excitation amplitude
is 0.002 m.

68

bottom: slope; motion: sway; amplitude: 1mm


10

only hexapod
hexapod with balance
balance
Boussinesq

Effort (N)

10
0.6

0.8

1.2
1.4
Period (s)

1.6

1.8

Figure 3.58: The comparisons of horizontal efforts with the case of inclined bottom tank. The excitation
amplitude is 0.001 m.

bottom: slope; motion: sway; amplitude: 2mm


14
only hexapod
hexapod with balance
balance
Boussinesq

12
10

Effort (N)

8
6
4
2
0
2
4
0.6

0.8

1.2
1.4
Period (s)

1.6

1.8

Figure 3.59: The comparisons of horizontal efforts with the case of inclined bottom tank. The excitation
amplitude is 0.002 m.

69

3.12

Conclusion

The computational models of sloshing motions in tanks are built based on extended Boussinesqtype model. Firstly, the flat bottom rectangular tank was considered with different filling levels.
For shallow water condition case, the Boussinesq model is validated by the results calculated
by Boussinesq-type equations in terms of depth-averaged fluid velocity published in Antuono
et al. [2012]. The linear damping term is considered here. The free surface elevations calculated by these two models coincide well. In finite depth conditions, finite element method,
finite difference method and analytical method are used for comparisons. A good results can
be found for the Boussinesq model with small excitation amplitude. When the free surface
elevations are too high, the slope of tangent of free surface profile close to the walls can not be
considered as zero. This numerical model can not work well with large excitation amplitude
motions.
For the intermediate water depth, the numerical results calculated by Boussinesq model
are compared to the solutions of linear modal approach firstly. For small excitation amplitude
1 mm and off-resonance excitation period, the linear theories still give a good agreement with
Boussinesq model which includes nonlinear free surface conditions. When excitation amplitude
is 2 mm, free surface profiles present nonlinearity compared with the linear modal approach.
The free surface elevations present steeper crests and flatter troughs in the Boussinesq model
based on nonlinear free surface conditions. The integral equation method is adopted with fully
nonlinear free surface conditions. The two nonlinear numerical methods present an overall
good agreement.
The flat and inclined bottom tank were compared to Boussinesq model. According to linear theory, there is no resonance on the second sloshing mode due to the symmetric shape.
The asymmetric shape, inclined bottom tank, was considered under horizontal harmonic motions. Resonance of the second sloshing mode appears in the numerical results calculated by
Boussinesq model.
The numerical results were validated by experimental results with flat bottom tank and
inclined bottom tank. An overall agreement is given between Boussinesq model and experimental results. The Boussinesq model can be used for the simulation of sloshing motions
with flat bottom tank and inclined bottom tank (small slope). To increase the bottom slope,
terms of order 2 must be taken into account. The efforts produced by the liquid in the tank
are compared between numerical results and experimental results, which presents an overall
agreement.

70

Chapter 4
The transformation of the equations of
motions between the frequency domain
and the time domain
In this chapter, the linearized boundary value problem of floating body, a rectangular
barge, will be introduced in the framework of potential theory. The problem is solved by a
semi-analytical method, eigen-function expansions. The fluid domain is divided into a succession of rectangular sub-domains where eigen-function expansions are used to express the
velocity potential. Matching of the velocity potentials and horizontal velocities at the successive boundaries gives the solutions of the problem. The hydrodynamic coefficients of the
barge, such as the diffraction loads, added masses and damping matrices, can be obtained and
then are substituted into the computations of the equations of motions in the time domain.
The response amplitude operators (RAOs) of the motions of the barge calculated in the time
domain are compared with the results in the frequency domain.

71

z
free surface

barge

1y = 2y
1

~n

2y = 3y

h
3

2 = 3

1 = 2
1z = 0

free surface

2z = 0

y1

y2

3z = 0

Figure 4.1: 2D sketch of a barge floats in waves with constant water depth.

4.1

The boundary value problem for a rectangular floating barge in the linear frequency domain theory

As shown in figure 4.1, the fluid domain is limited by the barge and the horizontal bottom.
The origin of the coordinate system oyz is located on the free surface at rest with z-direction
upwards. The linearized motions of the barge are excited by a harmonic incident wave propagating from the left-hand side semi-infinite sub-domain, along the y-axis. The velocity field in
time harmonic can be represented by:
(y, z, t) = {(y, z)eit }

(4.1)

where is the angular frequency of the incident wave. The velocity potential can be decomposed into the classical form:
(y, z) = i

X
Ag
(I (y, z) + D (y, z)) +
iXk Rk (y, z)

k=2,3,4

(4.2)

where A is the amplitude of the incident wave and g is the gravity acceleration. I denotes
the incident wave potential, and D is the diffraction wave potential due to the barge. The
last term Rk are the radiation potentials associated with the amplitudes of sway (X2 ), heave
(X3 ) and roll (X4 ).
The incident wave potential in the up-wave region (y < y1 ) is given:
I (y, z) =

cosh k0 (z + h) ik0 y
e
cosh k0 h

(4.3)

where h is the constant water depth and k0 satisfies the dispersion relationship 2 = gk0 tanh k0 h.
The governing equations for the velocity potential = I + D are given by:
= 0
~n = 0
z = 0
gz 2 = 0
72

fluid domain
wetted surface
z = h
z=0

(4.4)

And the radiation conditions should be added in the up-wave and down-wave regions. For
y y1 , it contains the incoming wave and a reflected wave; However a transmitted wave
should be included for y y2 .
The radiation potentials of the motions satisfy the following equations:
Rk = 0 fluid domain
Rk ~n = Nk wetted surface
Rkz = 0 z = h
2
gRkz Rk = 0 z = 0
+ propagation to the infinity

(4.5)

where Nk is the generalized normal vector N2 = ny , N3 = nz and N4 = (z zG )ny (y yG )nz .

4.2

Eigen-function expansions and Matching adjacent


regions

Following Cointe et al. [1990] and Liu [2010], the eigen-function expansions of the potential
in each region can be written as follows:
For the left-hand side fluid domain,
1 =

X
cosh k0 (z + h) ik0 y
(e
+ A0 eik0 y ) +
An ekn (yy1 ) cos kn (z + h)
cosh k0 h
n

where kn satisfy 2 = gkn tan kn h.


For the fluid domain under the barge,
X
y y2
y y1
2 = B0
+ C0
+
cos n (z + h)[Bn en (yy1 ) + Cn en (yy2 ) ]
y1 y2
y2 y1
n

(4.6)

(4.7)

n
and d is the draft of the barge.
hd
And for the right-hand side fluid domain,

where n =

3 =

X
cosh k0 (z + h)
D0 eik0 y +
Dn ekn (yy2 ) cos kn (z + h)
cosh k0 h
n

(4.8)

The A0 represents the reflected wave amplitude coefficient and D0 denotes the transmitted
wave amplitude coefficients. The coefficients An and Dn are the amplitude functions for the
evanescent modes at the boundaries, which decay exponentially with distance from the boundary.
At the boundaries, based on the principle of continuity on y = y1 and y = y2 , the following
equations should be satisfied with h z d.
1 = 2
2 = 3

1y = 2y
2y = 3y

(4.9)
(4.10)

And in the range of d z 0, the horizontal derivatives of the potentials 1y and 3y are
equal to zero on y = y1 and y = y2 respectively.
73

At the connection y = y1 , the orthogonal functions cosh k0 (z + h) and cos kn (z + h) are used
for multiplying the equations 4.9 respectively and then the integration along the z-direction
are taken:
Z d
Z d
1 cosh k0 (z + h) dz =
2 cosh k0 (z + h) dz
(4.11)
h

1y cosh k0 (z + h) dz =

1 cos kn (z + h) dz =

1y cos kn (z + h) dz =

2y cosh k0 (z + h) dz

(4.12)

2 cos kn (z + h) dz

(4.13)

2y cos kn (z + h) dz

(4.14)

h
d

Based on the equation 4.11 and 4.13, a linear system of n+1 equations can be given as follows:
~ 1 + MA1 A
~ = MB1 B
~ + MC1 C
~
R

(4.15)

And a similar linear system could be found with the equations 4.12 and 4.14:
~ 2 + MA2 A
~ = MB2 B
~ + MC2 C
~
R

(4.16)

~ B
~ and C
~ are the unknown coefficient vectors. The coefficient matrices MA1 , MA2 ,
where the A,
~ 1, R
~ 2 are known. Eliminating the unknown vector A,
~ two linear
MB1 , MB2 , MC1 , MC2 and R
systems 4.15 and 4.16 can be integrated:
1
~ + (MC2 MA2 M 1 MC1 )C
~ =R
~ 2 MA2 M 1 R
~
(MB2 MA2 MA1
MB1 )B
A1
A1 1

(4.17)

At the connection y = y2 , the orthogonal functions cosh k0 (z + h) and cos kn (z + h) are


also used for multiplying the equations 4.10 respectively and then the integration along the
z-direction are taken:
Z d
Z d
2 cosh k0 (z + h) dz =
3 cosh k0 (z + h) dz
(4.18)
h

2y cosh k0 (z + h) dz =

2 cos kn (z + h) dz =

2y cos kn (z + h) dz =

3y cosh k0 (z + h) dz

(4.19)

3 cos kn (z + h) dz

(4.20)

3y cos kn (z + h) dz

(4.21)

h
0

Based on the equation 4.18 and 4.20, a linear system of n+1 equations can be given as follows:
~ + MC3 C
~ = MD1 D
~
MB3 B

(4.22)

And the similar linear system is given with the equations 4.19 and 4.21.
~
~ = MD2 D
~ + MC4 C
MB4 B
74

(4.23)

~ is the unknown coefficient vector and MB3 , MB4 , MC3 and MC4 are the coefficient
where D
~ two linear systems can be integrated:
matrices. After eliminating the vector D,
1
~ + (MC3 MD1 M 1 MC4 )C
~ =0
(MB3 MD1 MD2
MB4 )B
D2

(4.24)

~ and C
~ can be given and then A,
~ D
~
Together equations 4.17 and 4.24, the values of vector B
can also be solved.
The similar decompositions and matching process is used for the radiation potentials. Due
to the differences of the governing equations between the diffraction potential and radiation potential, a particular solution should be added in the eigen-function expressions. The particular
solution of the radiation potential under the barge are given by

0,
k = 2, sway,

(z+h)2 y2
,
k = 3, heave,
Rk =
(4.25)
2(hd)

(yy
)
(yy
)
G

{ 3G (z + h)2 }, k = 4, roll.
2(hd)

With the values of potentials and Rk in the each fluid domain, the hydrodynamic loads
are given by:
Z
Fk = gA Nk ds
(4.26)
s

where Nk is the generalized normal vector on the wetted surface of the barge.
Moreover, the hydrodynamic coefficients are calculated as follows:
Z
2
2
Rl Nk ds
akl + ibkl =

(4.27)

where akl is the added mass matrix and bkl is damping coefficients matrix.
The equations of motions of the barge are given as follows:
( 2 (M + a22 ) + ib22 )X2 ( 2 a24 + ib24 )X4 = F2
( 2 (M + a33 ) + ib33 gB)X3 = F3
r
8
2
2
( a42 + ib42 )X2 ( (I + a44 ) + i(b44 +
BQ4 ) C44 )X4 = F4
X4

(4.28)
(4.29)
(4.30)

where X2 , X3 and X4 are the sway, heave and roll motions. M and I are the mass and the
moment of inertia of the barge.
M = Bd

I = Bd(B/3)2

(4.31)

B is the beam of the barge and C44 is the restoring coefficient.


C44 = g(

d
B3
Bd( + zg ))
12
2

(4.32)

where zg is the center of gravity of the barge.


X 4 is the standard deviation of the roll velocity which can be solved through iterative
method for a given irregular sea state defined by its spectrum S(). The coefficient BQ4 equals
to Cd B 4 /2.
75

According to the experiments carried out within a Club pour Les Actions de Recherches
sur les Ouvrages en Mer (CLAROM) project on roll resonance of barge described in Molin
and Remy [2001], the drag coefficient Cd is taken around 0.1 - 0.15 for rectangular barges of
somewhat lower beam over draft ratios (B/d = 5), and with lower centers of gravity (at the
free surface level). Due to the lower draft and higher center of gravity, higher value of Cd
coefficient is taken to 0.2. F2 , F3 and F4 are the hydrodynamic loads.

4.3

The results of hydrodynamic coefficients and loads

A two dimensional rectangular barge, with a beam B=1.0 m and a draft d=0.108 m, floats
in waves. The water depth h=1.0 m. The origin of the coordinate system is set on the free
surface at rest. The center of the gravity locates on the middle of the barge with vertical
position on the free surface. The hydrodynamic loads and coefficients will be calculated in the
frequency domain with the angular frequencies from 0.2 radian per second to 20 radian per
second.
The diffraction forces and moment are presented in figure 4.2. The added masses and
moment of inertia produced by the sway, heave and roll motions are given in figure 4.3. In
figure 4.4, the damping terms can be found.
1800

10000

1600

9000
8000

1400

7000

Force (zaxis)

1000
800
600

5000
4000

2000

200
0

6000

3000

400

1000

8
10
12
14
Angular frequency (rad/s)

16

18

20

8
10
12
14
Angular frequency (rad/s)

800
700
600

500

Moment

Force (yaxis)

1200

400

300
200

100
0

8
10
12
14
Angular frequency (rad/s)

16

18

20

Figure 4.2: The hydrodynamic forces and moment.

76

16

18

20

45

14
a24

40

12

35

10

30

Added mass (kg m/m)

Added mass (kg m/m)

a22

25
20
15

6
4
2

10

8
10
12
14
Angular frequency (rad/s)

16

18

20

8
10
12
14
Angular frequency (rad/s)

16

18

20

550
a33

Added mass (kg m/m)

500

450

400

350

300

8
10
12
14
Angular frequency (rad/s)

14

16

18

20

20
a42

a44

12

19

18

Added inertia (kg m 2/m)

Added mass (kg m/m)

10
8
6
4
2

16

15
14

13

2
4

17

8
10
12
14
Angular frequency (rad/s)

16

18

12

20

8
10
12
14
Angular frequency (rad/s)

16

18

20

Figure 4.3: These figures present the added masses produced by the sway (a22), heave (a33) and roll (a44)
motions and their coupling terms (a24 and a42).

77

250

70
b22

b24
60

200
50

40

Damping

Damping

150

30

100

20
10

50
0

8
10
12
14
Angular frequency (rad/s)

16

18

10

20

8
10
12
14
Angular frequency (rad/s)

16

18

20

1600
b33
1400
1200

Damping

1000
800

600

400
200

8
10
12
14
Angular frequency (rad/s)

70

16

18

20

25
b42

b44

60
20
50

15

Damping

Damping

40
30

10

20

10
5
0

10

8
10
12
14
Angular frequency (rad/s)

16

18

20

8
10
12
14
Angular frequency (rad/s)

16

18

20

Figure 4.4: These figures present the damping produced by the sway (b22), heave (b33) and roll (b44) motions
and their coupling terms (b24 and b42).

78

4.4

The equations of motions in the time domain

The general equations of motion in the time domain can be written as follows:
6
X
j=1

{(Mij + mij )
xj +

Kij (t )x j ( )d + Cij xj } = Fi (t) i = 1, ..., 6

(4.33)

where :
Mij is the inertia matrix;
mij is the infinite frequency added mass matrix;
Kij is the retardation function;
Cij is the matrix of hydrostatic restoring stiffness;
Fi is the time varying exciting forces, which may include all kinds of external influences;
xj is the motions of the floating body.
The problem we consider here is a two dimensional case, so only the sway, heave and roll
motions are calculated. Next, all the terms in the equations of motions will be discussed.
1. The inertia matrix

m 0 0
Mij = 0 m 0
0 0 I
where m is the mass of the barge and I is the moment of inertia which is related to the center
of the gravity.
2. The hydrostatic restoring force coefficients

0 0
0
C = 0 C33 0
0 0 C44

The coefficient in the sway motion is zero and the heave motion C33 is taken to gB. The
B3
d
roll motion C44 is given by g(
Bd( + zg )) where zg is the z coordinate of the center of
12
2
gravity.
3. The external exciting forces
The time history of the forcing functions in regular waves can be generated easily based on
the absolute value of the forces calculated in the frequency domain.
Fi (t) = ramp(t) fi () cos(t)

(4.34)

where ramp(t) is the ramp function used for the beginning of the motions.
A quadratic damping term, nonlinear viscous roll damping contributions, is added in the
right-side of the roll equation of motion, as an external moment, under the form:
1

Cv4 = CdB 4 L||
2

(4.35)

where Cd is a drag coefficient, the same as the frequency domain. B is the beam of the barge
and L is the length along x-axis. is the roll velocity.

79

4. The retardation function


The time-domain equations can describe motions of any kind, also harmonic motions. Comparing the equations of motions in the frequency domain and time domain in a simple harmonic
motion, the retardation functions are related to the damping coefficients in the frequencydomain description.
Z
2
Kij ( ) =
bij () cos( )d
(4.36)
0
In practice, the integration is carried out in a limited frequencies range.
Z
2 2
Kij ( ) =
bij () cos( )d
(4.37)
1
The details of numerical computations
Z
2 2
bij () cos( )d
Kij ( ) =
1
Z
2 2
sin( )
=
bij ()(
) d
1

2
2
= (bij (2 ) sin(2 ) bij (1 ) sin(1 )) +

bij ()(

cos
) d

N
2
2 X
bn
= (bij (2 ) sin(2 ) bij (1 ) sin(1 )) + 2
(cos(n+1 ) cos(n ))

n=1
(4.38)

With
N = 2 1 ; bn : bij (n+1) bij (n ); : n+1 n .

For = 0, the values of the retardation functions should be calculated by following equation.
Z
2 2
Kij ( = 0) =
bij ()d
(4.39)
1
5. The damping term
The lower limit of integration in the damping term shows the equation of motion is laid
on the past of the motions of floating body (Memory effect). Due to the characteristics of the
retardation function, only the very near past T of the floating body motion is important. The
limits can be obtained using a numerical test with a series of different limits.
Z t
Z t
Kij (t )x j ( )d
Kij (t )x j ( )d
(4.40)
0

tT

6. The infinite frequency added mass


Based on a similar process of derivation as the retardation function, the frequency-independent
added masses are related to the added masses and damping coefficients in the frequency-domain
description.
Z
1
mij = aij () +
Kij (t) sin(t)dt
(4.41)
0
where is an arbitrarily chosen value of the frequencies.
80

Theoretically, the values of mij are independent with the frequencies . Due to the various
approximations involved, these computations did not yield exactly the same value. An effective
method is to calculate a series of values of frequency independent added mass with different
frequencies and calculate the mean value.
The details of numerical computations present as follows:
Z t1
Z t1
cos(t)
K(t) sin(t)dt =
K(t)(
) dt

0
0
Z t1
cos(t)
1

dt
K (t)
= (K(0) K(t1 ) cos(t1 )) +

0
Z
sin(t)
1 t1
1
K (t)(
) dt
= (K(0) K(t1 ) cos(t1 )) +

nt
1
1 X K
= (K(0) K(t1 ) cos(t1 )) + 2
(sin(tn+1 ) sin(tn )) (4.42)

n=1 t
where K = K(tn+1 ) K(tn ); t nt = t1 .

4.5

Comparisons between frequency domain and time


domain

30

1400

29

1380

28

1360

27

1340

Constant added masses

Constant added masses

A three dimensional barge in beam waves is considered here with the length of the barge
L=2.55 m. The diffraction loads, added masses and radiation damping, given by the former
model, are just multiplied by the length of the barge. Theoretically, the added masses mij
(figure 4.5) should be constant for all the frequencies. Here the values for a range of frequencies
are given and the mean values will be taken as the constant added masses.

26
25
24
23
22

1300
1280
1260
1240

21
20

1320

1220

8
9
Frequency (rad/s)

10

11

1200

12

40

8
9
Frequency (rad/s)

10

11

12

8
9
Frequency (rad/s)

10

11

12

39
6.5

37

Constant added masses

Constant added masses

38

36
35
34
33
32

5.5

4.5

31
30

8
9
Frequency (rad/s)

10

11

12

Figure 4.5: This figure presents the constant added mass coefficients of sway(up left), heave(up right) and
roll(down left) motions and coupling motions(down right) between sway and roll motions.

81

The RAOs of the sway, heave and roll motions are compared in the frequency domain and
the time domain in figures 4.6, 4.7 and 4.8.
4.5
time
frequency

4
3.5

RAO

3
2.5
2
1.5
1
0.5
0

4
Period (s)

Figure 4.6: The comparison of RAO in sway motion between the frequency domain and time domain.

1.4
time
frequency
1.2

RAO

0.8

0.6

0.4

0.2

4
Period (s)

Figure 4.7: The comparison of RAO in heave motion between the frequency domain and time domain.

82

4
time
frequency

3.5

RAO

2.5
2

1.5
1
0.5
0

4
Period (s)

Figure 4.8: The comparison of RAO in roll motion between the frequency domain and time domain.

The RAOs in frequency domain are calculated by the equations of motions given by equation
4.28, 4.29 and 4.30. For the time domain computations, the values of ordinates are obtained
from the amplitudes of sway, heave and roll motions. The time varying external forces are
calculated according to equation 4.34. The quadratic damping term is considered with drag
coefficient 0.2. In figure 4.6, the sway motions of the barge are compared and the results in
the time domain coincide well with those of the frequency domain. The heave motion is shown
in figure 4.7. With an increasing incident wave period, RAO of heave motion tends to one
which means the barge floats up and down with the same amplitude as incident waves. The
resonance in roll motion can be found in figure 4.8 clearly. When the excitation period is close
to the roll natural period, the amplitudes of roll motions are much higher. The numerical
results calculated in the time domain coincide well with those of the frequency domain.

83

84

Chapter 5
The comparisons of coupled motions
between numerical results and
experimental results
In this chapter, the motions of a floating body, the rectangular tank supported by the
barge, will be discussed for comparisons between experimental and numerical results. The
experiments were carried out in the BGO-FIRST wave tank, at la Seyne sur mer. The tank is
partially filled with liquid and the barge is floating in waves with constant water depth. Both
the frequency domain analysis and the time domain analysis are carried out for the coupled
motions of the floating body which is limited by the internal and external fluid. In the frequency
domain, the linear theory is adopted for the liquid in the tank and eigen-functions expansion
method is used for the external fluid. A quadratic damping term is added in the equation
of the roll motion. In the time domain, the Bousinesq-type approach is used for the internal
fluid. The hydrodynamic coefficients needed for the equations of motions in the time domain
are based on the eigen-functions expansions method. The quadratic damping term in roll is
also taken into account in the coupled equations of motions. Both the results in the frequency
domain and the time domain simulations are compared with the experimental results.

85

5.1

The description of the experiments

The experimental campaign was carried out at the BGO-FIRST (Bassin de Genie Oceanique
FIRST) which is located at la Seyne sur mer. BGO-FIRST is a multipurpose basin mainly
dedicated to maritime structures studies for the offshore oil industry and in coastal engineering. It is a basin of ocean engineering designed to perform tests in combined wave, current and
wind. It has a central well for offshore applications, a vertical moveable floor and a horizontal
mobile platform. The main features are:
length : 40 m; width : 16 m;
depth : 5 m
depth of central well : 10 m,
diameter : 5 m
wave period : from 0.6 to 4 s
wave height : max. 0.6 m
current velocity :
maximum 0.4 m/s for 3 m
of depth
maximum 1.2 m/s for 1.1
m of depth
maximum longitudinal velocity
of the platform : 1.2 m/s
maximum transverse velocity of
the platform : 0.8 m/s
maximum velocity of wind : 5
m/s.
The experimental model is shown as follows. Two rectangular tanks are fixed on the middle
of a rectangular barge whose bow and stern are slanted. The mass of the barge model with
the empty glass-tanks was 169 kg, with a center of gravity at 0.237 m above keel level and a
gyration radius of 0.414 m in roll. With a total height of water of 0.58 m in the two tanks (0.29
m + 0.29 m, or 0.19 m + 0.39 m), no additional weight was required to achieve the target
draft. At lower filling level (0.19 m + 0.19 m), additional masses (40 kg) were put on the deck,
with their centers of gravity approximately 0.28 m above keel level and 0.35 m away from the
longitudinal axis (transverse direction).
The model was installed in the basin in the middle of the testing section, with its longitudinal axis parallel to the wave-maker line (only beam waves were considered). It was anchored
to the walls of the tank, with 4 aerial lines. The stiffness of the mooring system resulted in a
natural period in sway of 20 seconds (model scale).

86

The parameters of the barge


length : 3 m; width : 1 m;
depth : 0.267 m
Mass without tanks : 127 kg;
Mass of tanks : 42 kg
Additional masses : 116 kg
Draft : 0.108 m
Volume of displacement :
0.285 m3
The parameters of the tanks
Attached objects : 5 kg
Mass of empty tanks : 37 kg
Internal length : 0.8 m;
Internal width : 20.25 m;
Height : 0.6 m.

The motion of the barge was measured through the optical tracking system KRYPTON
RODYM DMM 6D. Its accuracy is very good. The reference point for the translation components was at the deck level (0.267 m above keel level). The water motion inside the tanks
was measured by five ORCA capacitive probes. Three were installed in tank 2, which always
had the higher filling level, at 25, 180 and 350 mm from the wall closet to the wave-maker.
Two were installed in tank 1, at 25 and 180 mm from the corresponding wall. The probes are
numbered from 1 to 5, with numbers 1, 2, 3 in tank 2, probe 1 being close to the middle and
probe 3 by the wall, and numbers 4 and 5 in tank 1, probe 4 being by the wall.

87

5.2

The frequency domain computations

The motions of the floating body (Figure 5.1) are affected by the internal fluid in the tank
and the external fluid which the barge floats. The forces and moments given by the internal
fluid are calculated based on the linear theories which have been described in chapter 2. The
forces and moments given by the external fluid are calculated based on the eigen-function
expansions method which has been described in chapter 4. In the framework of linear theory,
only the sway and roll motions are discussed in our two dimensional model.
According to the theories in chapter 2, the free surface elevation in the tank is written as:
X
(y, t) =
An (t)n (y)
(5.1)
n

where
n (y) =

n
n (y, 0)
g

(5.2)

The total velocity potential in the tank is given by


XZ t
X
z, t) =
(y,
( An ( )d )n n (y, z) +
X j (t)j
n

(5.3)

j=2,4

where j is the infinite frequency radiation potential. n (y, z) are the potentials of the natural
modes of oscillation in the tank. X 2 and X 4 are the velocity of the sway and roll motions.
The forces and moments produced by the internal liquid are given as follows:
F~ =

X
n

~
An (t)f~n Ma ()X

(5.4)

~
where An (t) are the amplitudes of sloshing modes and X(t)
= [y(t), (t)] are the motions of
~
the tank. The expressions of the coefficients fn = [fn2 , fn4 ] and infinite frequency added mass
matrix Ma () can be found in chapter 2.

liquid
free surface

barge

bottom
Figure 5.1: 2D sketch of a floating body in waves with constant water depth.

88

The forces and moments, together with the diffraction loads calculated by the eigen-function
expansions method, are taken in the sway and roll equations of motions of the floating body.
X
(M + m22 )
y + m24
+ b22 y + b24 = F2 +
fn2 An
(5.5)
n

X
(I44 + m44 )
+ m24 y + b44 + b24 y + BQ4 |
|
+ C44 = C4 +
[fn4 (zC zG )fn2 ]An (5.6)
n

Here M and I44 are the mass and the moment of inertia. m22 , m24 and m44 are the added
masses which has accounted the internal fluid Ma () and the external fluid. The sway and roll
motions are represented by y and . The bij terms are the outer radiation dampings and BQ4
the viscous damping in roll. C44 is the restoring coefficient. The F2 and C4 are the diffraction
loads. Here the term zC zG is due to the distance between the bottom of the tank and the
center of gravity of the floating body.
These equations are coupled with the equations of the sloshing modes
An + B1n A n + n2 An = Dn2 [
y
(zC zG )] + Dn4

(5.7)

The coefficient B1n is the linear damping term. The coefficients Dn2 and Dn4 are coupling
coefficients which are obtained from the n potentials.
In the frequency domain computations, the equations 5.5, 5.6 and 5.7 can be reformed as
follows:


(M + m22 ) 2 ib22 y
X


+ m24 2 ib24
fn2 An = F2 ()
(5.8)
n

"

(I44 + m44 ) 2 i b44 +

8
BQ4

+ C44

X


+ m24 2 ib24 y
[fn4 (zC zG )fn2 ]An = C4 ()

(5.9)

2Dn2 y + 2 [Dn4 (zC zG )Dn2 ]




+ 2 iB1n + n2 An = 0

(5.10)

where is the angular frequecy in waves and n are the natural angular frequencies of sloshing
modes. The last equation is repeated for each sloshing mode, for n = 1, N, where N is the
total number of modes.
For the damping term in the roll motion, the quadratic damping is introduced in the roll
equation of motion as follows:
1
Cv4 = Cd B 4 L||

(5.11)
2
If the stochastic linearization is applied, the damping moment is written as
r
2
Cd B 4 L
(5.12)
Cv4 =

where B is the beam of the barge, L its length along x-axis. Cd is a drag coefficient which is
taken to 0.2. is the standard deviation of the roll velocity. This procedure ensures identical
energy dissipations for Gaussian processes.
89

5.3

The time domain computations

As shown in the chapter 4, the equations of motions are given like this
Z t
X
{(Mij + mij )
xj +
Kij (t )x j ( )d + Cij xj } = Fi (t)

(5.13)

The coefficients in the equations of motions are the same as before, but the forces and moments
in the right side of the equations include both the internal and external effects. Following the
theory in chapter 3, the internal forces and moments are calculated by integrating the pressure
on the wetted surface s.
Z
Z
~
~
F = p ~n ds; M = p ~r ~n ds
(5.14)
s

where the pressure p is calculated by the Boussinesq model in the time domain.
The time history of the external forcing functions in irregular waves can be generated easily
when the forces in regular waves are known as a function of the wave frequency.
The general representation of the surface elevation corresponding to a particular energy
spectrum is given by:
X
(t) =
Ai cos(i t + i )
(5.15)
i

with i are random numbers. The amplitudes of the harmonic wave components are defined
by the spectral density of the waves:
A2i = 2S(i )i

(5.16)

where S(i ) is the spectral density.


Now the time history of the force are given by:
X
F (t) = R{
f (j )Aj ei(j t+j ) }

(5.17)

where f (j ) are the forces in the frequency domain with unit wave amplitude.

The equations of motions

Motions

The Boussinesq-type models

Forces
Figure 5.2: The iterative computational process in the time domain.

Initially, the free surface elevations and the velocity potential in the tank are given. The
motions of the liquid in the tank will be calculated by the Boussinesq model and the forces
and moments produced by the internal fluid in the tank can be output. The external forces
and moments are calculated by the equation 5.17. Both the internal and external forces and
moments are taken into the equations of motions of the floating body. And then the motions
of the floating body can be obtained, which are again used for the computations of the tank.
As shown in figure 5.2, the two-way coupling computations can be carried out with the time
stepping. And then the time histories of the free surface elevations in the tank and the motions
of the floating body can be obtained.
90

5.4

Validations of the coupled motions based on linearized free surface conditions

Firstly, the linear theory will be adopted in both the frequency domain and the time domain
for the validations of the coupled motions. The forces and moments produced by the liquid
in the tank are calculated by the linear theory in the time domain. According to the theories
in chapter 2, the pendulum equation for the sloshing motions in the tank still works when the
motions of the tank is not harmonic since the coefficients Dnj do not involve the frequency .
Some experiments of a floating body in a wave tank are present as follows.
The water depth of the external fluid is 1 meter and beam waves are considered here which
is parallel to the x-axis. The equivalent length of the barge along x-axis is 2.55 m. The draft
of the barge is 0.108 meter and the beam of the barge is 1 meter.
The length of the rectangular tank is 0.8 meter and the width is 0.5 meter. These values
refer to the internal dimensions. The filling level of the fluid in the tank is set to 0.19 meter.
The mass of the barge with empty glass tank was 169 kg, with a center of gravity at 24 cm
above keel level and a gyration radius of 41.4 cm in roll. For the case 19 cm filling level, the
additional mass is 40 kg and their centers of gravity are approximately 28 cm above keel line
and 35 cm away from the longitudinal axis.
The stiffness of the mooring system resulted in a natural period in sway of 20 seconds.
For the irregular waves, the produced wave spectra were of JONSWAP type, with = 2.
In figure 5.3, the spectral density of the irregular waves is given with the peak period Tp =1.6
second and the significant wave height Hs =0.066 meter. The signal of the incident waves used
in the time domain computations during 640 seconds is given in figure 5.4.
In the frequency domain computations, the response amplitude operator (RAO) of the sway
and roll motions should be calculated firstly. And then the spectral densities of the sway and
roll motions can be calculated based on the following equation:
p
RAO = Sout /Sin
(5.18)
where Sin is the spectral density of the incident wave and Sout is the spectral density of the
sway and roll motions.
The stochastic linearization method is applied in the computations of the standard deviation
in the equation of the roll motion. The iterative process is used for solving the standard
deviation (Molin [2002]). In jth iterative step, the equation of the roll motion should be firstly
solved with all the frequencies.
"
!
#
r
8
(j1)
(I44 + m44 ) 2 i b44 +

BQ4 + K44 (j)



X


+ m24 2 ib24 y
[fn4 (zC zG )fn2 ]An = C4 ()
(5.19)
n

Then the RAOs in the roll motion are integrated for the spectral density S().
(j)
( )2

(j) ()((j) ())S() 2d


1

The updated standard deviation will be used for the next iterative step.
91

(5.20)

In the time domain computations, the quadratic damping term is considered as an external
moments.
1
Cv4 = Cd B 4 L||

(5.21)
2
where the drag coefficient Cd is the same as that in the frequency domain.
The time histories of the sway and roll motions are transformed based on the FFT with
Welchs method (Welch [1967]). The spectral densities of the sway and roll motions calculated
in the time domain are compared to the results in the frequency domain.
Due to the free surface in the tank, the restoring coefficients of the floating body should be
decreased by:
g
(5.22)
bl3
12
where b = 0.5 m is the internal width of the tank and l = 0.8 m is the internal length.
The center of the bottom of the tank is not coincide with the center of gravity, and some
coefficients needed to be corrected are this distance z.
Dn4 Dn4 z Dn2
fn4 fn4 z fn2
Ma24 () Ma24 z Ma22 ()
Ma44 () Ma44 + z z Ma22 ()

(5.23)
(5.24)
(5.25)
(5.26)

There are three positions where free surface elevations are measured by capacitive probes in
the experiments. The probes 1, 2 and 3 are located at the positions 350 mm, 180 mm and 25
mm away from the left wall.
The results in the following figures are based on the JONSWAP type irregular waves with
the peak period Tp =1.6 second and significant wave height Hs =0.066 meter. The spectral
analysis in the time domain are based on the signals during 640 seconds.
In figures 5.5 and 5.7, the spectral densities of the sway and roll motions are given. The
peak angular frequency of the incident waves is 3.927 rad/s. At 19 cm filling level, the angular
frequency of the first sloshing mode, given by the linear dispersion relationship, is 4.95 rad/s.
It is close to the roll natural frequency 4.8 rad/s when the water in the tank is frozen. The
coupling effects can be found in these two figures around 6 rad/s. The figures 5.6 and 5.8 give
the time histories of the motions during 400 seconds.
The spectral densities of the relative free surface elevations at gauge 3, 2 and 1 are given
in figures 5.9, 5.11 and 5.13. The coupling effects can be found around the first sloshing mode.
Due to the linear theories in the tank, the second sloshing mode 8.2 rad/s can not be present in
these figures. But the third sloshing mode 10.6 rad/s appears. The position of the gauge 3 is
close to an anti-node in the third sloshing mode and the spectral density of the third sloshing
mode in gauge 3 is larger than in gauge 2 and gauge 1. The time histories of the relative free
surface elevations of the gauge 3, 2 and 1 are given in figures 5.10, 5.12 and 5.14.
The response amplitude operators (RAOs) of the roll motion and the relative free surface
elevations are shown. In figure 5.15, the RAO of the roll motion of the floating body in the
frequency domain is compared with the results in the time domain. The RAO in the time
domain is obtained from the cross-spectrum between the considered signal and the incident
waves. The RAOs at the positions of gauge 3, 2 and 1 are shown respectively in figure 5.16,
figure 5.17 and figure 5.18.

92

1.6

x 10

1.4

Spectral density

1.2
1
0.8

0.6
0.4
0.2
0

8
10
12
Frequency (rad/s)

14

16

18

20

elevations (m)

elevations (m)

elevations (m)

Figure 5.3: The peak period is Tp =1.6 second and the significant wave height is Hs =0.066 meter. The range
of the frequencies are from 0.25 times peak frequency to 2.5 times peak frequency.

0.1
0
0.1

20

40

60

80

100
time (s)

120

140

160

180

200

220

240

260

280

300
time (s)

320

340

360

380

400

0.1
0
0.1
200
0.1
0
0.1
400

450

500

550

600

time (s)

Figure 5.4: The incident wave adopted in the time domain with peak period Tp =1.6 s and significant wave
height Hs =0.066 m.

93

1.5

sway

x 10

Time
Frequency

Spectral density

0.5

5
6
7
Frequency (rad/s)

10

11

Figure 5.5: Comparison between frequency domain and time domain for the sway motion.

0.04
0.03

Amplitude (m)

0.02

0.01
0

0.01
0.02

0.03
0.04

20

40

60

80

100
Time (s)

120

140

160

180

200

220

240

260

280

300
Time (s)

320

340

360

380

400

0.04
0.03

Amplitude (m)

0.02

0.01

0
0.01
0.02

0.03
0.04
200

Figure 5.6: Time history of the sway motion during 400 seconds.

94

1.5

roll

x 10

Time
Frequency

Spectral density

0.5

5
6
7
Frequency (rad/s)

10

11

Figure 5.7: Comparison between frequency domain and time domain for the roll motions.

0.1
0.08
0.06

Amplitude (radian)

0.04
0.02
0
0.02
0.04
0.06
0.08
0.1

20

40

60

80

100
Time (s)

120

140

160

180

200

220

240

260

280

300
Time (s)

320

340

360

380

400

0.1
0.08
0.06

Amplitude (radian)

0.04
0.02
0
0.02
0.04
0.06
0.08
0.1
200

Figure 5.8: Time history of the roll motion during 400 seconds.

95

1.2

gauge 3

x 10

Time
Frequency
1

Spectral density

0.8

0.6

0.4

0.2

6
8
Frequency (rad/s)

10

12

14

Figure 5.9: Comparison between frequency domain and time domain for the capacitive probe 3. The probe is
located 25 mm away from the left wall.

0.1
0.08
0.06

Amplitude (m)

0.04
0.02
0
0.02
0.04
0.06
0.08
0.1

20

40

60

80

100
Time (s)

120

140

160

180

200

220

240

260

280

300
Time (s)

320

340

360

380

400

0.1
0.08
0.06

Amplitude (m)

0.04
0.02
0
0.02
0.04
0.06
0.08
0.1
200

Figure 5.10: Time history of the free surface elevation for probe 3 during 400 seconds.

96

gauge 2

x 10

Time
Frequency
5

Spectral density

6
8
Frequency (rad/s)

10

12

14

Figure 5.11: Comparison between frequency domain and time domain for the capacitive probe 2. The probe is
located 180 mm away from the left wall.

0.08
0.06

Amplitude (m)

0.04

0.02
0

0.02
0.04

0.06
0.08

20

40

60

80

100
Time (s)

120

140

160

180

200

220

240

260

280

300
Time (s)

320

340

360

380

400

0.08
0.06

Amplitude (m)

0.04

0.02

0
0.02
0.04

0.06
0.08
200

Figure 5.12: Time history of the free surface elevation for probe 2 during 400 seconds.

97

gauge 1

x 10

Time
Frequency
5

Spectral density

6
8
Frequency (rad/s)

10

12

14

Figure 5.13: Comparison between frequency domain and time domain for the capacitive probe 1. The probe is
located 350 mm away from the left wall.

0.02
0.015

Amplitude (m)

0.01

0.005
0

0.005
0.01

0.015
0.02

20

40

60

80

100
Time (s)

120

140

160

180

200

220

240

260

280

300
Time (s)

320

340

360

380

400

0.02
0.015

Amplitude (m)

0.01

0.005

0
0.005
0.01

0.015
0.02
200

Figure 5.14: Time history of the free surface elevation for probe 1 during 400 seconds.

98

3
Time
Frequency
2.5

RAO

1.5

0.5

5
6
frequency (rad/s)

Figure 5.15: The RAO of the roll motion of the floating body.

3.5
Time
Frequency
3

2.5

RAO

1.5

0.5

5
6
frequency (rad/s)

Figure 5.16: The RAO of the relative free surface elevation at gauge 3.

99

3
Time
Frequency
2.5

RAO

1.5

0.5

5
6
frequency (rad/s)

Figure 5.17: The RAO of the relative free surface elevation at gauge 2.

0.8
Time
Frequency

0.7

0.6

RAO

0.5
0.4

0.3
0.2
0.1

5
6
frequency (rad/s)

Figure 5.18: The RAO of the relative free surface elevation at gauge 1.

100

5.5

Comparisons between Boussinesq model and linear


theory

The Boussinesq model and the linear theory in the time domain are compared with the same
floating body but with a different incident wave amplitude. The peak period of the incident
wave is Tp =1.6 second and the significant wave height is Hs =0.017 m. The spectral density of
the incident wave is given in figure 5.19. In figure 5.20, the time history of the incident waves
during 640 seconds are shown in the time domain computations. The filling level of the tank
is also 19 cm and the water depth is 1 m.
The spectral densities of sway and roll motions are presented in figures 5.21 and 5.23 and
the time histories of these two motions are compared in figures 5.22 and 5.24. Due to the
nonlinear free surface conditions adopted in the Boussinesq model, the significant difference of
the spectral densities can be found in the second sloshing mode in gauge 3 which is located
close to the wall (figure 5.25). The spectral density around the first sloshing mode coincides
well. The reason is the anti-node of the second sloshing mode present at gauge 3. The situation
is similar for the gauge 1 which is located in the middle of the tank. The peak around the
second sloshing mode can be found in the Boussinesq model (figure 5.29). But the amplitude
around the first sloshing mode at gauge 1, a node for the first sloshing mode, is much smaller
than for the gauge 3.
The gauge 2 which is located around a quarter of the length of the tank, This location
corresponds to a node for the second sloshing mode and no significant peak of the spectral
density can be found in the Boussinesq model (figure 5.27).
The time histories of the free surface elevations at gauge 3, 2 and 1 are compared in figures
5.26, 5.28 and 5.30. These figures show the phases informations between the Boussinesq model
and the linear theory. The phases of the incident waves adopted in these two numerical models
are the same. Due to the acceleration of the tank used in the linear theory and the velocity of
the tank used in the Boussinesq model, the phase discrepancy (around half of the peak period)
may be observed for the 3 wave gauges.
The RAOs of the roll motion of the floating body are compared between the Boussinesq
model in the time domain and the linear theory in the frequency domain (figure 5.31). In figure
5.32, the RAO of the relative free surface elevation at gauge 3 is given and the difference for
the second sloshing mode is significant.
The difference between the figure 5.33 and figure 5.34 corresponds to 2 different measurement positions of the free surface. In figure 5.33, the RAO of the relative free surface elevations
at gauge 2 is given, which is located 0.18 m away from the left wall. Here the point is not
exactly the node for the second sloshing mode. The peak can be found around 8.2 rad/s corresponding to the small peak in figure 5.27. In figure 5.34, the measured point is located at
0.2 m away from the left wall. Due to the node of the second sloshing mode, there is no peak
around 8.2 rad/s. The RAO of the free surface elevation at gauge 1 is shown in figure 5.35.

101

1.2

x 10

Spectral density

0.8

0.6

0.4

0.2

8
10
12
Frequency (rad/s)

14

16

18

20

elevations (m)

elevations (m)

elevations (m)

Figure 5.19: The peak period is 1.6 second and the significant wave height is 0.017 m.

0.02
0
0.02

20

40

60

80

100
time (s)

120

140

160

180

200

220

240

260

280

300
time (s)

320

340

360

380

400

0.02
0
0.02
200
0.02
0
0.02
400

450

500

550

600

time (s)

Figure 5.20: The incident wave adopted in the time domain with peak period 1.6 s and significant wave height
0.017 m.

102

1.4

sway

x 10

linear
Boussinesq
1.2

Spectral density

0.8

0.6

0.4

0.2

5
6
Frequency (rad/s)

10

Figure 5.21: The spectral densities of the sway motions of the floating body.

0.01
linear
Boussinesq

amplitude (m)

0.005

0.005

0.01

0.015

10

20

30
time (s)

40

50

60

Figure 5.22: The time histories of the sway motions of the floating body during 60 seconds.

103

roll

x 10

linear
Boussinesq

Spectral density

5
6
Frequency (rad/s)

10

Figure 5.23: The spectral densities of the roll motions of the floating body.

0.04
linear
Boussinesq

0.03

amplitude (radian)

0.02
0.01
0

0.01
0.02
0.03

0.04

10

20

30
time (s)

40

50

60

Figure 5.24: The time histories of the roll motions of the floating body during 60 seconds.

104

1.5

gauge 3

x 10

linear
Boussinesq

Spectral density

0.5

5
6
Frequency (rad/s)

10

Figure 5.25: The spectral densities of the relative free surface elevations at gauge 3.

0.04
linear
Boussinesq
0.03

amplitude (m)

0.02

0.01

0.01

0.02

0.03

10

20

30
time (s)

40

50

60

Figure 5.26: The time histories of the relative free surface elevations at gauge 3 during 60 seconds.

105

gauge 2

x 10

linear
Boussinesq
6

Spectral density

5
6
Frequency (rad/s)

10

Figure 5.27: The spectral densities of the relative free surface elevations at gauge 2.

0.025
linear
Boussinesq

0.02
0.015

amplitude (m)

0.01
0.005
0
0.005
0.01
0.015
0.02
0.025

10

20

30
time (s)

40

50

60

Figure 5.28: The time histories of the relative free surface elevations at gauge 2 during 60 seconds.

106

gauge 1

x 10

linear
Boussinesq
5

Spectral density

5
6
Frequency (rad/s)

10

Figure 5.29: The spectral densities of the relative free surface elevations at gauge 1.

x 10

linear
Boussinesq

6
4

amplitude (m)

2
0
2
4
6
8
10
12

10

20

30
time (s)

40

50

60

Figure 5.30: The time histories of the relative free surface elevations at gauge 1 during 60 seconds.

107

7
Boussinesq
Frequency
6

RAO

6
7
frequency (rad/s)

10

Figure 5.31: The RAO of the roll motion of the floating body.

20
Boussinesq
Frequency

18
16
14

RAO

12
10
8
6
4
2
0

6
7
frequency (rad/s)

Figure 5.32: The RAO of the relative free surface elevation at gauge 3.

108

10

3.5
Boussinesq
Frequency
3

2.5

RAO

1.5

0.5

6
7
frequency (rad/s)

10

Figure 5.33: The RAO of the relative free surface elevation at gauge 2, i.e. 0.18 m from the left wall.

3.5
Boussinesq
Frequency
3

2.5

RAO

1.5

0.5

6
7
frequency (rad/s)

10

Figure 5.34: The RAO of the relative free surface elevation at 0.2 m away from the left wall.

109

18
Boussinesq
Frequency

16
14

RAO

12
10
8
6
4
2
0

6
7
frequency (rad/s)

10

Figure 5.35: The RAO of the relative free surface elevation at gauge 1.

5.6

Comparisons between numerical results and experimental results

A series of experiments have been carried out based on different peak periods and significant
wave heights. In the following comparisons, the peak period of the incident wave is 1.6 second
and the significant wave height is 0.034 m. The walls of the tanks are smooth and the filling
levels of both tanks are 19 cm. In figure 5.36, the comparison between the theoretical values
and the measured values of the JONSWAP incident waves is given.
Based on the signal of the incident waves measured in the wave tank, the values of the
phases can be obtained:
Z Tf
2
ij
(t)eij t dt
(5.27)
Aj e =
Tf 0
where t is the free surface elevation in the wave tank. j are the angular frequencies. Aj and
j are the amplitudes and phases. These phases will be taken into the irregular signals of the
external forces and moments of the floating body. The figure 5.37 presents the comparisons of
the free surface elevations between measured results and reconstructed solutions. The reconstructed free surface elevations are obtained using equation (5.15). The amplitudes and phases
come from equation 5.27. The discrepancy is due to limited range of frequencies.
Figure 5.38 5.40 present the time histories of the sway, heave and roll motions during 40
seconds. The amplitudes and phases are compared between the Boussinesq model (red line)
and the experimental results (blue line). The lower frequencies signals (mooring system) of
the sway motion has been filtered. Figure 5.41 and figure 5.42 give the time histories of the
free surface elevations at probe 2 which is located 180 mm away from the left wall and probe
3 which is located 350 mm away from the left wall.

110

4.5

Tp=1.6 second; Hs=0.034 m

x 10

JONSWAP
Experiment

4
3.5

spectral density

3
2.5
2
1.5
1
0.5
0

5
6
7
pulsation (rd/s)

10

11

Figure 5.36: The comparison of the spectral density between theoretical solution and experimental results.

0.02
reconstructed
measured

0.015

free surface elevation (m)

0.01
0.005
0

0.005
0.01
0.015

0.02

10

15

20
time (s)

25

30

35

40

Figure 5.37: The comparisons of free surface elevations between measured results and reconstructed solutions.

111

sway
0.015
Experiment
Boussinesq

Amplitude (m)

0.01

0.005

0.005

0.01

10

15

20
Time (s)

25

30

35

40

Figure 5.38: The comparison of the time histories of the sway motion between the Boussinesq method and the
experimental results.

heave
0.015

0.01

Amplitude (m)

0.005

0.005

0.01

0.015

10

15

20
Time (s)

25

30

35

40

Figure 5.39: The comparison of the time histories of the heave motion between the Boussinesq method and the
experimental results.

112

roll
0.04
Experiment
Boussinesq
0.03

Amplitude (rad)

0.02

0.01
0

0.01

0.02
0.03

0.04

10

15

20
Time (s)

25

30

35

40

Figure 5.40: The comparison of the time histories of the roll motion between the Boussinesq method and the
experimental results.

gauge2
0.03

0.02

Amplitude (m)

0.01

0.01

0.02

0.03

10

15

20

25

30

35

Time (s)

Figure 5.41: The comparison of the time histories of the free surface elevation at probe 2 between the Boussinesq
method and the experimental results.

113

gauge3
0.04
Experiment
Boussinesq

0.03

0.02

Amplitude (m)

0.01

0.01

0.02

0.03

0.04

10

15

20

25

30

35

Time (s)

Figure 5.42: The comparison of the time histories of the free surface elevation at probe 3 between the Boussinesq
method and the experimental results.

sway

x 10

Boussinesq
Linear
Experiment

3.5

spectral density

3
2.5
2
1.5
1
0.5
0

6
pulsation

10

Figure 5.43: The spectral densities of the sway motion of the barge.

114

11

heave

x 10

Boussinesq
Experiment
3.5

spectral density

2.5
2

1.5

1
0.5

6
pulsation

10

11

Figure 5.44: The spectral densities of the heave motion of the barge.

roll

x 10

Boussinesq
Linear
Experiment

spectral density

5
4

3
2

1
0

6
pulsation

10

Figure 5.45: The spectral densities of the roll motion of the barge.

115

11

4.5

gauge 3

x 10

Boussinesq
Linear
Experiment

4
3.5

spectral density

3
2.5
2
1.5
1
0.5
0

6
pulsation

10

11

Figure 5.46: The spectral densities of the relative free surface at gauge 3.

gauge 2

x 10

spectral density

Boussinesq
Linear
Experiment

6
pulsation

10

Figure 5.47: The spectral densities of the relative free surface at gauge 2.

116

11

gauge 1

x 10

Boussinesq
Linear
Experiment

spectral density

6
pulsation

10

11

Figure 5.48: The spectral densities of the relative free surface at gauge 1.

4.5
Boussinesq
linear
experiment

4
3.5

RAO

3
2.5
2
1.5
1
0.5
0

6
7
frequency (rad/s)

Figure 5.49: The RAO of the roll motion.

117

10

10
Boussinesq
linear
experiment

9
8
7

RAO

6
5
4
3
2
1
0

6
7
frequency (rad/s)

10

Figure 5.50: The RAO of the relative free surface at gauge 3.

4
Boussinesq
linear
experiment

3.5
3

RAO

2.5
2
1.5
1
0.5
0

6
7
frequency (rad/s)

Figure 5.51: The RAO of the relative free surface at gauge 2.

118

10

9
Boussinesq
linear
experiment

8
7

RAO

6
5
4
3
2
1
0

6
7
frequency (rad/s)

10

Figure 5.52: The RAO of the relative free surface at gauge 1.

Figure 5.43 5.45 present the sway, heave and roll motions of the floating body. The
spectral densities are calculated by the Fast Fourier Transformation (FFT) based on the signal
of the motions during 640 seconds. With the filling level 19 cm in the tanks, the angular
frequency of the first sloshing mode based on the linear dispersion equation is around 4.95
rd/s. And the roll natural frequency is around 4.8 rd/s (when the water in the tanks is frozen).
They are close to the peak frequency of the wave spectrum 3.93 rd/s. The coupled effect can
be found around the frequency 6.0 rd/s in the sway and roll motions. Due to the symmetric
shape of the rectangular tank, the second sloshing mode (8.2 rd/s) is not present in the sway
and roll motions.
In figure 5.46 5.48, the spectral densities of the relative free surface elevations are given.
For the first and second sloshing mode, the positions of the wave amplitude and the node can
be represented well based on the positions of the probes. In figure 5.46, the probe 3 is close
to the wall of the tank. From the figure, the peaks of the spectral density of the first and the
second sloshing modes can be found in the nonlinear Boussinesq method and the experiments.
According to the theories in chapter 2, the second sloshing mode can not be predicted by the
linear free surface conditions.
In figure 5.47, probe 2 is located a quarter of the length of the tank. It is a node for the
second sloshing mode but there are certain motions on the free surface due to the first sloshing
mode. So the values on the first sloshing mode is lower than the results of probe 3. And there
is no difference of the numerical results in the second sloshing mode between the linear free
surface conditions and the nonlinear free surface conditions. The coupled motions can also be
found around 6.0 rd/s.
In figure 5.48, the significant difference between the nonlinear and the linear theories can
be found in the second sloshing mode because the anti-node is located here. The nonlinear free
119

surface conditions in Boussinesq type approach coincides well with the experimental results.
Due to the position of a node, the amplitude of the first sloshing mode is smaller than the
values of probe 2 and 3.
The RAO of the roll motion of the floating body between Boussinesq model, linear theory
and experimental results is given in figure 5.49. The cross-spectrum between the considered
signal and the incident waves is used for the RAO of the roll motion. The figures 5.50, 5.51
and 5.52 show the RAOs of the relative free surface elevations at gauge 3, gauge 2 and gauge
1 respectively.

5.7

The stronger nonlinear case

The higher significant wave height Hs=0.066 m and the same peak period Tp=1.6 s is
used for the computations. In figure 5.53, the spectral density is compared between theoretical
value and experimental results. The parameters of the floating body are the same as the former
case. In figure 5.54, 5.55 and 5.56, the free surface elevations at gauge 1, 2 and 3 are shown
respectively. From the figures, the maximum amplitude of the free surface elevations close to
wall is almost 15 cm, with a peak greater than 25 cm.
The numerical results based on Boussinesq model are given in figures 5.57 and 5.58. As
shown in figures, the time domain simulations can not continued after 243 s. It corresponds to
very high value of the run-up on the wall. Due to the inappropriate damping terms included in
the Boussinesq model, the code becomes unstable and breaks. Furthermore, in our numerical
model, the derivative of the free surface profile on the wall is equal to zero. Actually, this condition can not work with strong nonlinearity. Maybe the problem can be solved by introducing
moving boundary algorithm. Future studies should be taken.
4

1.6

Tp = 1.6 s; Hs = 0.066 m

x 10

Theory
Experiment

1.4

spectral density

1.2
1
0.8

0.6
0.4
0.2

5
6
7
frequency (rad/s)

10

11

12

Figure 5.53: The spectral density of JONSWAP waves with peak period 1.6 s and significant wave height 0.066
m.

120

gauge 1
0.1
0.08

free surface elevation (m)

0.06
0.04

0.02
0

0.02
0.04

0.06

100

200

300
400
time (s)

500

600

700

Figure 5.54: The time history of the measured free surface elevation at gauge 1.

gauge 2
0.1
0.08

free surface elevation (m)

0.06
0.04
0.02
0
0.02
0.04
0.06
0.08

100

200

300
400
time (s)

500

600

700

Figure 5.55: The time history of the measured free surface elevation at gauge 2.

gauge 3
0.3
0.25

free surface elevation (m)

0.2
0.15

0.1
0.05

0
0.05

0.1

100

200

300
400
time (s)

500

600

Figure 5.56: The time history of the measured free surface elevation at gauge 3.

121

700

The free surface elevations on the wall


0.3
0.25

Amplitude (m)

0.2

0.15
0.1

0.05

0
0.05

0.1
200

205

210

215

220

225
230
Time (s)

235

240

245

250

Figure 5.57: A section of the time history of the free surface elevations at the left wall.

The free surface elevations on the wall


0.15

Amplitude (m)

0.1

0.05

0.05

0.1
200

205

210

215

220

225
230
Time (s)

235

240

245

250

Figure 5.58: A section of the time history of the free surface elevations at the right wall.

122

5.8

Conclusion

The coupled motions of the floating body are discussed based on wave tank experiments,
linear theory analysis and nonlinear Boussinesq model. The motions around the first sloshing
mode predicted by linear theory and Boussinesq model coincide well with experimental results.
Due to fully nonlinear free surface conditions, the free surface elevations calculated by Boussinesq model present a good agreement with experimental results on the second sloshing mode,
which is equal to zero for the linear theory. For violent incident waves, the Boussinesq model
can not work due to the limits of damping and boundary conditions on the walls.

123

124

Chapter 6
Concluding remarks and perspectives
The problem discussed in this thesis is the sloshing phenomenon. One part of the thesis
presents results for sloshing motions in a rectangular tank with forced motions. The other part
shows the coupled motions of a floating body, rectangular tank fixed on rectangular barge,
which is affected by the fluid in the tank and the external fluid around the barge. Both
linear theory and nonlinear numerical models are used for simulating sloshing motions. And
experiments were carried out in order to verify the numerical solutions.
In the framework of linearized free surface conditions, the free surface elevations in a rectangular tank can be considered simply as superposing a series of natural sloshing modes. The
sloshing motions can be related to the forced motions of tank through pendulum equation.
Although the derivation process is carried out in the frequency domain, the excitation motions are not limited to harmonic motions. From the time domain computations, the beating
phenomenon modulation of the envelope of the signal can be found in the beginning and the
transient gradually disappears due to the damping term in the equations of motions. The odd
natural sloshing modes can be solved by the linear modal approach. But the resonances of the
even natural sloshing modes can not be predicted due to its linearized free surface conditions.
Due to the limitations, the extended Boussinesq-type equations with fully nonlinear free
surface conditions are considered. The Boussinesq model is firstly validated by a flat bottom
tank. According to different filling levels in the tank, the regimes of nonlinear sloshing motions
can be divided roughly into shallow water theory, intermediate water depth and finite water
depth conditions. The nonlinear free surface elevations in shallow water are compared with
Antuono et al. [2012]. Firstly, a rectangular tank length l = 1.175m and filling level h = 0.06m
is considered. The free surface elevations on the wall with different excitation frequencies
are compared, which shows an overall good agreement. The second tank with length l =
0.6095m and filling level h = 0.06m is also considered. The ratio h/l = 0.098 is close to the
intermediate water depth regime. Two waves regime and single bump-wave-regime appear in
the tank with excitation frequencies /r = 0.97 and /r = 1.04 respectively. When the
excitation frequency is very close to natural frequency /r = 0.99, a nonlinear subharmonic
free surface elevation motion can be found. The results given by Boussinesq model present a
good agreement with results in the literature.
The finite water depth condition is validated by Wu et al. [1998] and Frandsen [2004]. The
2D rectangular tank with l/h = 2 is considered. From the comparisons of free surface profiles
at two instants, nonlinear free surface calculated by Boussinesq model and integral equation
method coincide well with the results in Wu et al. [1998]. Compared with linear solutions,
the steeper crests and flatter troughs can be found clearly. The time histories of free surface
125

elevations on the wall are compared with results of Frandsen [2004]. The results of Boussinesq
model present an overall good agreement with results of Frandsen [2004]. When the excitation
motions are more violent, the free surface elevations on the wall are higher. The slope of
tangent of the free surface profile close to walls is not equal to zero any more. The assumption
of boundary conditions on the wall is not suitable. So the comparisons between Boussinesq
model and results in Faltinsen et al. [2000] are not good. The moving boundary algorithm
should be applied.
For the intermediate water depth, the solutions of Boussinesq model are compared with
linear theory and integral equation method. The Boussinesq model is firstly validated through
small excitation amplitude (1 mm) sway motions. Due to small-amplitude external excitation,
the free surface elevations calculated by Boussinesq model give a good agreement with linear
theories in the off-resonance excitation periods. Larger excitation amplitude (2 mm) is also
considered for comparisons and the discrepancies between linear theory and nonlinear model
appear. The free surface profiles present steeper crests and flatter troughs with the Boussinesq
model due to its fully nonlinear free surface conditions. The fully nonlinear free surface conditions are also adopted in the integral equation method. As shown in the comparisons, the free
surface elevations calculated by Boussinesq model coincide well with integral equation method.
The resonance of second natural sloshing mode can be found with fully nonlinear free surface
conditions and asymmetric shape of tank. An asymmetric rectangular tank with an inclined
bottom is considered. The slope of the inclined bottom is small. A symmetric rectangular
tank is also used for comparisons with the same length of tank and same mass of internal
liquid. Compared with the solutions of symmetric rectangular tank, the natural periods of
asymmetric tank in the first and third sloshing modes are a little higher. The resonance in the
second sloshing mode appears in the solutions of asymmetric tank.
The numerical solutions of both flat and inclined bottom tank are compared with experimental results. Due to no damping term, the numerical amplitudes are a little higher than
experimental results. But the time histories of free surface elevations still present an overall
agreement with experimental results. A linear damping term is tested for the flat bottom
tank. The numerical amplitudes decrease and present an overall agreement with experimental
results. The forces calculated by Boussinesq model present an overall agreement with experimental results.
The sloshing motions in the tank are coupled with the motions of a barge. The diffraction
and radiation problem of the rectangular barge is solved by eigen-function expansions in the
framework of linearized free surface conditions. The comparisons of sway, heave and roll
motions in the frequency domain and time domain are given. Considering the length of barge
1.0 m and the draft 0.108 m, the drag coefficient in the quadratic damping term is taken equal
to 0.2.
The coupled motions of the floating body in the JONSWAP type beam waves is considered
in both numerical and experimental methods. The experiments were carried out at the BGOFIRST which is a multipurpose basin. Both the frequency domain and the time domain
computations are carried out for comparing experimental results. The linear theory is used in
the frequency domain and the Boussinesq approach is used in the time domain. The tank with
filling level 19 cm and internal length 0.8 m is considered. The natural frequencies at the first
sloshing mode and second sloshing mode are 4.95 rad/s and 8.2 rad/s respectively according
to linear dispersion relationship. The peak period of incident JONSWAP wave is 1.6 second.
The significant wave heights are 0.017 m, 0.034 m and 0.066 m respectively.
The comparison between Boussinesq model and linear theory is carried out with significant
126

wave height 0.017 m. Both the spectral densities and the time histories of motions and free
surface elevations are given. Due to the peak period close to the first sloshing mode and roll
natural period, the coupled effects can be found in the spectral densities of sway and roll
motions of floating body around the first sloshing mode. Three different positions on the free
surface of the tank are measured which are located at the probe 3 (close to wall), probe 2
(around a quarter of length) and probe 1 (close to middle of tank). The anti-node of the first
and second sloshing modes appear on the position of probe 3. As shown in the spectral density
of free surface elevations at probe 3, the coupled motions present around the first sloshing mode
coincide well for both linear theory and Boussinesq model. The resonance around the second
sloshing mode only appears in Boussinesq model due to its nonlinear free surface conditions.
The situation is different for the probe 2 corresponding to a node for the second sloshing mode.
There is no resonance in both Boussinesq model and linear theory. The anti-node of the second
sloshing mode is located at the position of probe 1. The significant difference can be found in
the spectral density around the second sloshing mode. The values of spectral density around
the first sloshing mode is much lower than in probe 3 because it is a node for the first sloshing
mode.
The experiment was carried out with significant wave height 0.034 m. The nonlinear free
surface conditions in Boussinesq model present an overall good agreement with experimental
results. The resonance of the second sloshing mode appear in both the Boussinesq model and
experimental results. For higher significant wave height 0.066 m, the Boussinesq model can
not work due to violent free surface. As shown in the experiments, the free surface elevations
in the tank has reached up to 25 cm. Damping and assumption of the Boundary condition
should be revised for stronger nonlinear simulations.

Perspectives
The slopes of the free surface profile on the walls are null, which is used in the thesis. For
the stronger nonlinear sloshing motions in the tank, the limitations are significant. A moving
boundary algorithm can be used to solve this problem. The free surface in the tank is extended
to the outside of the tank. In figure 6.1, the red domain is the extended fluid domain (dry
domain) and the black domain is the fluid domain in the tank (wet domain).
L


6= 0
y R


6= 0
y L

Figure 6.1: The sketch of the moving boundary algorithm.

A series of discretized points will be set in both the wet and dry free surface. The external
free surface is the linearly extended free surface profile which means that the red dashed line is
127

the tangent of the internal free surface. With the time stepping, the internal wet free surface
will be updated by the Boussinesq model and the external dry free surface will be extended
linearly. The discretized points close to the wet-dry boundary locate either in the wet domain
or the dry domain with the motions of the free surface. If the points locate in the wet domain,
the motions of the points will be calculated by the Boussinesq model. Otherwise, the positions
of the points will be determined by the linear extrapolation of the free surface displacement.

128

Bibliography
H. N. Abramson. Dynamic Behavior of Liquid in Moving Containers. Appl. Mech. Revs., 16:
501 506, 1963.
H. N. Abramson. The Dynamic Behavior of Liquids in Moving Containers, With Applications
to Space Vehicle Technology. Technical report, NASA, 1966.
Y. Agnon, P. A. Madsen, and H. A. Schaffer. A new approach to high order Boussinesq models.
J. Fluid Mech., 399:319 333, 1999.
M. Antuono, B. Bouscasse, A. Colagrossi, and C. Lugni. Two-dimensional modal method for
shallow-water sloshing in rectangular basins. J. Fluid Mech., 700:419 440, 2012.
H. F. Bauer, T.-M. Hsu, and J. T.-S. Wang. Liquid sloshing in elastic containers. Technical
report, NASA contractor report, 1976.
H. B. Bingham, P. A. Madsen, and D. R. Fuhrman. Velocity potential formulations of highly
accurate Boussinesq-type models. Coast. Eng., 56:467 478, 2009.
J. Boussinesq. Theorie des ondes et des remous qui se propagent le long dun canal rectangulaire horizontal, en communiquant au liquid contenu dans ce canal des vitesses sensiblement
pareilles de la surface au fond. J. Math. Pures Appl., 2nd Series 17:55 108, 1872.
C. A. Brebbia and J. Dominguez. Boundary elements: an introductory course. WIT Press,
Computational Mechanics Publications, 1992.
L. Brosset, Z. Mravak, M. L. Kaminski, S. Collins, and T. Finnigan. Overview of sloshel
project. In the 19th ISOPE, Osaka, Japan, 2009.
T. Bunnik and R. Huijsmans. Large-scale LNG sloshing model tests. International Journal of
Offshore and Polar Engineering, 19:814, 2009.
T. Bunnik and A. Veldman. MODELLING THE EFFECT OF SLOSHING ON SHIP MOTIONS. In OMAE, Shanghai, China, 2010.
B.-F. Chen and R. Nokes. Time-independent finite difference analysis of full non-linear and
viscous fluid sloshing in a rectangular tank. Journal of Computational Physics, 209:47 81,
2005.
G. F. Clauss, D. Testa, and F. Sprenger. Coupling effects between tank sloshing and motions
of a LNG carrier. In Proceeding of the 29th International Conference on Ocean, Offshore
and Arctic Engineering, Shanghai, China, 2010.
129

G. F. Clauss, S. A. Mavrakos, F. Sprenger, D. Testa, and M. Dudek. Hydrodynamic considerations for flng concepts. In Proceeding of the 30th International Conference on Ocean,
Offshore and Arctic Engineering, Rotterdam, The Netherlands, 2011.
G. F. Clauss, F. Sprenger, M. Dudek, and D. Testa. sloshing: From theory to offshore operations. In Proceeding of the 31th International Conference on Ocean, Offshore and Arctic
Engineering, Rio de Janeiro, Brazil, 2012.
R. Cointe, P. Geyer, B. King, B. Molin, and M.-P. Tramoni. Nonlinear and linear motions of
a rectangular barge in a perfect fluid. In Proc. 18th symposium naval hydrodynamics, Ann
Arbor, 1990.
J. W. Cooley and J. W. Tukey. An Algorithm for the Machine Calculation of Complex Fourier
Series. Math. Comp., 19:297301, 1956.
L. Delorme, A. Colagrossi, A. Souto-Iglesias, R. Zamora-Rodriguez, and E. Botia-Vera. A set
of canonical problems in sloshing, Part I: Pressure field in forced roll-comparison between
experimental results and SPH. Ocean Engineering, 36:168 178, 2009.
M. Ebrahimian, M. A. Noorian, and H. Haddadpour. A successive boundary element model for
investigation of sloshing frequencies in axisymmetric multi baffled containers. Engineering
Analysis with Boudnary Element, 37:383 392, 2013.
O. M. Faltinsen and A. N. Timokha. Sloshing. Cambridge University Press, 2009.
O. M. Faltinsen, O. F. Rognebakke, I. A. Lukovsky, and A. N. Timokha. Multidimensional
modal analysis of nonlinear sloshing in a rectangular tank with finit water depth. J. Fluid
Mech., 407:201 234, 2000.
P. Ferrant and D. Le Touze. Simulation of Sloshing Waves in a 3D Tank Based on a PseudoSpectral Method. In 16th International Workshop on Water Waves and Floating Bodies,
Hiroshima, Japan, 2001.
R. D. Firouz-Abadi, M. Ghamesi, and H. Haddadpour. A modal approach to second-order
analysis of sloshing using boundary element method. Ocean Engineering, 38:11 21, 2010.
B. Fornberg. A practical guide to pseudo spectral methods. 1995.
J. B. Frandsen. Sloshing motions in excited tanks. Journal of Computational Physics, 196:53
87, 2004.
Y. Fujino, L. Sun, B. Pacheco, and P. Chaiseri. Tuned Liquid Damper (TLD) for Suppressing
Horizontal Motion of Structures. J. Eng. Mech., 118:20172030, 1992.
A. Gedikli and M. E. Erg
uven. Evaluation of sloshing problem by variational boundary element
method. Engineering Analysis with Boundary Elements, 27:935943, 2003.
R. A. Gingold and J. J. Monaghan. Smoothed particle hydrodynamics: theory and application
to non-spherical stars. Mon. Not. R. Astron. Soc., 181:375 389, 1977.
M. F. Gobbi, J. T. Kirby, and G. Wei. A fully nonlinear Boussinesq model for surface waves.
Part 2. Extension to O(kh)4 . J. Fluid Mech., 405:181 210, 2000.
130

S. E. Hirdaris, N. J. White, N. Angoshtari, M. Johnson, Y. Lee, and N. Bakkers. Wave loads


and flexible fluid-structure interactions: current developments and future directions. Ships
and Offshore Structures, 5(4):307325, 2010.
C. W. Hirt and B. D. Nichols. Volume of Fluid (VOF) Method for the Dynamics of Free
Boundaries. Journal of Computational Physics, 39:201 225, 1981.
C. H. Hu and M. Kashiwagi. A CIP-based method for numerical simulations of violent free
surface flow. Journal of Marine Science and Technology, 9:143 157, 2004.
C. H. Hu, K.-K. Yang, and Y.-H. Kim. 3d numerical simulations of violent sloshing by cip-based
method. Journal of Hydrodynamics, 22:253258, 2010.
R. A. Ibrahim. Liquid Sloshing Dynamics: Theory and Applications. Cambridge University
Press, 2005.
E. Jamois, D. R. Fuhrman, Bingham H. B., and B. Molin. Wave-structure interactions and
nonlinear wave processes on the weather side of reflective structures. Coast. Eng., 53:929
945, 2006.
G.-G. Jes
us, A. G.-M. Carlos, Jose L. C., and Leo G. Two phase analysis of sloshing in a
rectangular container with Volume of Fluid (VOF) methods. Ocean Engineering, 73:208
212, 2013.
Y. Kim. A numerical study on sloshing flows coupled with ship motion - the anti-rolling tank
problem. Journal of Ship Research, 46:52 62, 2002.
Z. R. Kishev, C. H. Hu, and M. Kashiwagi. Numerical simulation of violent sloshing by a
cip-based method. Journal of Marine Science and Technology, 11:111 122, 2006.
J. F. Kuo, R. E. Sandstroem, T. W. Yung, M. N. Greer, and M. A. Danaczko. LNG tank
sloshing assessment methodology - The new generation. International Journal of Offshore
and Polar Engineering, 19:241 253, 2009.
Y. Liu. Effect of variable bathymetry on the linear and slow-drift wave responses of floating
bodies. PhD thesis, Ecole Centrale Marseille, 2010.
L. B. Lucy. A numerical approach to the testing of the fission hypothesis. J. Astron., 82:1013
1024, 1977.
P. A. Madsen and Y. Agnon. Accuracy and convergence of velocity formulations ofr water
waves in the framework of Boussinesq theory. J. Fluid Mech., 477:285 319, 2003.
P. A. Madsen and H. A. Schaffer. Higher order Boussinesq-type equations for surface gravity
waves-derivation and analysis. Phil. Trans. R. Soc. Lond., A 356:159, 1998.
P. A. Madsen, H. B. Bingham, and H. Liu. A new boussinesq method for fully nonlinear waves
from shallow to deep water. J. Fluid Mech., 462:1 30, 2002.
A. R. Mitchell and D. F. Griffiths. The Finite Difference Method in Partial Differential Equations. John Wiley & Sons, 1980.
131

N. N. Moiseyev. On the theory of nonlinear vibrations of aliquid of finite volume. Applied


Mathematics and Mechanics (PMM), 22(5), 1958.
B. Molin. Hydrodynamique des structures offshore. Editions Technip, Paris, 2002.
B. Molin and F. Remy. Model tests on a barge with two tanks. Technical report, 2001.
B. Molin, F. Remy, Rigaud S., and De Jouette C. LNG-FPSOs: frequency domain, cuopled
analysis of support and liquid cargo motions. In Proc. 10th IMAM Conf., Rethymnon, 2002.
B. Molin, F. Remy, A. Ledoux, and N. Ruiz. Effect of Roof Impact on Coupling Between Wave
Response and Sloshing in Tanks of LNG-Carriers. In Proceedings of the 27th International
Conference on Offshore Mechanics and Arctic Engineering, 2008.
J. J. Monaghan. smoothed particle hydrodynamics. Rep. Prog. Phys., 68:1703 1759, 2005.
J. N. Newman. Wave effects on vessels with internal tanks. In Proceeding of the 20th International Workshop on Water Waves and Floating Bodies, Longyearbyen, Svalbard, 2005.
D. H. Peregrine. Long waves on a beach. J. Fluid Mech., 27:815 827, 1967.
M. Peric, T. Zorn, O. Moctar, T. E. Schellin, and Y.-S. Kim. Simulation of sloshing in LNGtanks. Journal of Offshore Mechanics and Arctic Engineering, 131:03110111, 2009.
A. Rafiee, K. P. Thioagarajan, and J. J. Monaghan. SPH Simulation of 2D Sloshing Flow in
a Rectangular Tank. In ISOPE, Osaka, Japan, 2009.
A. Rafiee, S. Cummins, M. Rudman, and K. Thiagarajan. Comparative study on the accuracy
and stability of sph schemes in simulating energetic free-surface flows. European Journal of
Mechanics B/Fluids, 36:116, 2012.
J. N. Reddy. An Introduction to the Finite Element Method (Third ed.). McGraw-Hill, 2005.
W. J. Rider and D. B. Kothe. Reconstructing Volume Tracking. Journal of Computational
Physics, 141:112 152, 1998.
E. F. G. Van Daalen, J. Gerrits, G. E. Loots, and A. E. P. Veldman. Free surface anti-roll
tank simulations with a volume of fluid based navier-stokes solver. In Proceeding of the 15th
International Workshop on Water Waves and Floating Bodies, Ceaserea, Israel, 2000.
H. K. Versteeg and W. Malalasekera. An Introduction to Computational Fluid Dynamics: The
Finite Volume Method (second edition). Prentice Hall, 2007.
C. Z. Wang and B. C. Khoo. Finite element analysis of two-dimensional nonlinear sloshing
problem in random excitations. Ocean Engineering, 32:107 133, 2005.
P. D. Welch. The Use of Fast Fourier Transform for the Estimation of Power Spectra: A
Method Based on Time Averaging Over Short, Modified Periodograms. IEEE Trans. Audio
Electroacoustics, AU-15:70 73, 1967.
C.-H. Wu and B.-F. Chen. Sloshing waves and resonance modes of fluid in a 3d tank by a
time-independent finite difference method. Ocean Engineering, 36:500 510, 2009.
132

C.-H. Wu, O. M. Faltinsen, and B.-F. Chen. Numerical study of sloshing liquid in tanks with
baffles by time-independent finite difference and fictitious cell method. Computers & Fluids,
63:9 26, 2012.
G. X. Wu, Q. W. Ma, and R. E. Taylor. Numerical simulation of sloshing waves in a 3D tank
based on a finite element method. Applied Ocean Research, 20(6):337 355, 1998.
G. H. Xu, C. Xu, Q. I. En-rong, and X. K. Gu. Experimental Investigation of Sloshing Loads
and Structural Dynamic Responses. Ship Mechanics, 15, 2012.
T. Yabe and T. Aoki. A universal solver for hyperbolic equations by cubic-polynomial interpolation. I. one-dimensional solver. Computer Physics Communications, 66:219 232,
1991.
T. Yabe, T. Ishikawa, P. Y. Wang, T. Aoki, Y. Kadota, and F. Ikeda. A universal solver
for hyperbolic equations by cubic-polynomial interpolation. II. two- and three- dimensional
solvers. Computer Physics Communications, 66:219 232, 1991.
V. E. Zakharov. Stability of periodic waves of finite amplitude on the free surface of a deep
fluid. J. Appl. Mech. Tech. Phys., 9:190 194, 1968.

133

You might also like