Elliptic and hyperelliptic functions describing the particle motion beneath small-amplitude water waves with constant vorticity

We provide analytic solutions of the nonlinear differential equation system describing the particle paths below small-amplitude periodic gravity waves travelling on a constant vorticity current. We show that these paths are not closed curves. Some solutions can be expressed in terms of Jacobi elliptic functions, others in terms of hyperelliptic functions. We obtain new kinds of particle paths. We make some remarks on the stagnation points which could appear in the fluid due to the vorticity.


Introduction
The present work is confined to two-dimensional water waves with constant vorticity. To gain insight into the motion in such waves we determine analytically the trajectories beneath small-amplitude water waves in constant vorticity flows. In order to get this Lagrangian feature of the flow, that is, the evolution of individual water particles, we will firstly find a solution of the Eulerian system of equations within the framework of small amplitude waves and then we will calculate the solutions of the nonlinear differential equations system which describes the particle motion.
Many of the theoretical results concerning waves on water make the initial assumption of irrotational flow. There are circumstances in which this is well justified but there are cases where it is inappropriate. Waves with vorticity are commonly seen in nature, for example, in shear currents. Tidal flow is a well-known example when constant vorticity flow is an appropriate model (see Da Silva T. A. and Peregrine D. H. [16]). For a discussion of the physical relevance of flows with constant vorticity see also [5]. In 1802, Gerstner [20] constructed an explicit example of a periodic travelling wave in water of infinite depth with a specific non-constant vorticity 1 . The fact that this flow is very special is confirmed also by the fact that this is the only steady flow satisfying the constraint of constant pressure along the streamlines cf. [33]. Gerstner's wave is a two-dimensional wave which adopts the Lagrangian viewpoint, describing the evolution of individual water particles. Its surface profile is symmetric [3]. Beneath Gerstner's wave it is possible to have a motion of the fluid where all particles describe circles with a depth-dependent radius [3], [24].
In 1934, Dubreil-Jacotin [18] considered the problem of the existance of steady periodic water waves with general vorticity. She proved the existence of large classes of small-amplitude water waves with vorticity. For large-amplitude water waves with vorticity, in 2004, Constantin and Strauss [13] proved that, for an arbitrary vorticity distribution and for a given c > 0 and relative mass flux p 0 , there is a global continuum of steady periodic waves travelling at speed c in water of finite depth and such that the horizontal component of the velocity u < c throughout the fluid. The continuum contains waves with u arbitrarily close to the wave speed c. The existence of global continua of smooth solutions for the related problem of periodic waves of infinite depth was proved by Hur [25]. The vorticity does not destroy the symmetry. The construction in [13] assumes that the wave profiles are symmetric. Constantin and Escher [8], [9] proved that the symmetry of the wave profile is not a hypothesis but rather a conclusion when the wave profile is monotone between crest and trough and the vorticity is positive and non-increasing with greater depth. Free from restrictions on the vorticity but requiring a quite precise knowledge of all streamlines in the fluid, Hur proved in [26] that if the wave profile is monotone near the trough and every streamline has a single minimum per wavelength located below the trough, then the steady periodic water waves of finite depth are symmetric. Constantin, Ehrnstöm and Wahlén [7] showed that for an arbitrary vorticity distribution, a steady periodic water wave with a profile that is monotone between crests and troughs has to be symmetric.
Another remarkable feature of the rotational steady waves is that they could contain stagnation points. If (u, v) denotes the velocity field and c the constant speed of the wave, then a point where u = c and v = 0 is called stagnation point. There are interesting problems related to the so-called extreme waves: these are waves with the stagnation points at their crests. Varvaruca [40] proved for a certain class of vorticity functions, the existence of extreme waves and showed that at such a stagnation point the profile of the wave has either a corner of 120 0 or a horizontal tangent. For a recent survey of different aspects of the theory of steady water waves with vorticity see [39]. This paper has interest in finding information about the flow below water waves with constant vorticity, more precisely, we will investigate how the presence of vorticity influences the particle paths. Throughout the hydrodynamics literature, it has been quite common to assume that beneath an irrotational periodic two-dimensional travelling water wave, the particles trace closed, circular or elliptic, orbits. (see for example [34], [35], [17], [32]). While in this first approximation all particle paths appear to be closed, Constantin and Villari showed in [12], using phase-plane considerations for the nonlinear system describing the particle motion, that in linear irrotational periodic gravity water waves no particles trajectory is actually closed, unless the free surface is flat. Similar results hold for the particle trajectories in irrotational deep-water (see Constantin, Ehrnström and Villari [6]), and in irrotational shallow water (see Ionescu-Kruse [27] and [28], Section 5.1). Ionescu-Kruse [27], [28] obtained the exact solutions of the nonlinear differential equation system which describes the particle motion in small-amplitude shallow water waves and showed that there does not exist a single pattern for all particles: depending on the strength of the underling uniform current, some particle trajectories are undulating curves to the right, or to the left, others are loops with forward drift, and others are not physically acceptable, in the last case it seems necessary to study the full nonlinear problem. For the full nonlinear problem, Constantin proved in [4], by analyzing a free boundary problem for harmonic functions in a planar domain, that all water particles in Stokes waves display a forward drift. For an extension of the investigation in [4] to deep-water Stokes waves see Henry [21]. In a very recent paper [14], Constantin and Strauss recovered the results in [4] by a simpler approach and they also investigated the effect of an underlying current on the paths of the particles. While in periodic waves within a period each particle experiences a backward-forward motion with a forward drift, Constantin and Escher showed in [10] that in a solitary water wave there is no backward motion: all particles move in the direction of wave propagation at a positive speed, the direction being upwards or downwards if the particle precedes, respectively, does not precede the wave crest.
There have also been some studies of particle paths for rotational waves. Within the linear theory, by using phase-plane considerations for the nonlinear system describing the particle motion, Ehrnström and Villari [19] found that for positive constant vorticity, the behavior of the streamlines is the same as for the irrotational waves, though the physical particle paths behave differently if the size of the vorticity is large enough. For negative vorticity they showed that in a frame moving with the wave, the fluid contains a cat's-eye vortex (see [36], Ex. 2.4). The paper [42] by Wahlén which contains an existence result for smallamplitude solutions, based on local bifurcation theory, showed also that the predictions for negative vorticity [19] in the linear theory are true. We mention that an alternative approach to the existence result in [42] for small-amplitude steady waves with constant vorticity was very recently proposed by Constantin and Varvaruca [15]. Beside the phase-plane analysis, the exact solutions of the nonlinear system describing the particle motion, allow a better understanding of the dynamics. For small-amplitude shallow-water waves with vorticity and background flow Ionescu-Kruse found in [28] the exact solutions and showed that depending on the relation between the initial data and the constant vorticity some particles trajectories are undulating curves to the right, or to the left, others are loops with forward drift, or with backward drift, others can follow peculiar shapes (see [28], Fig. 7e).
Removing the shallow-water restriction, in the present paper we provide explicit solutions for the nonlinear system describing the motion of the particles beneath small-amplitude gravity waves which propagate on the surface of a constant vorticity flow. In Section 2 we recall the governing equations for gravity water waves. In Section 3 we present their nondimensionalisation and scaling. We present two different scalings, in one the constant vorticity ω 0 is scaled whereas in another one ω 0 remains unscaled. We choose x and z the space coordinates, thus, the sign of the constant vorticity ω 0 is opposite to the sign of the constant vorticity considered if x and y are chosen the space coordinates. In Section 4 we obtain the periodic travelling solutions of the considered linearized problems (see (27), respectively (44)), and the speed of propagation of the linear wave c (see (26), respectively (43)). The solutions are also written in the original physical variables: see (28), (29), respectively (45), (46). We observe that the speed of the wave and the pressure have different expressions in the two linearizations. In Section 5 we find the solutions of the nonlinear differential equation systems (51) and we describe the possible particle trajectories beneath constant vorticity water waves. In the study of the system (51) it is interesting to observe that, for the first linearization, that is, the one made around still water in which the constant vorticity ω 0 is scaled, the sign of the wave speed c will influence the sign of the parameter A which appears in the components u and v of the velocity field. Thus, if we consider left-going waves c < 0, we get A < 0 and if we consider right-going waves c > 0, we obtain A > 0. For the second linearization, that is, the one made around a laminar flow characterized by u = ω 0 z + α, v = 0, α being a constant, we obtain that, independent of the sign of ω 0 , the sign of A depends on the sign of c − h 0 ω 0 − √ gh 0 α, where h 0 is the finite depth and g the constant gravitational acceleration. The expression c − h 0 ω 0 − √ gh 0 α could be regarded as "the speed" of a wave which can be left-going or right-going. In the study of the system (51) a peakon-like trajectory (72) comes up (see also Ionescu-Kruse [31]). This solution contains the arctanh(·) function, having a vertical asymptote in the positive direction ( Figure 2). For this solution u = c and a stagnation point in the fluid appear only for t → ±∞, where the path of the particle has a horizontal tangent. The other solutions of the system (51) are given by (87). We show these solutions are not closed curves. Some of these solutions can be expressed with the aid of the Jacobi elliptic functions, others are expressed with the aid of the hyperelliptic functions. We draw some of the curves obtained for different values of the parameters (see Figure 3, Figure 4, Figure 5). At the end, we make some remarks on the stagnation points inside the fluid.

The water wave problem
The two-dimensional gravity waves on constant vorticity water of finite depth are described by the following boundary value problem: is the velocity field of the water -no motion takes place in the y-direction, p(x, z, t) denotes the pressure, g is the constant gravitational acceleration, p 0 being the constant atmospheric pressure and ω 0 is the constant vorticity. The water moves in a domain with a free upper surface at

Non-dimensionalization and scaling
We non-dimensionalize the set of equations (1) using the undisturbed depth of the water h 0 , as the vertical scale, a typical wavelength λ, as the horizontal scale, and a typical amplitude of the surface wave a (for more details see [32]). Thus, we define the set of non-dimensional variables where, to avoid new notations, we have used the same symbols for the nondimensional variables x, z, η, t, u, v, p on the right-hand side. The partial derivatives u z and v x will be then replaced by and the natural scaling for the vorticity is where we have used the same symbol for the non-dimensional ω 0 on the righthand side. Therefore, in non-dimensional variables (2), (4), the water-wave problem (1) becomes: where we have introduced the amplitude parameter ǫ = a h0 and the shallowness parameter δ = h0 λ . After the non-dimensionalization of the system (1) let us now proceed with the scaling transformation. First we observe that, on z = 1 + ǫη, both v and p are proportional to ǫ. This is consistent with the fact that as ǫ → 0 we must have v → 0 and p → 0. We can consider the following scaling of the non-dimensional variables where we avoided again the introduction of a new notation. For this scaling of u and v, we also get ω 0 → ǫω 0 The water-wave problem (1) writes in non-dimensional scaled variables (2), (4), (6), (7), as By letting ǫ → 0, δ being fixed, we obtain a linear approximation of the problem (8), that is, This linearization is used in [27], [28] for irrotational and constant vorticity shallow water waves, in [29], [30] for capillary-gravity waves and in [31] for constant vorticity gravity waves.
For constant vorticity flows all the streamlines are real-analytic as proved recently in [11]. Thus, for travelling water waves with constant vorticity one can get the analytic validity of the linearization (9). In the case of irrotational water waves, the rigorous analysis of the validity and relevance of the linearizations around some reference states is carried out in [1]. Instead of the scaling (6), we can also consider the following one with α constant. In this case ω 0 remains unscaled. Thus, the water-wave problem (1) writes in non-dimensional scaled variables (2), (4), (10), as By letting ǫ → 0, δ being fixed, we obtain a linear approximation of the problem (11), that is, We observe that the forth equation in the system (12), which represents the vorticity equation, becomes in these scaled variables the vorticty equation for an irrotational flow. The linearization (12) is used in [12] for irrotational gravity waves, in [6] for the corresponding deep-water waves, in [22], [23] for capillary-gravity waves and in [19] for gravity waves over finite depth with constant vorticity. This linearization is around a laminar flow. Such shear flows are characterized by the flat surface, z = 1, corresponding to η = 0, p = 0, v = 0 and u = ω 0 z + α.

Solutions of the linearized problems
Let us solve the linearized systems (9) and (12) and compare their solutions. From the third equation and the forth equation in (9), we obtain that Applying the method of separation of variables, we seek the solution of the equation (13) Substituting (14) into the equation (13), separating the variables and taking into account the expressions of v on the boundaries, that is, the fifth equation and the last equation in (9), we find where K ≥ 0 is a constant that might depend on time. Taking into account (15) and the fourth equation of the system (9), we obtain where F (x, t) is an arbitrary function. The components u and v of the velocity have to fulfill also the third equation in (9), hence, in view of (15) and (16), we get The above relation must hold for all values of x ∈ R, and 0 ≤ z ≤ 1, thus, it follows and Seeking periodic travelling wave solutions for the equation (19), we take K = 2π (20) and we choose the following solution where c represents the non-dimensional speed of propagation of the linear wave and is to be determined. From (18) the function F (x, t) is independent of x, therefore we will denote this function by F (t). Thus, the components of the velocity filed are We return now to the systems (9) in order to find the the expressions of the pressure. Taking into account the first two equations in (9) and the expressions (22) of the velocity field, we obtain On the free surface z = 1 the pressure (23) has to fulfill the sixth equation of the system (9). Hence, in view of (21), we get The above relation must hold for all values x ∈ R, therefore, we get and we provide the non-dimensional speed of the linear wave Summing up, the solution of the linear system (9) is: The solution of this equation can be written into the form c 1 , c 2 being integration constants. We observe that to the linearized problem (9) we can get solutions different from the usual one η(x, t) = cos(K(x − ct)).
with c given by (26). Taking into account (2), (4), (6), (7), we return to the original physical variables. The speed of the wave (26) and the solution (27) become: where is the wave number. The sign minus in (28) indicates a left-going wave. Let us look now at the linearized system (12). From the third equation and the forth equation in (12), we obtain again the equation (13). Applying the method of separation of variables, we seek the solution of this equation in the form (14). Substituting (14) into the equation (13), separating the variables and taking into account the expressions of v on the boundaries, that is, the fifth equation and the last equation in (12) where K ≥ 0 is a constant that might depend on time. Taking into account (31) and the fourth equation of the system (12), we obtain where F(x, t) is an arbitrary function. The components u and v of the velocity have to fulfill also the third equation in (12), hence, in view of (31) and (32), we get The above relation must hold for all values of x ∈ R, and 0 ≤ z ≤ 1, thus, it follows ∂F(x, t) ∂x = 0 and Seeking periodic travelling wave solutions for the equation (35), we take and we choose the following solution where c represents the non-dimensional speed of propagation of the linear wave and is to be determined. From (34) the function F(x, t) is independent of x, therefore we will denote this function by F(t). Thus, the components (32), (31) of the velocity filed are We return now to the systems (12) in order to find the the expressions of the pressure. Taking into account the first two equations in (12) and the expressions (38) of the velocity field, we obtain On the free surface z = 1 the pressure (39) has to fulfill the sixth equation of the system (12). Hence, in view of (37), we get The above relation must hold for all values x ∈ R, therefore, we get and the non-dimensional speed of the linear wave c satisfies the relation Solving this equation we find Summing up, the solution of the linear system (12) is given by with c from (43). Taking into account (2), (4), (10), we return to the original physical variables. The speed of the wave (43) and the solution (44) have in physical variables the following expressions: is the wave number. The solution (45), (46) with α = 0, c 0 = 0 was also obtained in [19]. Comparing (28), (29) with (45), (46) we observe that the speed of the wave and the pressure have different expressions in the two linearizations. The velocity field has in the two linearizations the form: where, for the linearization (9): and for the linearization (12):

Particle trajectories
Let (x(t), z(t)) be the path of a particle in the fluid domain, with location (x(0), z(0)) := (x 0 , z 0 ) at time t = 0. The motion of the particles below the small-amplitude water waves in constant vorticity flows with the velocity field (48), is described by the following differential system The values of A, B, C are either (49) or (50), depending on which linearization we consider. From (49), (50), for any ω 0 = 0, we get in the both cases For the first linearization (9), the sign of A depends on the sign of the wave speed c. Thus, if we choose in (49) the square root with minus, that is, we consider left-going waves, we have A < 0 and if we choose in (49) the square root with plus, that is, we consider right-going waves, we get A > 0.
For the second linearization (12), the sign of A depends on the sign of c − h 0 ω 0 − √ gh 0 α. Looking at the expression (50) of c − h 0 ω 0 − √ gh 0 α, we get that independent of the sign of ω 0 , if we choose in (50) the square root with minus then c − h 0 ω 0 − √ gh 0 α < 0, thus, A < 0, and if we choose in (50) the square root with plus then c − h 0 ω 0 − √ gh 0 α > 0, thus, A > 0. The expression c − h 0 ω 0 − √ gh 0 α could be regarded as "the speed" of a wave, which is left-going if we take the square root with minus in (50), and is right-going if we take the square root with plus in (50).
To study the exact solution of the system (51) it is more convenient to rewrite it in the following moving frame This transformation yields We write the second equation of this system in the form Integrating, we get If kA sin X(t) dt < 0 (57) then Z(t) = 2 arctanh exp kA sin X(t) dt (58) Taking into account the formula: and the expression (58) of Z(t), the first equation of the system (54) becomes where we have denoted by With (57) in view, we have 0 < w < 1 (62) Differentiating with respect to t this relation, we obtain From (63) we have furthermore Thus, taking into account (64), (65), the equation (60) becomes We make the following substitution A being different from zero (52). Differentiating with respect to t this relation, we get We replace (67), (68) into the equation (66) and we obtain the equation A solution of the equation (69) is which, in view of (67) and (63) implies Therefore, from (58) with the condition (57), and further from (53), a solution of the system (51) is const 1 and const 2 are constants determined by the initial conditions (x(0), z(0)) := (x 0 , z 0 ). This peakon-like solution was also presented in the paper [31]. The graph of the parametric curve (72) is drawn in the Figure 2. Calculating the derivatives of x(t) and z(t) with respect to t, we get Hence, for the solution (72) a stagnation point in the fluid, where x ′ (t) = c, z ′ (t) = 0, appear only for t → ±∞. We observe that at these points the path of the particle has a horizontal tangent.
The other solutions of the equation (69) satisfy The homogeneous equation: where θ is an integration constant. By the method of variation of constants, the general solution of the non-homogeneous equation (74) is given by where θ(w) is a continuous function which satisfies the equation The solution of the equation (78) is β being a constant. Therefore, the solution of the non-homogeneous equation (74) has the expression Taking into account (67), we get We separate the variables in (81): From (58), (61), we have Thus, (82) can be written as By (55) we obtain Further, from (53), we get another solution of the system (51):

Z(t) being the solution of the equation (85).
We observe that the solutions (87) are not closed curves.
Indeed, if there exists t 2 > t 1 such that Z(t 2 ) = Z(t 1 ), then, in view of (85), we . Thus, although in the moving frame we obtain in this case a closed curve with Z(t 1 ) = Z(t 2 ) and X(t 1 ) = X(t 2 ), in the fixed frame we get z(t 2 ) = z(t 1 ) and x(t 2 ) − x(t 1 ) = c(t 2 − t 1 ) = 0. If c > 0, the particles which follow these curves will have a forward drift, if c < 0, they will have a backward drift.
Let us now investigate more the equation (85). Using the formula: sinh 2 (x) = cosh(2x)− 1 2 , the equation (85) can be written in the form: Taking into account the expression of cosh(x) as Taylor series: we get under the square root in (88) the following power series The constant A is different from zero (52), thus, the power series (90) contains for sure powers of Z higher than four. A partial sum of the above series is a polynomial of degree higher than four. Thus, considering only a partial sum of this series, the solution of the equation (85) involves a hyperelliptic integral (for hyperelliptic integrals see, for example, [2], page 252). Its inversion would lead to a hyperelliptic function.
There are special cases when a hyperelliptic integral can be reduced to an elliptic one and thus, its inversion will contain the Jacobi elliptic functions sn, cn, sc, etc. If in (90) we have C = c (by choosing appropriate constants c 0 , c 0 in (49), respectively (50)) and we consider powers of Z till six, the solution of the equation (85) involves the following hyperelliptic integral ± dZ We consider the substitution and thus, the left-hand side in (91) reduces to an elliptic integral of the first kind: This elliptic integral of the first kind may by reduced to the Legendre normal form.
Case 1: all the zeroes of the cubic polynomial under the square root in (93) are real and distinct. We denote them byẐ 1 <Ẑ 2 <Ẑ 3 . Because the leading coefficient of this cubic polynomial is smaller then zero and its constant term is greater then zero, we have either Case 1a: the condition (94) is fulfilled. Then we introduce the variable ϕ by (see [38] Ch. VI, §4, page 602) and we get where the constant 0 < k 2 1 < 1 is given by Therefore we obtain the Legendre normal form of the integral in (93): the constant factor in front of the integral being equal to The inverse of the integral in (98) is sn (the Jacobi elliptic function sine amplitude, see, for example, [2]) sn (C 1 t; k 1 ) = sin ϕ In view of the notations (92) and (96) we get cn being the Jacobi elliptic function cosine amplitude (see, for example, [2]). Taking into account the expressions for the derivatives of sn and cn (see, for example, [2]), that is, where dn (t; k) := 1 − k 2 sn 2 (t; k), we obtain We introduce (101) and (102) in (87) and we get x(t) and z(t) explicitly. We remark that, if k 2 A 2 3 − B 2 4 < 0 and k 2 A 2 + 2βB > 0 (this can happen, for example, for a small enough A and for (B > 0 & β > 0) or (B < 0 & β < 0)), the coefficients of the cubic polynomial in (93) have alternating signs. Thus, by Descartes' rule of signs, if all the roots are real, the situation (94) occurs.
Case 1b: the condition (95) is fulfilled. Z 1 ,Ẑ 2 ,Ẑ 3 being the zeroes of the real cubic polynomial under the square root in (93), this polynomial has the unique decomposition −4β 2 (Ẑ −Ẑ 1 )(Ẑ −Ẑ 2 )(Ẑ − Z 3 ). A suitable change of variable transforms the integral (93) to the Legendre normal form (98) up to a constant. In general, the elliptic functions can have complex arguments (for example, if the constant factor in front of the integral (98) is a complex number, then the obtained sine amplitude function will depend on a complex variable) but here we are interested only in the real case. We are also looking for a real Z, so,Ẑ introduced by (92) has to be greater then zero. With the change of variable (96), which brings the integral (93) to the Legendre normal form (98), because nowẐ 2 < 0, we end up with aẐ which can be positive, negative or zero. Thus, in this case, we get the expression (101) of Z(t) only if t satisfieŝ Z 2 sn 2 (C 1 t; k 1 ) +Ẑ 3 cn 2 (C 1 t; k 1 ) > 0 (103) that is, by sn 2 (C 1 t; k 1 ) + cn 2 (C 1 t; k 1 ) = 1, only if t satisfies where 0 <Ẑ 3 Z3−Ẑ2 < 1. For a very small positive solutionẐ 3 → 0, the set of t's which fulfill the above inequality (104) tends to the empty set. Thus, in this case, the hyperelliptic integral in (91) can not be reduced to an elliptic one and the solution can not be expressed with the aid of the Jacobi elliptic functions. The solution will be expressed with the aid of a hyperelliptic function obtained by the inversion of the integral in (91).
Case 2: the cubic polynomial under the square root in (93) has only one real solution denotedẐ 0 . Because the leading coefficient of this cubic polynomial is smaller then zero and its constant term is greater then zero, we havê We denote by p and q the real coefficients such that We introduce the variable ψ by (see [38] Ch. VI, §4, page 602) and we get where the constant 0 < k 2 2 < 1 is given by Therefore we obtain the Legendre normal form of the integral in (93): the constant factor in front of the integral being equal to The inverse of the integral in (109) is Taking into account (107), we get If t satisfies the following inequality < 1, then,Ẑ from (112) is greater than zero and we obtain by (92) the expression of Z(t). We can calculate the time derivative of Z(t) and by (87) we get x(t) and z(t) explicitly.
As in the case 1b, for very small positive solutionẐ 0 → 0, the set of t's which fulfill the inequality (113) tends to the empty set. Thus, in this case, the hyperelliptic integral in (91) can not be reduced to an elliptic one and the solution can not be expressed with the aid of the Jacobi elliptic functions. The solution will be expressed with the aid of a hyperelliptic function obtained by the inversion of the integral in (91).
We remark that, if k 2 A 2 3 − B 2 4 < 0 and k 2 A 2 + 2βB < 0 (this can happen, for example, for a small enough A and for (B < 0 & β > 0) or (B > 0 & β < 0)), the coefficients of the cubic polynomial in (93) have signs ---+. Thus, by Descartes' rule of signs, only one root is positive and we are in the case 1b or in the case 2. This positive root it will be close to zero.
Let us draw below some of the curves obtained for different values of the parameters, using Mathematica 3 . We consider k = 1, h 0 = 1, g = 9.8, ǫ = 0.1, α = 0, β = 1 and ω 0 = 2 > 0. Then, by (50), choosing the square root with sign plus, we get c = 4.07454 > 0, A = 0.176526 > 0, B = 2 > 0. We take c 0 such that C = c. In this case, all the roots of the cubic polynomial under the square root in (93) are real and we get Z(t) in the form (101). The graph of the curve obtained is drawn in Figure 3. Using the same values for k, h 0 , g, ǫ, α and β but taking ω 0 = 20 > 0, we get by (50), choosing the square root with sign minus, c = 4.29294 > 0, A = −1.33654 < 0, B = 20 > 0. We take c 0 such that C = c. Then all the roots of the cubic polynomial under the square root in (93) are real and we get Z(t) in the form (101). The graph of the curve obtained is depicted in Figure 4. For the same k, h 0 , g, ǫ, α and β as above, with ω 0 = 2 > 0 but choosing the square root with sign minus in (50), we get c = −1.59773 < 0, A = −0.306137 < 0, B = 2 > 0. We take c 0 such that C = c. Then all the roots of the cubic polynomial under the square root in (93) are real and we get Z(t) in the form (101). The graph of the curve obtained is presented in Figure 5. We choose now ω 0 = −20 < 0, k, h 0 , g, ǫ α and β having the same values as above. We get by (50), choosing the square root with sign plus, c = −4.29294 < 0, A = 1.33654 > 0, B = −20 < 0. We take c 0 such that C = c. Then the cubic polynomial under the square root in (93) has only one real root, that is,Ẑ 0 = 0.000798. The real coefficients p and q from (106) have the values: p = 9.55422, q = 24.8588 and the right hand side in (113) has the value 0.99968. We conclude that in this case the solution can not be expressed through Jacobi elliptic functions.
We would like to make some remarks on the stagnation points inside the fluid. Calculating the derivatives with respect to t of x(t) and z(t) from (87), we get where, taking into account (85),  .