Table of Contents Author Guidelines Submit a Manuscript
Journal of Applied Mathematics
Volume 2012, Article ID 686824, 15 pages
http://dx.doi.org/10.1155/2012/686824
Research Article

A Time-Domain Boundary Element Method for Wave Diffraction in a Two-Layer Fluid

State Key Laboratory of Coastal and Offshore Engineering, Dalian University of Technology, Liaoning, Dalian 116024, China

Received 21 January 2012; Revised 6 April 2012; Accepted 25 April 2012

Academic Editor: Ioannis K. Chatjigeorgiou

Copyright © 2012 Ying Gou et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Abstract

A time-domain numerical model is established based on the higher-order boundary element method (HOBEM) to simulate wave diffraction problem in a two-layer fluid of finite depth. There are two possible incident wave modes (surface-wave mode and internal-wave mode) exist in the incident wave for a prescribed frequency in a two-layer fluid. For surface-wave mode, the hydrodynamic characters of fluid particles are similar to single-layer fluid. For the internal-wave mode, through the definition of a new function respected to velocity potentials of upper and lower fluid on the interface by using matching condition, a single set of linear equations is set up to compute the time histories of wave forces and wave profiles by using a fourth-order Runge-Kutta method. An artificial damping layer is adopted both on the free surface and interface to avoid the wave reflection. Examinations of the accuracy of this time-domain algorithm are carried out for a truncated cylinder and a rectangular barge, and the results demonstrate the effectiveness of this method.

1. Introduction

In the classical view of ocean hydrodynamics, the fluid under consideration is usually assumed to be of constant density. In the real ocean, however, density of sea water changes due to the variations in temperature and salinity in the water depth direction. In the ocean with density stratification, we usually simplify this complex case as a two-layer fluid model. In this model there is a thin layer called the pycnocline, the density of fluid above and below this pycnocline is approximately constant. In this paper we also focus on this two-layer fluid.

The internal waves will be generated on the interface between the two fluid layers with different types such as incident solitary wave and periodic wave. For internal solitary wave, the wave length is much longer than the characteristic length of structure. So the internal solitary wave forces acting on structures can be simulated only by the Morison formula [1, 2]. For the other typical harmonic internal waves, the wave length is over a wide range. The diffraction/radiation theory should be adopted when the characteristic length of structure is relative large.

Lamb [3] and Landau and Lifshitz [4] discussed some of the types of linear wave motion in a two-layer fluid. It was shown that incident waves in a two-layer fluid can propagate with two different wave numbers for a given frequency, corresponding to the surface-wave mode and internal-wave mode, respectively. Yeung and Nguyen [5] derived the Green functions in a two-layer fluid of finite depth for 3D problems. Ten and Kashiwagi [6] and Kashiwagi et al. [7] studied a two-dimensional radiationdiffraction problem for a general body floating in a two-layer fluid of finite depth by means of a boundary integral-equation method. Linton and McIver [8] studied the horizontal cylinders interacted with waves in two-layer fluids by multipole expansion method. By using the same method, Sturova [9] solved the problem of wave motions of a heavy fluid driven by the oscillations of a cylinder (radiation problem) and the scattering of an internal wave by a fixed cylinder (diffraction problem). Das and Mandal [10] also used this method to study the problems of radiation of water waves by a submerged sphere in deep water as well as in finite depth water with an ice-cover. You et al. [11] calculated the radiation and diffraction of water waves by a floating circular cylinder in a two-layer fluid of finite depth by using analytical method, the wave exciting forces for a floating circular cylinder due to incident waves of both surface- and internal-wave modes were presented.

To authors knowledge, most of the studies of the interaction between internal waves and structures are based on analytic method or frequency-domain approach based on the use of wave Green function. Because of the limitation of analytic method, it couldnot be used in practical engineering. For the wave Green function approach whether in frequency-domain or in time-domain approach, the complex wave Green function in two-layer fluid should be calculated. Therefore, developing a time-domain approach easy to implement to study the wave interactions between internal waves and structures has important significance.

In this present work, a time-domain higher-order boundary element method (THOBEM) was developed for internal wave diffraction from a 3D body located in the upper layer fluid. Integral equations in the upper and lower fluid domains are derived by applying the Green’s second identity to simple Green function and velocity potential in each layer, respectively. By the construction of a function on the interface, a single set of linear equations is set up to compute the time histories of wave forces and wave profiles by using a fourth-order Runge-Kutta method. An artificial damping layer is adopted to dissipate the scattering waves on both free surface and interface. The internal wave force and moment on a floating body are calculated and compared with published frequency-domain solution and analytic method, and a relatively good agreement was found.

2. Numerical Simulations

A fixed Cartesian coordinate system 𝑂𝑥𝑦𝑧 is adopted with the origin 𝑂 located at the undisturbed free surface, where the 𝑧-axis is positive upwards and the body is located in the upper layer fluid, as shown in Figure 1. Denote the densities and depths of the upper and lower fluid layers by 𝜌1, 1 and 𝜌2, 2, respectively. The ratio between the upper and lower layer density can be defined as 𝛾=𝜌1/𝜌2 and the total water depth is=1+2.

686824.fig.001
Figure 1: The coordinate system.
2.1. Incident Wave Potential and Dispersion Relation

The expressions of incident wave velocity potential of upper and lower domain in a two-layer fluid had been listed by some researchers [5, 6], and the detail of the present derivative process is listed below.

The fluid in each layer is assumed to be inviscid and incompressible, and the flow irrotational. The velocity potentials Φ(1)(𝑥,𝑦,𝑧,𝑡) and Φ(2)(𝑥,𝑦,𝑧,𝑡) in the fluid domain Ω1 and Ω2 satisfy Laplace equation: 2Φ(𝑚)=0,(2.1) where 𝑚=1,2 denote the upper and lower fluid domains. The linearized boundary conditions for free surface and interface are described as: Φ(1)𝑡𝑡+𝑔Φ𝑧(1)Φ=0,𝑧=0,𝑧(1)=Φ𝑧(2),𝑧=1,𝛾Φ(1)𝑡𝑡+𝑔Φ𝑧(1)=Φ(2)𝑡𝑡+𝑔Φ𝑧(2),𝑧=1.(2.2) A rigid sea bottom satisfies the no-penetration condition in the lower domain is given as: Φ𝑧(2)=0,𝑧=12.(2.3) The velocity potential is defined as the following equation: Φ(𝑚)𝜙(𝑥,𝑦,𝑧,𝑡)=Re(𝑚)(𝑥,𝑦,𝑧)𝑒𝑖𝜔𝑡.(2.4) Inserting (2.4) into boundary conditions of (2.2)~(2.3): 𝜕𝜙(1)(𝑥,𝑦,𝑧)=𝜔𝜕𝑧2𝑔𝜙(1)(𝑥,𝑧)=𝑘0𝜙(1)(𝑥,𝑦,𝑧),𝑧=0,(2.5)𝜕𝜙(1)(𝑥,𝑧)=𝜕𝑧𝜕𝜙(2)(𝑥,𝑧)𝜕𝑧,𝑧=1,𝛾(2.6)𝜕𝜙(1)(𝑥,𝑦,𝑧)𝜕𝑧𝑘0𝜙(1)=(𝑥,𝑦,𝑧)𝜕𝜙(2)(𝑥,𝑦,𝑧)𝜕𝑧𝑘0𝜙(2)(𝑥,𝑦,𝑧),𝑧=1,(2.7)𝜕𝜙𝑧(2)(𝑥,𝑧)𝜕𝑧=0,𝑧=12,(2.8) where 𝑘0=𝜔2/𝑔. The solutions of (2.5), (2.6), and (2.8) can be written as the following forms: 𝜙(1)(𝑥,𝑦,𝑧)=𝑍1(𝑧)𝑒𝑖𝑘(𝑥cos𝛽+𝑦sin𝛽),𝑍1(𝑧)=𝑎1cosh𝑘(𝑧+)cosh𝑘+𝑏1sinh𝑘𝑧sinh𝑘1𝜙(2)(𝑥,𝑦,𝑧)=𝑍2(𝑧)𝑒𝑖𝑘(𝑥cos𝛽+𝑦sin𝛽),𝑍2(𝑧)=𝑎2cosh𝑘(𝑧+),cosh𝑘(2.9) where 𝑘 denotes the wave number, 𝑎1, 𝑏1, and 𝑎2 are the unknown coefficients. Inserting the expression of 𝜙(𝑚)(𝑥,𝑦,𝑧) into the free surface boundary conditions (2.5), (2.6), and (2.8), we can get 𝑍1(𝑧)=𝑎1cosh𝑘(𝑧+)+𝑘cosh𝑘0𝑘𝑍tanh𝑘sinh𝑘𝑧2(𝑧)=𝑎1𝑘1+0𝑘tanh𝑘cosh𝑘1cosh𝑘sinh𝑘2cosh𝑘(𝑧+).cosh𝑘(2.10) The wave elevations 𝜂(1) and 𝜂(2) have the relations as following: Φ𝑡(1)=𝑔𝜂(1)Φ,𝑧=0,(2.11)𝑧(2)=𝜂𝑡(2),𝑧=1.(2.12) Assuming that 𝜂(2)=𝐴𝑒𝑖(𝑘𝑥cos𝛽+𝑘𝑦sin𝛽𝜔𝑡), the coefficients 𝑎1 can be obtained based on (2.12): 𝑎1=𝑖𝑤𝐴cosh𝑘𝑘sinh𝑘2𝑘1+0/𝑘tanh𝑘cosh𝑘1cosh𝑘/sinh𝑘2,(2.13) we can then obtain the incident wave potentials and wave elevations as follows: Φ𝑖(1)=𝑘𝜔𝐴0/𝑘sinh𝑘𝑧+cosh𝑘𝑧𝑘0cosh𝑘1𝑘sinh𝑘1Φsin(𝑘𝑥cos𝛽+𝑦sin𝛽𝜔𝑡)intheupperdomain𝑖(2)=𝜔𝐴cosh𝑘(𝑧+)𝑘sinh𝑘2𝜂sin(𝑘𝑥cos𝛽+𝑦sin𝛽𝜔𝑡)inthelowerdomain𝑖(1)=𝐴cosh𝑘11𝑘/𝑘0tanh𝑘1𝜂cos(𝑘(𝑥cos𝛽+𝑦sin𝛽)𝜔𝑡).𝑖(2)=𝐴cos(𝑘(𝑥cos𝛽+𝑦sin𝛽)𝜔𝑡),(2.14) where 𝐴 denotes the incident wave amplitude on the interface.

By inserting the expressions of 𝜙(𝑚)(𝑥,𝑦,𝑧) into (2.7), the dispersion equation in two-layer fluid is obtained as (𝛾+1)tanh𝑘1tanh𝑘2𝑘2+𝑘0tanh𝑘1tanh𝑘2𝑘0𝑘20(𝛾+1)tanh𝑘1tanh𝑘2=0.(2.15) It can be seen that three are two possible roots of (2.15), that is, traveling waves with a prescribed frequency 𝜔 in a two-layer fluid can propagate with two possible different wavenumbers 𝑘1 and 𝑘2. In order to find the physical significance of 𝑘1 and 𝑘2 more explicitly, on the assumption that the depth of each layer fluid goes to infinite, (2.15) can be simplified to the following version: 𝑘𝑘02𝑘(1𝛾)2𝑘0+(𝛾+1)=0.(2.16) The two roots of (2.16) are 𝑘𝑘0𝑘=1,𝑘0=(1+𝛾).(1𝛾)(2.17) It can be seen obviously that the former root is corresponding to the surface wave mode, and the other is corresponding to the internal wave mode. In this paper we only concern the internal wave mode.

2.2. Diffraction Problem

Under the assumptions of linear water wave theory, the total velocity potential in upper and lower fluid domain can be divided into a known incident potential and an unknown diffractive potential, while the wave elevation on the free surface and interface can be divided into a known incident wave elevation and an unknown diffractive wave elevation as the following: Φ(1)=Φ𝑖(1)+Φ𝑑(1),Φ(2)=Φ𝑖(2)+Φ𝑑(2),𝜂(1)=𝜂𝑖(1)+𝜂𝑑(1),𝜂(2)=𝜂𝑖(2)+𝜂𝑑(2).(2.18)

For the diffraction potential, it must satisfy the following governing equations and the boundary conditions: 2Φ𝑑(𝑚)=0,𝑚=1,2𝜕Φ(1)𝑑𝑧𝜕𝑛=𝜕Φ(1)𝑖𝑧𝜕𝑛on𝑆𝐵Φ(2)𝑑𝑧=0on𝑧=12,𝜂(2.19)(1)𝑑𝑡=Φ(1)𝑑𝑧Φ(a)on𝑧=0,(1)𝑑𝑡=𝑔𝜂𝑑(1)𝜂(b)(2.20)(2)𝑑𝑡=Φ(2)𝑑𝑧Φ(c)(1)𝑑𝑧=Φ(2)𝑑𝑧(d)on𝑧=1𝛾Φ(1)𝑑𝑡+𝑔𝜂𝑑(2)=Φ(2)𝑑𝑡+𝑔𝜂𝑑(2)(e).(2.21) For simplicity, we define a function 𝜑=𝛾Φ𝑑(1)Φ𝑑(2) on 𝑧=1, then the boundary condition (e) of (2.21) can be rewritten as: 𝜑𝑡=(1𝛾)𝑔𝜂𝑑(2)on𝑧=1.(2.22) The velocity potential on interface is discontinuous, but the construction of this function 𝜑 makes it easier to solve this boundary value problem. The details are introduced in the next subsection.

In order to avoid the reflection of scatter waves on the outer computational boundary, we introduce an artificial damping layer to absorb the scattered wave energy. On the outer part of the free surface and interface, a damping term is added to the free surface boundary conditions 𝜕𝜂𝑑(𝑚)=𝜕𝑡𝜕𝜙𝑑(𝑚)𝜕𝑧𝜈(𝑟)𝜂𝑑(𝑚),𝑚=1,2,𝜕Φ𝑑(𝑚)𝜕𝑡=𝑔𝜂𝑑(𝑚)𝜈(𝑟)Φ𝑑(𝑚),𝑚=1,𝜕𝜑𝜕𝑡=(1𝛾)𝑔𝜂𝑑(𝑚)𝜈(𝑟)𝜑,𝑚=2,(2.23) where 𝑚=1,2, and 𝜈(𝑟) denotes the damping coefficient is given by 𝑣(𝑟)=𝐶𝜔𝑟𝑟0𝛽𝜆2𝑟0𝑟𝑟00+𝛽𝜆𝑟<𝑟0,(2.24) where 𝐶 is the damping coefficient, 𝛽 is the beach breadth coefficient, and 𝜆 is the characteristic wave length. 𝑟0 and 𝑟0+𝛽𝜆 represent the inner and outer radius of the damping layer, respectively. In this paper values of 𝛼 and 𝛽 are equal to 1.0.

A ramping function is utilized when the simulation started. Its aim is to make the numerical simulation more stable and reach steady state faster. In this present simulation, a ramp function is applied to the first two waves by the following forms: 𝑅𝑚1(𝑡)=21cos𝜋𝑡2𝑇𝑡2𝑇1𝑡2𝑇,(2.25) where 𝑇 is the wave period.

2.3. Integral Equation and Its Solution by Time-Domain Method

Based on Green’s theorem, the boundary integral equations in the upper and lower layer fluid can be obtained to solve the above boundary value problems. The integral equation in the upper fluid includes the integration over the body surface, the free surface, and the interface:𝛼𝜙1𝑑𝑆𝐵𝜙1𝑑𝜕𝐺1𝜕𝑛𝑑𝑠+𝑆𝐹𝐺1𝜕𝜙1𝑑𝜕𝑛𝑑𝑠+𝑆𝐼𝐺1𝜕𝜙1𝑑𝜕𝑛𝑑𝑠=𝑆𝐵𝐺1𝜕𝜙1𝑑𝜕𝑛𝑑𝑠+𝑆𝐹𝜙1𝑑𝜕𝐺1𝜕𝑛𝑑𝑠+𝑆𝐼𝜙1𝑑𝜕𝐺1𝜕𝑛𝑑𝑠,(2.26) where Green’s function 𝐺1=1/4𝜋𝑟 is the Rankine source. In our case the body located in the upper layer is considered, thus the integral equation in the lower fluid only includes the integration over the interface. The integral equation in the lower layer domain has the forms as 𝑆𝐼𝐺2𝜕𝜙2𝑑𝜕𝑛𝑑𝑠=𝑆𝐼𝜙2𝑑𝜕𝐺2𝜕𝑛𝑑𝑠𝛼𝜙2𝑑,(2.27) where 𝐺2=1/4𝜋𝑟1/4𝜋𝑟2 includes its images about the seabed.𝑟is the distance between the field and the source point, and 𝑟2 is the distance between the field point and the image of the source point about the sea bed. 𝛼 is the so-called solid angle coefficient and satisfies the following forms: 0𝛼=1(sourcepointwithinthedomain)(sourcepointoutsidethedomain)1solidangle4𝜋(sourcepointontheboundary).(2.28) In this paper the direct method is used to calculate the solid angle [12].

Compared with constant panel method the higher-order boundary element method [13] is more accurate and efficient. In this numerical method, quadratic isoparametric elements are used to discretize the surface over which the integral is performed. Curved quadrilateral and triangular elements are placed by eight and six nodes, respectively.

Through the discretization of (2.26) and (2.27), two sets of linear equations can be obtained as: 𝑎11𝑎12𝑎13𝑎21𝑎22𝑎23𝑎31𝑎32𝑎33𝜕Φ𝑑(1)𝜕𝑛𝑆𝐹Φ𝑑(1)𝑆𝐵𝜕Φ𝑑(1)𝜕𝑛𝑆𝐼=𝑠11𝑠12𝑠13𝑠21𝑠22𝑠23𝑠31𝑠32𝑠33Φ𝑑(1)𝑆𝐹𝜕Φ𝑑(1)𝜕𝑛𝑆𝐵Φ𝑑(1)𝑆𝐼[𝑏](2.29)𝜕Φ𝑑(2)𝜕𝑛𝑆𝐼=[𝑡]Φ𝑑(2)𝑆𝐼.(2.30) It is worth noting that the normal velocities of the upper and lower domains on the interface are in the opposite direction. It means that the normal derivatives of the potentials at the interface have the following relation: 𝜕Φ𝑑(1)𝜕𝑛=𝜕Φ𝑑(2).𝜕𝑛(2.31) Applying the above relation to (2.30) and taking advantage of the function 𝜑 constructed in (2.22), we can combine (2.29) and (2.30) to get a single set of linear equations as follows: 𝑎11𝑎12𝑎13+1𝛾𝑠13[𝑡]1[𝑏]𝑎21𝑎22𝑎23+1𝛾𝑠23[𝑡]1[𝑏]𝑎31𝑎32𝑎33+1𝛾𝑠33[𝑡]1[𝑏]𝜕Φ𝑑(1)𝜕𝑛𝑆𝐹Φ𝑑(1)𝑆𝐵𝜕Φ𝑑(1)𝜕𝑛𝑆𝐼=𝑠11𝑠121𝛾𝑠13𝑠21𝑠221𝛾𝑠23𝑠31𝑠321𝛾𝑠33Φ𝑑(1)𝑆𝐹𝜕Φ𝑑(1)𝜕𝑛𝑆𝐵{𝜑}𝑆𝐼.(2.32) All the coefficient matrices [𝑎],[𝑠],[𝑡], and [𝑏] are constant at all time steps due to the time-independent integral boundary.

A time-stepping integration method is employed to obtain the values of wave elevations, the velocity potentials over the body surface, and the derivative of interface. The initial conditions as Φ𝑑(1)||𝑆𝐹=0,𝜂𝑑(1)Φ=0,𝑑(2)||𝑆𝐼=0,𝜂𝑑(2)𝜑||=0,𝑆𝐼=0.(2.33) After solving the boundary value problem and obtaining the fluid velocity on the free surface and interface at each time step, the general forms of the dynamic and kinematic boundary conditions (2.23) are used as ordinary differential equations to be marched in time and can be rewritten as 𝛿𝜂𝑚𝑑Φ𝛿𝑡=𝑓𝑚𝑑,𝜂𝑚𝑑,𝑡𝑚=1,2𝛿Φ1𝑑Φ𝛿𝑡=𝑔1𝑑,𝜂1𝑑,𝑡𝛿𝜑2𝑑𝜑𝛿𝑡=2𝑑,𝜂2𝑑.,𝑡(2.34) A 4th-order Runge-Kutta method is adopted for integrating (2.34) with time.

When the unknown velocity potentials over body surface are obtained, the pressure can be derived from the Bernoulli equation, the forces and moments acting on the body can be calculated by integrating the pressure over the mean body surface.

3. Numerical Validation of the Time-Domain Method

3.1. A Truncated Cylinder in a Two-Layer Fluid of Finite Depth

The water depths of the upper and the lower layers are 1/=0.7, 2/=0.3. The densities of the fluids in the upper and the lower layers are 𝜌1=998.2kg/m3 and 𝜌2=1027.2kg/m3. That means 𝛾=0.97. The cylinder has a radius of 𝑎/=0.5, and a draft of 𝑇/=0.5. The truncated cylinder is located in the upper fluid, as shown in Figure 2. In the following calculations, the surge and heave exciting forces have been nondimensionalized by 𝜌1𝑔𝑎𝐴, while the dimensionless factor for pitch exciting moment about 𝑦-axis is 𝜌1𝑔𝑎2𝐴.

686824.fig.002
Figure 2: Definition of sketch for a truncated cylinder.

The structural meshes are used to discretize the body surface, the free surface, and the interface. If the mesh size or the time step is too large, the convergence will be difficult to achieve. Conversely, it will take too much computation time if a large number of meshes are used. The numerical convergence tests are desired to ensure the mesh size, and the time step is suitable. Based on the numerical convergence tests, the radius of the free surface should be two to three times of the wave length. The mesh size on the free surface and interface should be about ten percent of the wave length, and the sizes of elements should be uniformly distributed. The time step for time marching is taken as Δ𝑡=𝑇/100 (𝑇 is the wave period).

The time histories of the wave forces and the wave moments on the cylinder due to internal wave mode at the dimensionless wave number 𝑘=4.0 are plotted in Figures 3, 4, and 5. It can be seen that very steady time-history results can be obtained by this time-domain method.

686824.fig.003
Figure 3: Time histories of horizontal wave force.
686824.fig.004
Figure 4: Time histories of vertical wave force.
686824.fig.005
Figure 5: Time histories of pitch moment.

To examine the accuracy of the present model, the calculation of horizontal wave force, vertical wave force and 𝑦-axis wave moment is carried out by the present time-domain method and the analytic solutions of results of You et al. [11] are also shown to compare. The results are shown in Figures 6, 7, and 8. The dimensionless wave frequency varied from 0.1 to 0.4.

686824.fig.006
Figure 6: Horizontal wave force on cylinder.
686824.fig.007
Figure 7: Vertical wave force on cylinder.
686824.fig.008
Figure 8: Pitch moment on cylinder.

From these figures, it is can be seen that the present results are in close agreement with the results obtained by the analytical solutions. The magnitudes of the horizontal and vertical wave forces decrease as dimensionless wave frequency decreases within the calculated range. The 𝑦-axis wave moment is the combination of two components: due to horizontal wave force and vertical force, so the trend is different.

3.2. A Rectangular Box in a Two-Layer Fluid of Finite Depth

To further validate the present method, the wave diffraction of a rectangular box is considered. In order to make a comparison with the numerical results of Nguyen and Yeung [14] by the frequency-domain method, the same box is studied which the length 𝐿, width 𝐵, and draft 𝐷 are 90 m, 90 m, and 40 m, respectively. The box is also located in the upper fluid with no forward speed. The water depths of the upper and the lower layers are 1=48m and 2=16m. The dimensionless factor for wave forces are 𝜌1𝑔𝐿𝐷, and for moment is 𝜌1𝑔𝐿2𝐷. The results are shown in Figures 9, 10, and 11.

686824.fig.009
Figure 9: Horizontal wave force on box.
686824.fig.0010
Figure 10: Vertical wave force on box.
686824.fig.0011
Figure 11: Pitch moment on box.

Three types of density ratio 𝛾=0.9, 0.7, and 0.1 are adopted. From the comparisons, we can find that the present results by time-domain method also have a good agreement with the frequency-domain method. From these figures, it is clear to see that the density ratio has obvious effects on the hydrodynamic forces of internal-wave mode. As the density ratio decrease, for a given frequencies, the magnitudes of the horizontal and vertical forces increase. As we have mentioned before, the trends of the magnitude of the wave moment is different.

4. Conclusions

A time-domain higher-order boundary element method is developed and it can be used to simulate the wave diffraction problem in a two-layer fluid. Through the construction of a function on the interface, the two integral equations in the upper and lower layers are combined, and a single set of linear equations are set up. The wave diffraction problem are solved with a fourth-order Runge-Kutta method for time marching.

In order to verify the validity of the present approach, a truncated cylinder and a rectangular box are taken as examples. Compared with the analytic solution and the numerical results by frequency-domain approach, the computation in the present time-domain approach agrees well with them. So the present time-domain approach can be used to study the wave diffraction problem in two-layer fluid and it can be extend to simulate the interaction between internal wave and floating structure.

In the present study, we find that the density ratio has obvious effects on the hydrodynamic forces of internal-wave mode. As the density ratio decrease, for a given frequencies, the magnitudes of the horizontal and vertical forces increase. The 𝑦-axis wave moment is the combination of two components: horizontal wave force and vertical force, so the trend is different.

Acknowledgments

This work is financially supported by the National Research Program of China (973-Program, Grant no. 2011CB013703). This work is also partially supported by Projects of Natural Science Foundation of China (Grant no. 51009023, 11072052, 50921001, 51109022).

References

  1. Y. L. Cheng, J. C. Li, and Y. F. Liu, “The induced flow field by internal solitary wave and it action on cylindrical piles in the stratified ocean,” in Proceedings of the 4th International Conference on Fluid Mechanics, pp. 296–299, 2004.
  2. Z. J. Song, B. Teng, Y. Gou et al., “Comparisons of internal solitary wave and surface wave actions on marine structures and their responses,” Applied Ocean Research, vol. 33, no. 2, pp. 120–129, 2011. View at Google Scholar
  3. H. Lamb, Hydrodynamics, Cambridge University Press, Cambridge, UK, 1993. View at Zentralblatt MATH
  4. L. D. Landau and F. M. Lifshitz, Fluid Mechanics, Pergamon Press, Oxford, UK, 1959, (translation into Japanese by H. Takeuchi, published in 1970 by Tokyo-Tosho).
  5. R. W. Yeung and T. Nguyen, “Radiation and diffraction of waves in a two-layer fluid,” in Proceedings of the 22nd Symposium on Naval Hydrodynamics, pp. 875–891, Washington, DC, USA, 1999.
  6. I. Ten and M. Kashiwagi, “Hydrodynamics of a body floating in a two-layer fluid of finite depth. Part 1. Radiation problem,” Journal of Marine Science and Technology, vol. 9, pp. 127–141, 2004. View at Google Scholar
  7. M. Kashiwagi, I. Ten, and M. Yasunaga, “Hydrodynamics of a body floating in a two-layer fluid of finite depth. Part 2. Diffraction problem,” Journal of Marine Science and Technology, vol. 11, pp. 150–164, 2006. View at Google Scholar
  8. C. M. Linton and M. McIver, “The interaction of waves with horizontal cylinders in two-layer fluids,” Journal of Fluid Mechanics, vol. 304, pp. 213–229, 1995. View at Publisher · View at Google Scholar · View at Zentralblatt MATH
  9. I. V. Sturova, “Problems of radiation and diffraction for a circular cylinder in a stratified fluid,” Fluid Dynamics, vol. 34, pp. 512–533, 1999. View at Google Scholar
  10. D. Das and B. N. Mandal, “Water wave radiation by a sphere submerged in water with an ice-cover,” Archive of Applied Mechanics, vol. 78, no. 8, pp. 649–661, 2008. View at Publisher · View at Google Scholar
  11. Y. X. You, Q. Shi, G. Wei, and G. P. Miao, “The interaction of water waves with a vertically floating cylinder in a two-layer fluid,” Journal of Shanghai Jiao Tong University, vol. 41, no. 9, pp. 1460–1464, 2007. View at Google Scholar
  12. B. Teng, G. Ying, N. Dezhi et al., “A higher-order BEM for wave-current action on structures direct computation of free-term coefficient and CPV integrals,” China Ocean Engineering, vol. 20, no. 3, pp. 395–410, 2006. View at Google Scholar
  13. B. Teng and R. Eatock Taylor, “New higher-order boundary element methods for wave diffraction/radiation,” Applied Ocean Research, vol. 17, pp. 71–77, 1995. View at Google Scholar
  14. T. C. Nguyen and R. W. Yeung, “Unsteady three-dimensional sources for a two-layer fluid of finite depth and their applications,” Journal of Engineering Mathematics, vol. 70, no. 1–3, pp. 67–91, 2011. View at Publisher · View at Google Scholar