- Email: [email protected]
Long time evolution of standing gravity waves in deep water
ELSEVIER
Wave
Motion 23 ( 1996) 279-287
Long time evolution of standing gravity waves in deep water Makoto Okamura Research Institute for Applied Mechanics, Kyushu Universir?: Kasuga. Fukuoka 816. Japan
Received 9 March 1995; revised 14 November
1995
Abstract A numerical study of the long time evolution of standing gravity waves with a breaking model is presented. We propose two simple and heuristic breaking models. The models simulate numerically the frequency downshift phenomenon which Lake et al. (1977) have reported experimentally. Our results show a new recurrence after the frequency downshift. A new method to obtain a standing wave profile is also proposed.
1. Introduction Since the experiment by Lake et al. [l] many authors have tried to explain the frequency downshift phenomenon theoretically. Uchiyama and Kawahara [2] use the Dysthe equation, which is the fourth order nonlinear Schriidinger equation, with a damping term of induced mean flow. They show a permanent frequency downshift for certain cases. Kato and Oikawa [3] use the Dysthe equation with a dissipative term, which shows an interaction between waves and the mean flow, to simulate the frequency downshift phenomenon. Both results show a permanent frequency downshift but we cannot derive the additional terms from the basic equations for water waves. Trulsen and Dysthe [4] use the Dysthe equation with a heuristic damping term. The damping term has the effect that the amplitude A decreases the critical value A0 very rapidly if A > Ao. Their result also shows a permanent frequency downshift. Kharif [5] studied the nonlinear evolution of a uniform wave train with effects of viscosity and surface tension using a high order spectral method under the basic equations for water waves. His results show a permanent frequency downshift similar to the results by Trulsen and Dysthe [4] and Hara and Mei [6]. Melville and Rapp 171 investigated simultaneous measurements of the surface displacement and the horizontal velocity on the surface in detail. They found that the velocity change is much larger than the surface-displacement change in breaking. In other words the velocity decreases during breaking but breaking has little influence on the surface displacement. There are some studies of the frequency downshift phenomenon using the weak nonlinear equation with an additional term but very few under the basic equations for water waves. As the additional term is not based on the water wave equations it is difficult to extend the weak nonlinear equation to a full nonlinear equation. In this paper we use the basic equations for water waves including a wave breaking model. The essence of the model is the above finding by Melville and Rapp [7]. The subject of investigation is not a travelling wave but a standing wave because it is easier to make a breaking model. 0165.2125/96/$09.50 0 1996 Elsevier Science SSDIOl65-2125(95)00060-7
B.V. All rights reserved
280
M. Okamura/
Wave Motion 23 (1996) 279-287
2. Formulation of the problem We assume that fluid is incompressible and inviscid and we treat two-dimensional irrotational motion. The capillary effect on a free surface is neglected and the fluid depth is infinite. Basic equations for velocity potential 4 and surface elevation n are
a440
-
for y -+ -oo, ’
ay
where x and y are the horizontal and vertical coordinates, y being directed upwards, t is the time, g the gravitational acceleration and P the pressure. We define the x-axis as the free surface at rest. Eq. (4) corresponds to the kinematic boundary condition but it is different from the usual kinematic boundary condition. It is derived to investigate the tip of a breaking travelling wave by Longuet-Higgins [8]. The boundary condition expresses that the rate of change of the pressure following a particle on the free surface must vanish. The expression is more complicated than the usual kinematic boundary condition. The reason why we use the boundary condition (4) in this paper is that there is no space-derivative of the free surface in the expression. This is helpful in obtaining a standing wave profile with a Fourier series expansion. Here we have introduced the dimensionless variables using wave number K and wave frequency w of the standing wave as follows: Kx+x,
Ky+y,
$wP_ -&R.
wt-+t,
Then we treat a standing wave with 27~ period for time and space without loss of generality. Note that the new g depends on the unknown quantity w. We solve the above equations numerically to seek a standing wave 6(x, y, t), ij(x, t) under the following symmetry conditions:
ax,y, t) =4(x + *n, y,t),
(6)
bk
(7)
y9 t>
=4(-x, y, t),
6(X, y, t) = 6k 4<% Y?t) = -6k
y, f + 27r),
(8)
y, -t),
(9)
gx, y, t) = -qQ?T - x, y, ?T - I).
(10)
We can express the velocity potential as a Fourier series under the above symmetry conditions,
the basic equations
( 1) and (21, N-2 iG,YJ)
= c
Akj
cos
kx tzkq’sin j t,
(11)
k=O j=2-mod(k.2)
where mod(k, 2) denotes the remainder under integer when k is divided by 2 and N is the maximum expansion. The reason why the summation fork in (11) is up to N - 2 will be stated later.
order of the
M. Okamura/Wave Motion 23 (1996) 279-287
The relation between wave steepness E and the surface elevation
q(x, 0) at r = 0 is as follows.
F = $a, 0) - rl(n. 0)
(12)
? We use the collocation method to obtain a standing wave. At first we produce independent collocation points into (3) and (4). The collocation points are selected as: j-l [i = -lr N-l
j=
1,2 ,...,
i=l,2
,...,
’
281
N-l
.Xi=
’
2
i-l
N_2n,
i=
1,2 . . . . . N-l
equations by substituting
(13)
Note that i-l xi=-l-r. N-2
N-l 2
for *=T,
( 14)
because of the symmetry condition (10). The number of the collocation points is N(N - 1)/2. We substitute ( 11) into (3), (4) and the above collocation points into X, t in them to obtain N (N - 1) independent equations. The total number of the independent equations including the wave steepness condition (12) becomes N( N - 1) f 1. We next consider the number of unknown quantities: Akj, q(xi, tj) and g. Solving a standing wave up to the Nth order expansion for AII = O(S) <( 1.
(15)
the orders of the coefficients AiN = O(sN),
i=l,3,5
ANi =O(SN+‘),
=O(S
Nfl)
,...,
i = 1,3,5,.
Ai,N-1 = 0(8N-‘), AN-1.i
Akj for which k or j is equal to N or N - 1 are:
i = 0,2,4, 3
i=2,4
,...,
N-2; . . , N; . . . , N - 3; N-l.
(16)
The number of the unknown quantities Akj which are up to the Nth order of S is N(N - 1)/2. Here the coefficients ANi, AN-l,; are omitted because of the higher order than N. This is why the summation for k in (11) is up to N - 2. The number of other unknown quantities g, f(xi, tj) is N (N - 1)/2 + 1. So the total number of the unknown quantities is N (N - 1) + 1, and is the same as that of the independent equations. We solve the N (N - 1) + 1 nonlinear simultaneous equations by the Newton method. The initial solution of the iteration is a second order nonlinear standing wave as follows: Al1 = --E,
ij(x.r)
A02 = $e2,
otherwise
= ECOSXCOSt + ~(Cos2*cos2t
g = 1 + i&2.
Akj = 0,
(17)
+cos2X),
(18) (19)
We stop the iteration if the maximum of the differences between unknown quantities before an iteration and unknown quantities after the iteration is smaller than 10-15. All the results here will use .s = 0.15, 0.2 and N = 15. If q is also expanded with Fourier series there is a problem with convergence. Schwartz and Whitney [9] show that the convergence radius is about 0.3. Therefore the convergence of the solution is not good even for E = 0.2 in their method without Pad6 approximants. Our new method has no difficulty of the convergence in obtaining standing waves.
282
M. Okamura/Wave
Motion 23 (1996) 279-287
3. Time evolution of standing waves without modelling In this section we check accuracy of the 15th order standing wave profile obtained in Section 2 by performing time integration of the standing wave. We compare a result of time evolution of the standing wave including small disturbances with a result under the linear stability theory. Recurrence is also presented similar to the case for travelling waves. We use an efficient boundary-integral method developed by Dold and Peregrine [lo] for time evolution of water waves. At an initial time t = 0 there are 10 standing waves with 27r wavelength without motion. It is sufficient to consider the region 0 I x ( 10~ due to symmetry of standing waves. All cases here are calculated in the region even for the case without disturbances. Note that we re-introduce the different dimensionless variables from (5) in Sections 3 and 4. Kx-+x,
Ky+y,
&%+t,
where K denotes the wave number of the standing waves considered. We check the accuracy of the scheme by performing time evolution of 10 standing waves with wave steepness 0.2. The time series of the first harmonic of the surface elevation oscillates with a time period. The information required is the envelope of the time series. We calculate the time series of the harmonic at each time step (the difference is 0.008) then we pick up its peak to evaluate the envelope. The envelope of the time series should be constant forever for standing waves without disturbances. Fig. 1 shows the time evolution of the envelope of time series of the first harmonic component W of the carrier wave profile with wave steepness 0.2 for 0 -C t < 5000. The envelope is constant within 5 x 10e6. The total energy remains constant within the relative error 5 x 10-7. The first check on accuracy of the profile of standing waves and the time evolution scheme is satisfactory.
Fig. 1. The time evolution 0
of the envelope
of time series of the first harmonic
component
of the carrier wave with wave steepness
0.2 for
Fig. 2. The time evolution of the envelope of time series of the five Fourier components of surface elevation: cosx, cosO.9~. cosO.8~. cos 1.1x and cos 1.2.x. We express the five components as W, L1, L2, HI and Hz, respectively. E = 0.15, ‘~1 = cy2 = 10e3. The theoretical growth rate is 1.94 x lop3 which is related to L 1 and HI.
M. Okamura/Wave
Motion 23 (1996) 279-287
283
Fig. 2 shows the time evolution of the envelope of the time series of five Fourier components of surface elevation: cos X, cos 0.9x, cos 0.8x, cos 1. lx and cos 1.2x. We express the five components as W, L I, Lz, HI and Hz, respectively. Wave steepness E is chosen as 0.15. The initial condition is 4(x. 0) = 0,
(21)
77(X*0) = 6(x, 0) + ats(cos
(22)
1.1x + ~0~0.9~) + (~~E(COS1.2x + cosO.8~),
where {(x, 0) is the surface profile of a standing wave with wave steepness 0.15 at t = 0, and art and ~2 are the parameters which denote the amplitude of the disturbance and are set at 10-j. Okamura [ 1 l] shows that the disturbance related to ~1 is unstable and its growth rate is 1.94 x low3 for wave steepness E = 0.15 while the disturbance related to (~2 is neutrally stable. The result of the linear stability analysis is also shown in Fig. 2. The growth rate agrees well with the result from the time evolution and the time evolution simulates the neutral stability for the disturbance related to CYZ.The curves of the envelope of the time evolution oscillate slowly. The results remain unchanged even if we set the time step to be a half. One of the reasons why they oscillate slowly is that the added disturbance is not the eigenfunction related to the growth rate 1.94 x 10e3 but we do not check because it is enough to check accuracy of the profile of standing waves and the time evolution scheme. The above relsults of the checks are satisfactory. Fig. 3 shows recurrence which is very similar to the case of travelling waves. Wave steepness E is chosen as 0.15. The initial conditions are (2 1) and (22) for at = (~2 = 0.1. The total energy remains constant within the relative error 5 x 10P6. The meaning of W, L1,t2,HI and H2 is the same as that in Fig. 2. Fig. 4 shows the time evolution of the envelope of time series of the five Fourier components of surface elevation. Wave steepness is chosen as 0.2. The initial conditions are (21) and (22) for CXI= 0.1, CQ = 0.01. The computation breaks down at t x 623 because of the discontinuity of the derivative of the surface profile. The meaning of W, L 1,Lz,H1 and H2 is the same as that in Fig. 2. The computation is broken down before it recurs if the wave steepness is larger than a critical value.
Fig. 3. The time evolution
of the envelope
of time series of the five Fourier components
0 Fig. 4. The time evolution
200
400
of the envelope of time series of the five Fourier components
600
of surface elevation:
E = 0.15, ‘~1 = (~2 = 0.1.
t
of surface elevation:
E = 0.2. (Y, = 0.
I, (~2= 0.01.
M. Okamura/Wave Motion 23 (1996) 279-287
284
4. Two breaking models We introduce two models for wave breaking: the first one emphasizes Melville and Rapp [7] and the second is more realistic.
the essence of the experimental
result by
4.1. Model 1 The fundamental idea is based on the result of the experiment by Melville and Rapp who showed that a rapid increase and decrease of fluid velocity are observed during breaking but the signature of the breaking event in surface displacement is not clear. The essence of their result is that the kinematic energy decreases but the potential energy remains constant during breaking. One of the simple ways to model their result is a manipulation of only the velocity potential 4 at time to when the surface elevation modulates greatly. For example the relation between the new velocity potential &+,(x, fn) and the old velocity potential &td(X, to) is taken as $new(X, to) = @old(X, fO)(l - coso.1(x
(23)
-x0))/%
where xu denotes the place where the surface elevation modulates most. However, we do not manipulate the surface elevation n at tn. Almost all the kinematic energy near the maximum surface elevation is removed suddenly. The initial conditions are (21) and (22) with E = 0.2, ot = 0.1 and o2 = 0.01, which is the same as that in Fig. 4. The breaking time to is selected as 616. The surface profile is shown in Fig. 5 at t = 616, which suggests x0 = 0. We manipulate the velocity potential such as (23) only at the time. As a result there is a decrease of 90% in the kinematic energy at that time. A typical time evolution of the amplitude of the five Fourier components for surface elevation is shown in Fig. 6. The meaning of W, L 1, L2, HI and H2 is the same as that in Fig. 2. The amplitude of the carrier component W decreases and the amplitude of the side-band components increases for r < 616. The standing wave modulates greatly at I z 616. After the breaking the lower side-band L 1 remains
1.
rl
I
I
I 20
I
Ol
-‘0
10
Fig. 5. The surface profile at
Fig. 6. The time evolution
I
I1 X30
t = 616: F = 0.2, (~1 = 0.1. q = 0.01.
of the envelope of time series of the five Fourier components
of surface elevation:
E = 0.2. LYE= 0. I, q
= 0.01.
285
M. Okamura/ Wave Motion 23 (1996) 279-287
Fig. 7. The time evolution of the envelope of time series of the five Fourier components of surface elevation: A vertical bar at the top denotes the time at which we manipulate the velocity potential.
F = 0.2. LY1 =
0.1. ~2 = 0.0 I.
constant with the largest amplitude for 6 16 < t -C 1100. We can regard the behaviour for t > I 100 as a modulation of the carrier wave L 1 with the upper side-band W and the lower side-band L2. The new carrier wave L 1 is born at t = 616 because of the breaking. This very simple and rough model simulates the frequency downshift phenomenon. Note that we treat frequency in experiments but wave number in theories. Fig. 6 shows the discontinuity of the envelope of the time series of Fourier components at t = 616. Of course it does not mean the discontinuity of the time series of Fourier components of the surface elevation at t = 616. Because of the manipulation of velocity potential the time derivative of surface elevation has the discontinuity at t = 6 I6 and the envelope of the time series of the surface elevation has the discontinuity. At t = 616 nine tenth of the kinematic energy is removed and the total energy just after breaking is reduced to 60% of that just before breaking. As the energy reduction is much more slow in experiments this model is not realistic. 4.2. Model II In this section we present another model (model II) which is more realistic than model I. We assume that the breaking associated with frequency downshift is weak and the energy reduction is very small and the amplitude of the surface elevation remains constant during the breaking process, in other words the main effect of breaking is a change of velocity potential. Melville and Rapp [7] support this assumption by experiment. We describe the above assumption with a mathematical expression. We cut higher Fourier components than the Mth component of velocity potential 4 when the maximum surface inclination is larger than a critical angle &. We can see wave breaking near the highest crest in experiments. We cut higher Fourier components in order to model a wave breaking near the highest crest which is associated with higher Fourier components. If the velocity potential $,,td just before the manipulation for 10 standing waves is
(24) the velocity potential Gnew j ust after the manipulation
&vC~. to) =
eu,,
cos(nx/lO),
rt=o
is
(25)
M. Okamura/Wave
286
Motion 23 (1996) 279-287
where to is the time at which the maximum surface inclination is larger than the critical angle 0,. The maximum angle is checked at each time step whose difference is about 0.008. The energy reduction of a time step is very small but the manipulation number is large. For example the number is about 5000 in the case of Fig. 7. A typical result is shown in Fig. 7. The initial conditions are (21) and (22) with F = 0.2, CYI= 0.1 and (112= 0.01, which is the same as that in Fig. 4. The cut number M is 18 and the critical angle 0, is 2 lo, A vertical bar at the top in Fig. 7 denotes the time at which we manipulate the velocity potential. The total energy at t = 5000 is reduced to 59% of that at the initial time. The time evolution of the Fourier modes in this figure is the same as that in Fig. 4 for t < 500 because the maximum inclination angle of the free surface is smaller than the critical angle 0,. There are a lot of vertical bars for 700 < t < 900, which shows that breaking occurs during the time. Then the first lower side-band remains constant and to be the largest for 800 < t -c 1400. This is the frequency downshift phenomenon. W and L2 correspond to the first upper- and lower-side bands of the new carrier L 1 for 800 -C t. The figure shows recurrence of the carrier L 1. Note that the results for model I and II are similar. A typical surface profile is similar to that in Fig. 5 when the modulation is fairly large. Breaking occurs near the peak of the large amplitude in experiments. As the larger amplitude contains the higher Fourier modes of the velocity potential, cutting the higher Fourier modes plays a role in the reduction of the fluid velocity near the peak of the large amplitude.
5. Concluding
remarks
We have presented the two simple models to simulate the frequency downshift phenomenon. The basic idea is the experimental result by Melville and Rapp that a rapid increase and decrease of the fluid velocity are observed during breaking but the signature of the breaking event in the surface displacement is not clear. The idea is emphasized in model I. Model I, which is very simple but simulates frequency downshift well, states that a rapid energy damping is one of the main reasons for frequency downshift. At an initial time (t = 0) there are 10 waves in the region (0 5 x ( 20~). The number of waves reduces to nine when the waves modulate largely at t z 600 (see Figs. 4 and 5). Nine nearly standing waves restart without motion (4 = 0) at t = 616. In a sense it is trivial to simulate the downshift phenomenon in model I. Model II describes the above experimental result with a more accurate mathematical expression. There are two parameters: the cutting number M and the critical angle 0,. As we do not know M and 0,, they are set as large as possible within a stable time evolution. It is also possible to simulate travelling waves with model II. The mechanism of model II is as follows. Uniform standing waves with small disturbances modulate due to Benjamin-Feir instability. Then breaking occurs near the peak of the largest crest, which means that wave energy near the peak decreases. As the larger crest contains the higher Fourier modes, cutting the higher Fourier modes corresponds to energy reduction near the highest peak. The above discussion is similar to that in Trulsen and Dysthe [4]. The time evolution of standing waves without modelling shows a recurrence in Fig. 3. The lower side-band Lt is larger than the upper side-band Ht at t x 1000 when the waves modulate largely. The difference between L 1 and HI comes from the high nonlinearity (i.e. the higher Fourier modes) because we cannot see the difference under the nonlinear Schrbdinger equation but we can under the Dysthe equation. Afterwards the waves become uniform if their steepness is small enough. The most important factor is the high nonlinearity (the higher Fourier modes) for recurrence. So if the higher Fourier modes are reduced at t x 1000, a mechanism to recur may not work well. The lower side-band L 1 keeps largest for some time. These models are heuristic and give no deep understanding of breaking. However it is valuable for an understanding of breaking to simulate frequency downshift under the basic equations for water waves.
M. Okamura/Wave
Motion 23 (1996) 279-287
287
Acknowledgements The author is indebted to Dr. J.W. Dold and Professor D.H. Peregrine for permission scheme of water waves.
to use their time evolution
References [ I ] B.M. Lake, H.C. Yuen, H. Rungaldier and W.E. Ferguson, Nonlinear deep-water waves: Theory and experiment. Part 2. Evolution of a continuous wave train, J. Fluid Mech. 83,49-74 (1977). [2] Y. Uchiyama and T. Kawahara, A possible mechanism for frequency down-shift in nonlinear wave modulation, Wave Motion 20, 99-l IO (1994). [3] Y. Kato and M. Oikawa, Frequency down-shift of a nonlinear modulational wavetrain, Eng. Sci. Reports, Kyushu Univ. IS. 3053 1I (1993) [in Japanese]. [4] K. Trulsen and K.B. Dysthe, Frequency down-shift through self modulation and breaking, Water Wave Kinematics, 561-572 (1990). 151 C. Kharif, Subharmonic transition of a nonlinear short gravity wave train on deep water, Proc. Nonlinear Water Waves Workshop, University of Bristol, 54-61 (1991). [h] T. Hara and C.C. Mei, Frequency downshift in narrowbanded surface waves under the influence of wind, J. Fluid Mech. 230,429-477 (1991). 171 W.K. Melville and R.J. Rapp, The surface velocity field in steep and breaking waves, J. Fluid Mech. 189, 1-22 (1988). [S] MS. Longuet-Higgins, A technique for time-dependent free-surface flows, Proc. Roy. Sot. London, A371,441-451 (1980). [9] L.W. Schwartz and A.K. Whitney, A semi-analytic solution for nonlinear standing waves in deep water. J. Fluid Mech. 107. 147171 (1981). IO) J.W. Dold and D.H. Peregrine, Water-wave modulation, Coastal Eng. 163-175 (1986). I I] M. Okamura, Maximum wave steepness and instabilities of finite-amplitude standing waves, Fluid Dyn. Res. I. 201-214 (1986).
Wave
Motion 23 ( 1996) 279-287
Long time evolution of standing gravity waves in deep water Makoto Okamura Research Institute for Applied Mechanics, Kyushu Universir?: Kasuga. Fukuoka 816. Japan
Received 9 March 1995; revised 14 November
1995
Abstract A numerical study of the long time evolution of standing gravity waves with a breaking model is presented. We propose two simple and heuristic breaking models. The models simulate numerically the frequency downshift phenomenon which Lake et al. (1977) have reported experimentally. Our results show a new recurrence after the frequency downshift. A new method to obtain a standing wave profile is also proposed.
1. Introduction Since the experiment by Lake et al. [l] many authors have tried to explain the frequency downshift phenomenon theoretically. Uchiyama and Kawahara [2] use the Dysthe equation, which is the fourth order nonlinear Schriidinger equation, with a damping term of induced mean flow. They show a permanent frequency downshift for certain cases. Kato and Oikawa [3] use the Dysthe equation with a dissipative term, which shows an interaction between waves and the mean flow, to simulate the frequency downshift phenomenon. Both results show a permanent frequency downshift but we cannot derive the additional terms from the basic equations for water waves. Trulsen and Dysthe [4] use the Dysthe equation with a heuristic damping term. The damping term has the effect that the amplitude A decreases the critical value A0 very rapidly if A > Ao. Their result also shows a permanent frequency downshift. Kharif [5] studied the nonlinear evolution of a uniform wave train with effects of viscosity and surface tension using a high order spectral method under the basic equations for water waves. His results show a permanent frequency downshift similar to the results by Trulsen and Dysthe [4] and Hara and Mei [6]. Melville and Rapp 171 investigated simultaneous measurements of the surface displacement and the horizontal velocity on the surface in detail. They found that the velocity change is much larger than the surface-displacement change in breaking. In other words the velocity decreases during breaking but breaking has little influence on the surface displacement. There are some studies of the frequency downshift phenomenon using the weak nonlinear equation with an additional term but very few under the basic equations for water waves. As the additional term is not based on the water wave equations it is difficult to extend the weak nonlinear equation to a full nonlinear equation. In this paper we use the basic equations for water waves including a wave breaking model. The essence of the model is the above finding by Melville and Rapp [7]. The subject of investigation is not a travelling wave but a standing wave because it is easier to make a breaking model. 0165.2125/96/$09.50 0 1996 Elsevier Science SSDIOl65-2125(95)00060-7
B.V. All rights reserved
280
M. Okamura/
Wave Motion 23 (1996) 279-287
2. Formulation of the problem We assume that fluid is incompressible and inviscid and we treat two-dimensional irrotational motion. The capillary effect on a free surface is neglected and the fluid depth is infinite. Basic equations for velocity potential 4 and surface elevation n are
a440
-
for y -+ -oo, ’
ay
where x and y are the horizontal and vertical coordinates, y being directed upwards, t is the time, g the gravitational acceleration and P the pressure. We define the x-axis as the free surface at rest. Eq. (4) corresponds to the kinematic boundary condition but it is different from the usual kinematic boundary condition. It is derived to investigate the tip of a breaking travelling wave by Longuet-Higgins [8]. The boundary condition expresses that the rate of change of the pressure following a particle on the free surface must vanish. The expression is more complicated than the usual kinematic boundary condition. The reason why we use the boundary condition (4) in this paper is that there is no space-derivative of the free surface in the expression. This is helpful in obtaining a standing wave profile with a Fourier series expansion. Here we have introduced the dimensionless variables using wave number K and wave frequency w of the standing wave as follows: Kx+x,
Ky+y,
$wP_ -&R.
wt-+t,
Then we treat a standing wave with 27~ period for time and space without loss of generality. Note that the new g depends on the unknown quantity w. We solve the above equations numerically to seek a standing wave 6(x, y, t), ij(x, t) under the following symmetry conditions:
ax,y, t) =4(x + *n, y,t),
(6)
bk
(7)
y9 t>
=4(-x, y, t),
6(X, y, t) = 6k 4<% Y?t) = -6k
y, f + 27r),
(8)
y, -t),
(9)
gx, y, t) = -qQ?T - x, y, ?T - I).
(10)
We can express the velocity potential as a Fourier series under the above symmetry conditions,
the basic equations
( 1) and (21, N-2 iG,YJ)
= c
Akj
cos
kx tzkq’sin j t,
(11)
k=O j=2-mod(k.2)
where mod(k, 2) denotes the remainder under integer when k is divided by 2 and N is the maximum expansion. The reason why the summation fork in (11) is up to N - 2 will be stated later.
order of the
M. Okamura/Wave Motion 23 (1996) 279-287
The relation between wave steepness E and the surface elevation
q(x, 0) at r = 0 is as follows.
F = $a, 0) - rl(n. 0)
(12)
? We use the collocation method to obtain a standing wave. At first we produce independent collocation points into (3) and (4). The collocation points are selected as: j-l [i = -lr N-l
j=
1,2 ,...,
i=l,2
,...,
’
281
N-l
.Xi=
’
2
i-l
N_2n,
i=
1,2 . . . . . N-l
equations by substituting
(13)
Note that i-l xi=-l-r. N-2
N-l 2
for *=T,
( 14)
because of the symmetry condition (10). The number of the collocation points is N(N - 1)/2. We substitute ( 11) into (3), (4) and the above collocation points into X, t in them to obtain N (N - 1) independent equations. The total number of the independent equations including the wave steepness condition (12) becomes N( N - 1) f 1. We next consider the number of unknown quantities: Akj, q(xi, tj) and g. Solving a standing wave up to the Nth order expansion for AII = O(S) <( 1.
(15)
the orders of the coefficients AiN = O(sN),
i=l,3,5
ANi =O(SN+‘),
=O(S
Nfl)
,...,
i = 1,3,5,.
Ai,N-1 = 0(8N-‘), AN-1.i
Akj for which k or j is equal to N or N - 1 are:
i = 0,2,4, 3
i=2,4
,...,
N-2; . . , N; . . . , N - 3; N-l.
(16)
The number of the unknown quantities Akj which are up to the Nth order of S is N(N - 1)/2. Here the coefficients ANi, AN-l,; are omitted because of the higher order than N. This is why the summation for k in (11) is up to N - 2. The number of other unknown quantities g, f(xi, tj) is N (N - 1)/2 + 1. So the total number of the unknown quantities is N (N - 1) + 1, and is the same as that of the independent equations. We solve the N (N - 1) + 1 nonlinear simultaneous equations by the Newton method. The initial solution of the iteration is a second order nonlinear standing wave as follows: Al1 = --E,
ij(x.r)
A02 = $e2,
otherwise
= ECOSXCOSt + ~(Cos2*cos2t
g = 1 + i&2.
Akj = 0,
(17)
+cos2X),
(18) (19)
We stop the iteration if the maximum of the differences between unknown quantities before an iteration and unknown quantities after the iteration is smaller than 10-15. All the results here will use .s = 0.15, 0.2 and N = 15. If q is also expanded with Fourier series there is a problem with convergence. Schwartz and Whitney [9] show that the convergence radius is about 0.3. Therefore the convergence of the solution is not good even for E = 0.2 in their method without Pad6 approximants. Our new method has no difficulty of the convergence in obtaining standing waves.
282
M. Okamura/Wave
Motion 23 (1996) 279-287
3. Time evolution of standing waves without modelling In this section we check accuracy of the 15th order standing wave profile obtained in Section 2 by performing time integration of the standing wave. We compare a result of time evolution of the standing wave including small disturbances with a result under the linear stability theory. Recurrence is also presented similar to the case for travelling waves. We use an efficient boundary-integral method developed by Dold and Peregrine [lo] for time evolution of water waves. At an initial time t = 0 there are 10 standing waves with 27r wavelength without motion. It is sufficient to consider the region 0 I x ( 10~ due to symmetry of standing waves. All cases here are calculated in the region even for the case without disturbances. Note that we re-introduce the different dimensionless variables from (5) in Sections 3 and 4. Kx-+x,
Ky+y,
&%+t,
where K denotes the wave number of the standing waves considered. We check the accuracy of the scheme by performing time evolution of 10 standing waves with wave steepness 0.2. The time series of the first harmonic of the surface elevation oscillates with a time period. The information required is the envelope of the time series. We calculate the time series of the harmonic at each time step (the difference is 0.008) then we pick up its peak to evaluate the envelope. The envelope of the time series should be constant forever for standing waves without disturbances. Fig. 1 shows the time evolution of the envelope of time series of the first harmonic component W of the carrier wave profile with wave steepness 0.2 for 0 -C t < 5000. The envelope is constant within 5 x 10e6. The total energy remains constant within the relative error 5 x 10-7. The first check on accuracy of the profile of standing waves and the time evolution scheme is satisfactory.
Fig. 1. The time evolution 0
of the envelope
of time series of the first harmonic
component
of the carrier wave with wave steepness
0.2 for
Fig. 2. The time evolution of the envelope of time series of the five Fourier components of surface elevation: cosx, cosO.9~. cosO.8~. cos 1.1x and cos 1.2.x. We express the five components as W, L1, L2, HI and Hz, respectively. E = 0.15, ‘~1 = cy2 = 10e3. The theoretical growth rate is 1.94 x lop3 which is related to L 1 and HI.
M. Okamura/Wave
Motion 23 (1996) 279-287
283
Fig. 2 shows the time evolution of the envelope of the time series of five Fourier components of surface elevation: cos X, cos 0.9x, cos 0.8x, cos 1. lx and cos 1.2x. We express the five components as W, L I, Lz, HI and Hz, respectively. Wave steepness E is chosen as 0.15. The initial condition is 4(x. 0) = 0,
(21)
77(X*0) = 6(x, 0) + ats(cos
(22)
1.1x + ~0~0.9~) + (~~E(COS1.2x + cosO.8~),
where {(x, 0) is the surface profile of a standing wave with wave steepness 0.15 at t = 0, and art and ~2 are the parameters which denote the amplitude of the disturbance and are set at 10-j. Okamura [ 1 l] shows that the disturbance related to ~1 is unstable and its growth rate is 1.94 x low3 for wave steepness E = 0.15 while the disturbance related to (~2 is neutrally stable. The result of the linear stability analysis is also shown in Fig. 2. The growth rate agrees well with the result from the time evolution and the time evolution simulates the neutral stability for the disturbance related to CYZ.The curves of the envelope of the time evolution oscillate slowly. The results remain unchanged even if we set the time step to be a half. One of the reasons why they oscillate slowly is that the added disturbance is not the eigenfunction related to the growth rate 1.94 x 10e3 but we do not check because it is enough to check accuracy of the profile of standing waves and the time evolution scheme. The above relsults of the checks are satisfactory. Fig. 3 shows recurrence which is very similar to the case of travelling waves. Wave steepness E is chosen as 0.15. The initial conditions are (2 1) and (22) for at = (~2 = 0.1. The total energy remains constant within the relative error 5 x 10P6. The meaning of W, L1,t2,HI and H2 is the same as that in Fig. 2. Fig. 4 shows the time evolution of the envelope of time series of the five Fourier components of surface elevation. Wave steepness is chosen as 0.2. The initial conditions are (21) and (22) for CXI= 0.1, CQ = 0.01. The computation breaks down at t x 623 because of the discontinuity of the derivative of the surface profile. The meaning of W, L 1,Lz,H1 and H2 is the same as that in Fig. 2. The computation is broken down before it recurs if the wave steepness is larger than a critical value.
Fig. 3. The time evolution
of the envelope
of time series of the five Fourier components
0 Fig. 4. The time evolution
200
400
of the envelope of time series of the five Fourier components
600
of surface elevation:
E = 0.15, ‘~1 = (~2 = 0.1.
t
of surface elevation:
E = 0.2. (Y, = 0.
I, (~2= 0.01.
M. Okamura/Wave Motion 23 (1996) 279-287
284
4. Two breaking models We introduce two models for wave breaking: the first one emphasizes Melville and Rapp [7] and the second is more realistic.
the essence of the experimental
result by
4.1. Model 1 The fundamental idea is based on the result of the experiment by Melville and Rapp who showed that a rapid increase and decrease of fluid velocity are observed during breaking but the signature of the breaking event in surface displacement is not clear. The essence of their result is that the kinematic energy decreases but the potential energy remains constant during breaking. One of the simple ways to model their result is a manipulation of only the velocity potential 4 at time to when the surface elevation modulates greatly. For example the relation between the new velocity potential &+,(x, fn) and the old velocity potential &td(X, to) is taken as $new(X, to) = @old(X, fO)(l - coso.1(x
(23)
-x0))/%
where xu denotes the place where the surface elevation modulates most. However, we do not manipulate the surface elevation n at tn. Almost all the kinematic energy near the maximum surface elevation is removed suddenly. The initial conditions are (21) and (22) with E = 0.2, ot = 0.1 and o2 = 0.01, which is the same as that in Fig. 4. The breaking time to is selected as 616. The surface profile is shown in Fig. 5 at t = 616, which suggests x0 = 0. We manipulate the velocity potential such as (23) only at the time. As a result there is a decrease of 90% in the kinematic energy at that time. A typical time evolution of the amplitude of the five Fourier components for surface elevation is shown in Fig. 6. The meaning of W, L 1, L2, HI and H2 is the same as that in Fig. 2. The amplitude of the carrier component W decreases and the amplitude of the side-band components increases for r < 616. The standing wave modulates greatly at I z 616. After the breaking the lower side-band L 1 remains
1.
rl
I
I
I 20
I
Ol
-‘0
10
Fig. 5. The surface profile at
Fig. 6. The time evolution
I
I1 X30
t = 616: F = 0.2, (~1 = 0.1. q = 0.01.
of the envelope of time series of the five Fourier components
of surface elevation:
E = 0.2. LYE= 0. I, q
= 0.01.
285
M. Okamura/ Wave Motion 23 (1996) 279-287
Fig. 7. The time evolution of the envelope of time series of the five Fourier components of surface elevation: A vertical bar at the top denotes the time at which we manipulate the velocity potential.
F = 0.2. LY1 =
0.1. ~2 = 0.0 I.
constant with the largest amplitude for 6 16 < t -C 1100. We can regard the behaviour for t > I 100 as a modulation of the carrier wave L 1 with the upper side-band W and the lower side-band L2. The new carrier wave L 1 is born at t = 616 because of the breaking. This very simple and rough model simulates the frequency downshift phenomenon. Note that we treat frequency in experiments but wave number in theories. Fig. 6 shows the discontinuity of the envelope of the time series of Fourier components at t = 616. Of course it does not mean the discontinuity of the time series of Fourier components of the surface elevation at t = 616. Because of the manipulation of velocity potential the time derivative of surface elevation has the discontinuity at t = 6 I6 and the envelope of the time series of the surface elevation has the discontinuity. At t = 616 nine tenth of the kinematic energy is removed and the total energy just after breaking is reduced to 60% of that just before breaking. As the energy reduction is much more slow in experiments this model is not realistic. 4.2. Model II In this section we present another model (model II) which is more realistic than model I. We assume that the breaking associated with frequency downshift is weak and the energy reduction is very small and the amplitude of the surface elevation remains constant during the breaking process, in other words the main effect of breaking is a change of velocity potential. Melville and Rapp [7] support this assumption by experiment. We describe the above assumption with a mathematical expression. We cut higher Fourier components than the Mth component of velocity potential 4 when the maximum surface inclination is larger than a critical angle &. We can see wave breaking near the highest crest in experiments. We cut higher Fourier components in order to model a wave breaking near the highest crest which is associated with higher Fourier components. If the velocity potential $,,td just before the manipulation for 10 standing waves is
(24) the velocity potential Gnew j ust after the manipulation
&vC~. to) =
eu,,
cos(nx/lO),
rt=o
is
(25)
M. Okamura/Wave
286
Motion 23 (1996) 279-287
where to is the time at which the maximum surface inclination is larger than the critical angle 0,. The maximum angle is checked at each time step whose difference is about 0.008. The energy reduction of a time step is very small but the manipulation number is large. For example the number is about 5000 in the case of Fig. 7. A typical result is shown in Fig. 7. The initial conditions are (21) and (22) with F = 0.2, CYI= 0.1 and (112= 0.01, which is the same as that in Fig. 4. The cut number M is 18 and the critical angle 0, is 2 lo, A vertical bar at the top in Fig. 7 denotes the time at which we manipulate the velocity potential. The total energy at t = 5000 is reduced to 59% of that at the initial time. The time evolution of the Fourier modes in this figure is the same as that in Fig. 4 for t < 500 because the maximum inclination angle of the free surface is smaller than the critical angle 0,. There are a lot of vertical bars for 700 < t < 900, which shows that breaking occurs during the time. Then the first lower side-band remains constant and to be the largest for 800 < t -c 1400. This is the frequency downshift phenomenon. W and L2 correspond to the first upper- and lower-side bands of the new carrier L 1 for 800 -C t. The figure shows recurrence of the carrier L 1. Note that the results for model I and II are similar. A typical surface profile is similar to that in Fig. 5 when the modulation is fairly large. Breaking occurs near the peak of the large amplitude in experiments. As the larger amplitude contains the higher Fourier modes of the velocity potential, cutting the higher Fourier modes plays a role in the reduction of the fluid velocity near the peak of the large amplitude.
5. Concluding
remarks
We have presented the two simple models to simulate the frequency downshift phenomenon. The basic idea is the experimental result by Melville and Rapp that a rapid increase and decrease of the fluid velocity are observed during breaking but the signature of the breaking event in the surface displacement is not clear. The idea is emphasized in model I. Model I, which is very simple but simulates frequency downshift well, states that a rapid energy damping is one of the main reasons for frequency downshift. At an initial time (t = 0) there are 10 waves in the region (0 5 x ( 20~). The number of waves reduces to nine when the waves modulate largely at t z 600 (see Figs. 4 and 5). Nine nearly standing waves restart without motion (4 = 0) at t = 616. In a sense it is trivial to simulate the downshift phenomenon in model I. Model II describes the above experimental result with a more accurate mathematical expression. There are two parameters: the cutting number M and the critical angle 0,. As we do not know M and 0,, they are set as large as possible within a stable time evolution. It is also possible to simulate travelling waves with model II. The mechanism of model II is as follows. Uniform standing waves with small disturbances modulate due to Benjamin-Feir instability. Then breaking occurs near the peak of the largest crest, which means that wave energy near the peak decreases. As the larger crest contains the higher Fourier modes, cutting the higher Fourier modes corresponds to energy reduction near the highest peak. The above discussion is similar to that in Trulsen and Dysthe [4]. The time evolution of standing waves without modelling shows a recurrence in Fig. 3. The lower side-band Lt is larger than the upper side-band Ht at t x 1000 when the waves modulate largely. The difference between L 1 and HI comes from the high nonlinearity (i.e. the higher Fourier modes) because we cannot see the difference under the nonlinear Schrbdinger equation but we can under the Dysthe equation. Afterwards the waves become uniform if their steepness is small enough. The most important factor is the high nonlinearity (the higher Fourier modes) for recurrence. So if the higher Fourier modes are reduced at t x 1000, a mechanism to recur may not work well. The lower side-band L 1 keeps largest for some time. These models are heuristic and give no deep understanding of breaking. However it is valuable for an understanding of breaking to simulate frequency downshift under the basic equations for water waves.
M. Okamura/Wave
Motion 23 (1996) 279-287
287
Acknowledgements The author is indebted to Dr. J.W. Dold and Professor D.H. Peregrine for permission scheme of water waves.
to use their time evolution
References [ I ] B.M. Lake, H.C. Yuen, H. Rungaldier and W.E. Ferguson, Nonlinear deep-water waves: Theory and experiment. Part 2. Evolution of a continuous wave train, J. Fluid Mech. 83,49-74 (1977). [2] Y. Uchiyama and T. Kawahara, A possible mechanism for frequency down-shift in nonlinear wave modulation, Wave Motion 20, 99-l IO (1994). [3] Y. Kato and M. Oikawa, Frequency down-shift of a nonlinear modulational wavetrain, Eng. Sci. Reports, Kyushu Univ. IS. 3053 1I (1993) [in Japanese]. [4] K. Trulsen and K.B. Dysthe, Frequency down-shift through self modulation and breaking, Water Wave Kinematics, 561-572 (1990). 151 C. Kharif, Subharmonic transition of a nonlinear short gravity wave train on deep water, Proc. Nonlinear Water Waves Workshop, University of Bristol, 54-61 (1991). [h] T. Hara and C.C. Mei, Frequency downshift in narrowbanded surface waves under the influence of wind, J. Fluid Mech. 230,429-477 (1991). 171 W.K. Melville and R.J. Rapp, The surface velocity field in steep and breaking waves, J. Fluid Mech. 189, 1-22 (1988). [S] MS. Longuet-Higgins, A technique for time-dependent free-surface flows, Proc. Roy. Sot. London, A371,441-451 (1980). [9] L.W. Schwartz and A.K. Whitney, A semi-analytic solution for nonlinear standing waves in deep water. J. Fluid Mech. 107. 147171 (1981). IO) J.W. Dold and D.H. Peregrine, Water-wave modulation, Coastal Eng. 163-175 (1986). I I] M. Okamura, Maximum wave steepness and instabilities of finite-amplitude standing waves, Fluid Dyn. Res. I. 201-214 (1986).