XV Konferencja Naukowa - Korbielów 2003
 Metody Komputerowe w Projektowaniu i Analizie Konstrukcji Hydrotechnicznych
Standing waves of finite amplitude in water
of non-uniform depth
Kazimierz Szmidt1
1. INTRODUCTION
In analysing water structure interaction problems we have to solve boundary conditions
on a wetted surface of a structure which in general is a moving surface. The simplest case in
this field seems to be a rigid breakwater, not moving in space, with a wetted surface
changing in time. An example of the letter case is a breakwater constructed in the form of a
rigid vertical wall immersed in fluid of constant depth loaded with water waves. In many
practical cases we also meet cases with sloping bottom in the vicinity of the breakwater
structure, and thus, in analysis of the latter a change of water depth should be taken into
account. Frequently, in real conditions, one may observe that water waves approaching the
breakwater have a relatively big height, often changing in time. In describing the latter cases
a finite amplitude of the waves must be also taken into account. In describing the
aforementioned problems we usually apply the space variables and time (Eulerian
description). With these variables, the fundamental equation describing the irrotational
motion of an inviscid incompressible fluid assumes the form of the Laplace equation for a
velocity potential. The main difficulty of a solution of the equation is to solve the initial and
boundary conditions of the problem considered.
2
Z
piston generator
E D
H
a a 1
x(t)
H
~
b2
C
b1
H
2
a�
A
B
1
Z
L L
1 2
Fig. 1. Definition sketch.
1
Prof. dr hab. inż., Institute of Hydro-Engineering of PAS, ul. Kościerska 7, 80-953 Gdańsk
1
There are also regions of our interest however where a more preferable is to use the material
(Lagrangian) description, which allows us to simplify the solution of the boundary
conditions on a moving boundary of the fluid domain. In particular, the Lagrangian
variables are more convenient in calculating the hydrodynamic interaction of the fluid and
structure. In this paper we confine our attention to the material description of the initial
value problem of fluid being at rest and starting to move at a certain moment of time.
Following the problems mentioned above we will consider an approximate solution to the
problem of water waves of finite amplitude reflected from a rigid vertical boundary.
2. FORMULATION OF THE PROBLEM
The theoretical model considered is shown schematically in Fig.1. The waves are
generated by a piston type generator forming the left boundary of the domain. Due to
reflection from the right boundary (the vertical wall CD in the figure) we have the case of
standing waves in the fluid domain. The frequency of the generation was chosen to be equal
to the resonance frequency of the fluid domain. In this way a standing wave of finite
amplitude
growing in time on the right boundary is obtained. The case considered corresponds to
smooth waves sufficiently far off breaking. In what follows we confine our attention to the
plane problem of a fluid motion in the Euclidean space. Fundamental relations of the
problem considered may be found in Szmidt (2001). To make the discussion clear some of
the important results obtained in the work are summarized below. In order to describe the
fluid motion we introduce the Cartesian coordinate system in the actual
(zr , r =�1,2)
configuration. In the reference configuration, the Cartesian coordinates, corresponding to
l�
names of the fluid particles, are denoted by . Moreover, it is convenient to
(Z , l� =� 1,2)
introduce a common Cartesian coordinate system. The motion of the fluid is described by
the mapping of the names into the positions occupied by the material points at the time
t ł� 0
a� l� l�
zi (Z ,t) =� d�i Z +� wi (Z ,t) , (1)
l�
where d�i is the Kronecker s delta and are components of the displacement vector.
wi
l�
The Jacobian of the transformation is the determinant of the matrix of the transformation
gradient
J =� det[�z,ia�]�
, (2)
a�
,a�
where the symbol denotes the partial derivative with respect to . In a similar way,
Z
,t
,i
the symbol denotes the partial derivative with respect to , and the subscript means
zi
the partial derivative with respect to time. For the discussed, two  dimensional case, the
inverse of the matrix of the transformation gradient reads
��z,2 -� z12 ł�
1 2 ,
[�Z,a�]�=�
ę� ś� , (3)
i
J
ę�-� z,2 z1 ś�
1 ,1
�� ��
2
a�
where the subscripts mean the partial derivatives with respect to the material variables
Z
for a� =� 1,2 .
Knowing the above relations we may transform important formulae from the Eulerian
variables into the Lagrangian variables and vice versa. Thus, let us consider the potential
function expressed in terms of the material variables. In these variables the
F�(Zl� ,t)
Laplace equation assumes the following form
d�rs[�F�,b�Z,b�]� Z,a� =� 0
. (4)
s r
,a�
With respect to the potential, the velocity components are
l� b�
&�
wr (Z ,t) =� F�,b�Zr ,
(5)
where the dot denotes the material time derivative.
In the discussed case the fluid density is constant and the Jacobian of the mapping
l�
is equal to one. Accordingly, in what follows it is convenient to introduce
zi =� zi (Z ,t)
the  pressure function
p 1
l�
&�
P(Z ,t) =� =� -�h -� F� +� d�rsF�,b� Z,b�F�,g� Z,g� +� C(t).
r s
r� 2
(6)
C(t)
where is a  constant of the solution, and is the potential of the mass force due to
h
the gravitational field. In the spatial description the potential is given by the relation
h(zr ) =� gz3 , (7)
g
where the coordinates are chosen in such a way that acts vertically upwards and is
z3
the gravitational acceleration. For the considered two-dimensional problem one may
introduce the classical notation for the actual configuration and
z1 =� x, z2 =� y
2
for the reference configuration. With respect to this notation Eq (5)
Z1 =� X , Z =� Y
gives
&� &�
w1 =� u =� F�, X (1 +� v,Y ) -� F�,Y v, X ,
(8)
&� &�
w2 =� v =� -�F�, X u,Y +� F�,Y (1+� u,X ).
Having the velocity we may calculate the displacement components
t
l� a� l�
&�
u(Z ,t) =�
��u(Z , x�)dx� +� u(Z , t =� 0),
0
(9)
t
l� a� l�
&�
v(Z ,t) =�
��v(Z , x�)dx� +� v(Z , t =� 0) .
0
3
In order to describe the initial and boundary conditions, let us consider the case shown in
Fig. 1. The motion of the fluid is induced by the piston-type wavemaker (the rigid wall AE
in the figure) starting to move at a certain moment of time. For the case shown in Fig 1 the
boundary conditions are:
1 1
a) u(Z =� 0, t) =� x0 (t), b) u(Z =� L,t) =� 0,
,
2 2
c) v(Z =� 0, t) =� 0, d) P(Z =� H, t) =� const.,
(10)
x0 (t)
where describes the horizontal displacement of the wall AE in the figure, and the
constant in Eq (10d) will be assumed equal to zero. In the further part we will consider the
initial value problem of the fluid motion starting from rest.
3. GENERATOR MOTION
Let us consider the piston type generator starting to move at a certain moment of time.
It is assumed a smooth beginning of the fluid motion for which not only the velocity, but
also the acceleration field disappear at the initial moment of time. The motion of the
generator is assumed in the following form (Wilde and Wilde, 2001)
x0 (t) =� a �� A(t�) cosw�t +� D(t�)sin w�t
, (11)
w�
where is the angular frequency, and
a =� s-�3
1
A(t�) =� t3 exp(-�t�),
3!
ć�1+� 1 1
D(t�) =� 1 -� t� +� t�2 +� t�3 ��exp(-�t�), t� =� h�t,
�� ��
2! 3!
Ł� ł�
(12)
h� is a parameter.
where
One may check that the displacement together with its first and second derivatives are equal
to zero at the starting point. Moreover, with passing time, the generation goes
asymptotically to harmonic displacement with constant amplitude. In the further discussion
we confine our attention to generation described by the latest formulae. The non-linear
problem at hand has no closed analytical solution and therefore, in order to find a solution
of it, we have to approximate the fundamental equations by ones which are more tractable.
One of the methods of approximation is the infinitesimal - wave approximation based on a
perturbation procedure (Wehausen and Laitone, 1960) in which functions entering the
problem are expanded into power series in a small parameter.
4. SMALL PARAMETER REPRESENTATION OF THE FUNDAMENTAL
RELATIONS
Let us consider the equation of motion together with appropriate boundary and initial
e�
conditions. The functions entering the equations are expanded in a small parameter .
4
Thus, we have
l�
F�(Z , t) =� e�f�1 +� e�2f�2 +� e�3f�3 +� �� ����,
(13)
l� 2 2 3 l�
u(Z ,t) =� e�u1 +� e� u +� e�3u +� �� ����, v(Z ,t) =� e�v1 +� e�2v2 +� e�3v3 +� ���� �� .
where , and for
f�i
ui vi i =�1,2,3,������ are  components of the solutions.
In what follows we confine our attention to the second order expansion. Substituting these
e�
expressions into Eq (4) and collecting terms with the same powers in , one finds
e� �� f�111 +� f�122 =� 0,
, ,
2
e� �� f�2 +� f�2 +�
,11 ,22
+� 2[�f�111v1 -� f�1 (�u,1 +� v1 )�+� f�122u1 ]� =� 0,
, ,2 ,12 2 ,1 , ,1
(14)
where the terms up to the second order are displayed.
It is seen the linear component of the expansion results in the Laplace equation for the
l�
velocity potential identical to that one in the space variables, while for the
f�1(Z ,t)
higher component we have the Poisson s equation. Therefore, in constructing the first order
l�
solution for the velocity potential
we do not need to distinguish the co-ordinate
f�1(Z ,t)
systems. In a similar way, the expansion of the velocity components reads
l�
&�
u(Z ,t) =� e� �� f�11 +�
,
2
+� e� ��[�f�2 +� f�11v,1 -� f�12v1 ]� ,
,1 , 2 , ,1
(15)
l�
&�
v(Z ,t) =� e� �� f�12 +�
,
2
+� e� ��[�f�22 +� f�12u1 -� f�11u12]�.
, , ,1 , ,
Knowing that
l� l� 3 l�
h(Z ,t) =� gz3(Z ,t) =� g[�Z +� v(Z ,t)]�, (16)
we may assume C =� gH in Eq (6) and write
1
l� 3 l� l�
&�
P(Z ,t) =� g(H -� Z ) -� gv(Z ,t) -� F�(Z ,t) +� d�rsF�,b� Z,b�F�,g�Z,g�
.
r s
2
(17)
The first term on the right hand side of the equation means the hydrostatic pressure
3
. The pressure function may be thus rewritten as
P0 =� g(H -� Z )
, (18)
P =� P0 +� e�P1 +� e�2P2
5
where
1
&� &�
P1 =� -�gv1 -� f�1, P2 =� -�gv2 -� f�2 +� [�(f�1 )2 +� (f�12 )2 ]�.
,1 ,
2
(19)
With respect to the expansion (18) we can write the sequence of the dynamic boundary
conditions on the upper boundary. From the first of Eq (19) it follows
&�
gv1(Z1, H,t) +� f�1(Z1, H,t) =� 0 . (20)
Calculating the partial time derivative of the equation one obtains
&�&�1 +� gf�12 =� 0 .
f�
(21)
, 2
Z =�H
A similar procedure for the square term (the second equation in the relations 19) gives
&�&�2 2 , 1 , 1 , &� &�
e�2 �� f� +� gf�,2 +� g(�f�12u,1 -� f�11u,2)�-�(�f�11f�11 +� f�12f�12)� =� 0.
, , , 2
Z =�H
(22)
2
The boundary condition on the bottom A-B ( ) leads to the result
Z =� 0
&�
v1 Z2 =� 0 =� 0, �� f�12 Z =� 0 =� 0,
, 2
&�2 Z =�0
v =� 0, �� f�22 +� f�12u1 -� f�1 u12 Z =�0 =� 0.
2 , , ,1 ,1 , 2
(23)
On the slope B-C the normal component of the velocity equals zero, and thus
-� f�1 sin a� +� f�12 cos a� =� 0,
,1 ,
BC
1
1
-� f�,2 sin a� +� f�22 cos a� -� f�11
1 , , ,x�h�
��f� dt =� 0,
cos a�
t
BC
(24)
x�,h�
where the co-ordinate axes ( ) are defined by the line B-C and the normal vector of it.
On the vertical wall we have the condition, that the normal component of the
Z1 =� 0
&�
x0 (t)
velocity field is equal to the velocity of the generator. From the first of Eq. (15) it
follows
2 1 1
&�
f�11 Z1=�0 =� x0 (t), f�,1 +� f�11v,2 -� f�12v,1 Z1 =�0 =� 0.
, , ,
(25)
6
In a similar way, on the right boundary we have
Z1 =� L
2 1 1
f�11 Z1=�L =� 0, f�,1 +� f�11v,2 -� f�12v,1 Z1=�L =� 0.
, , ,
(26)
5. FD FORMULATION OF THE PROBLEM
In order to find a solution of the problem discussed we resort to discrete formulation of
it by means of the finite difference method. With this method a finite number of nodal
points of an assumed net is considered. For the domain of fluid shown in Fig.1 a non-
homogeneous spacing of points is assumed in which the spacing of vertical lines is equal to
a
while the spacing of horizontal lines is equal to b1 in the lower part, and b2 in the upper
part of the fluid domain, respectively. The assumed discrete model is shown in Fig. 1. The
differential equations for the components of the velocity potential are substituted by the
finite difference equations at all nodal points of the assumed net. To save the place we shall
a�
f�
use the notation for the first order velocity potential . For a typical nodal
f�1(Z ,t)
(i, j)
point within the lower part the fluid domain, where i means the number of a vertical
j
line and denotes the number of a horizontal line, the finite difference representation of
the Laplace equation reads
-� e�1f�i-�1, j -� f�i, j-�1 +� K1f�i, j -� f�i, j+�1 -� e�1f�i+�1, j =� 0
,
(27)
where
2
e�1 =� (�b1 a)� and K1 =� 2(1+� e�1) (28)
(i, j)
For a typical point within the upper part of the domain the relevant equations are
similar, but now, instead of b1 one should introduce b2 , and at the same time, e�2 and
K2 , respectively. In order to write the equations for points where the vertical dimension of
the mesh are changed the Taylor expansion has been applied.
The equations (27) are written for all nodal points of an assumed net, including
boundary points. At the same time, it should be noted that the boundary conditions of the
problem at hand involve not only the values of the potential but also the first and second
time derivatives as well as the derivatives with respect to the material co-ordinates of the
function. In order to perform integration of the equations in the time domain we resort to di-
screte description of the time, i.e. instead of the continuous time we introduce a sequence of
time steps with the increment . In the discrete time domain in constructing an ap-
D�t >� 0
proximate solution to the boundary value problem it is convenient to use the Wilson me-
q�
thod which enables us to establish algebraic equations of the problem mentioned on a com-
mon level of the discrete time. The method is based on a linear approximation of an acce-
leration vector at every point of the discrete time (Bathe, 1982). In the discussed problem of
7
the potential motion of the fluid, in place of the acceleration we are dealing with the second
time derivatives of the potential function.
6. THE FIRST ORDER SOLUTION
The first order solution of the problem mentioned is similar to a linear solution of the
potential flow in the Eulerian variables. In both cases we have to solve the classical Laplace
equation within the fluid domain together with appropriate initial and boundary conditions.
In constructing the solution we do not need to distinguish the co-ordinate system. Thus, we
have the following system of equations:
- the Laplace equation
2
ś�2f�1 ś� f�1
�2f�1 =� +� =� 0
, (29)
2 2
ś�x ś�y
- the boundary conditions
ś�f�1 ś�f�1 ś�f�1
=� x0(t), =� 0, &�&�1 +� g =� 0,
&� f�
ś�x ś�n ś�y
x=�0 A-� B -�C -� D D-� E
(30)
and initial conditions at that the displacement, velocity and acceleration field in the
t =� 0+�
fluid domain equal zero. In order to find a set of solutions of the problem at a sequence of
time steps we apply the formulae derived in the previous section.
7. THE SECOND ORDER SOLUTION
Having the first order solution at a certain moment of time (at time t =� const and
, where ) we can construct the second order solution of the
t =� const +� DT DT =� q�D�t
problem considered. To do this we have to solve the Poisson s equation together with
proper boundary and initial conditions. As compared to the first order solution the problem
becomes more complicated. Let us consider now the second order velocity potential which
should satisfy the Poisson s equation (the second of Eq.14)
�2f�2 +� RA =� 0 , (31)
where depends solely on the first order solution
RA
ć� ��
RA =� -�4��f�1XX 1 dt +� f�1XY 1XYdt
, ,XX , ,
��f� ��f� �� . (32)
�� ��
Ł� t t ł�
In the finite difference formulation, the last quantity must be calculated at each point of the
assumed net. From the last formula it is seen that in the numerical solution we have to
calculate approximate values of the second order derivatives of the first order potential with
respect to the material co-ordinates. A similar difficulty emerges in describing the
aforementioned free surface boundary condition on the upper boundary of the fluid domain.
The procedure described above enables us to find the two components of the solution for
the potential function at each level of time.
8
8. CONCLUDING REMARKS
Following the procedure mentioned above some numerical computations have been
performed. The frequency of the generation was chosen to be equal to the resonance
frequency of the fluid domain (the length of the generated wave was equal twice the
distance ED in Fig.1.). Some of the results of computations are illustrated in Fig.2, where
the free surface elevation at the right boundary (the vertical displacement of the point D in
the figure) and the resultant of pressure forces acting on the boundary, and its moment
relative to point C in the figure, are depicted in the subsequent plots. Since we are dealing
with generation of resonance frequency, the calculated quantities grow in time. From the
plots it is seen that the linear component of the solution behaves like a linear oscillator in
the resonance range i.e. it grows in time according to the formula . The frequency
t sin w�t
of the second component of the solution equals twice of the first component. The main
feature of the second component is that it grows in time in an exponential manner, and thus,
the solution obtained is valid only in a limited range of time. One should remember that in
constructing the approximate solution it has been assumed that the second component of the
solution should be sufficiently small as compared to the first order one. The results obtained
allow us to examine main features of the discrete formulation of the non-linear problem
considered.
9
E D
a d
d
=� 1 , =� 5 , de =� 1 , D�t =� 0.0541s
C
b c
c
A B
a b
5
4
Generator motion.
3
2
1
0
-1
-2
-3
-4
-5
0 50 100 150 200 250
Time step
15
10
Surface elevation at "D".
5
0
-5
-10
-15
0 50 100 150 200 250
Time step
15
10
First order term.
5
0
-5
-10
-15
0 50 100 150 200 250
Time step
15
10
Second order term.
5
0
-5
-10
-15
0 50 100 150 200 250
Time step
Fig. 2. Generation of standing waves in fluid domain with sloping bottom.
10
Displacement [cm]
Displacement [cm]
Displacement [cm]
Displacement [cm]
1.50
f
1
x
Resultant of pressure on CD = r� �gd2 f
1.25 2
1.00
0.75
0.50
0 50 100 150 200 250
Time step
2.0
f 1
x
Moment of pressure on CD = r� �gd3 f
6
1.5
1.0
0.5
0.0
0 50 100 150 200 250
Time step
Fig. 2. Continued.
REFERENCES
1. Bathe K. J. (1982) Finite Element Procedures in Engineering Analysis, Prentice-Hall,
Inc., Englewood Cliffs, New Jersey.
2. Szmidt J.K. (2001) Zbadanie dyskretnego rozwiązania potencjalnego ruchu falowego
cieczy w zmiennych Lagrange a, oprac. wewn. IBW-PAN, str. 1-31.
3. Wehausen J. V. and Laitone E. V. (1960) Surface Waves in Encyclopaedia of Physics
ed. by Flugge S., vol. IX, Fluid Dynamics III, Springer Verlag, Berlin.
4. Wilde P. and Wilde M. (2001) On the Generation of Water Waves in a Flume,
Archives of Hydro-Engineering and Environmental Mechanics, vol. 48, No.4, pp.69-
83.
Fale stojące o skończonej amplitudzie w wodzie
o zmiennej głębokości
STRESZCZENIE
W pracy dyskutuje się rozwiązanie dyskretne zagadnienia narastającej fali stojącej o
skończonej amplitudzie. W analizowanym modelu jest to wymuszenie ruchu cieczy w ob-
szarze skończonym za pomocą generatora typu sztywnego tłoka. Częstość wymuszenia od-
powiada częstości rezonansowej obszaru cieczy. Cały układ rozpoczyna ruch ze stanu spo-
koju. Rozwiązanie zbudowano za pomocą metody małego parametru, przy ograniczeniu
rozwiązania do drugiego rzędu przybliżenia. Rozważania zilustrowano przykładem licz-
bowym dotyczącym pionowego falochronu w wodzie o liniowo zmiennej głębokości.
11
Force
Moment