A Numerical Study of Linear Long Water Waves over Variable Topographies using a Conformal Mapping

. In this work we present a numerical study of surface water waves over variable topographies for the linear Euler equations based on a conformal mapping and Fourier transform. We show that in the shallow-water limit the Jacobian of the conformal mapping brings all the topographic effects from the bottom to the free surface. Implementation of the numerical method is illustrated by a MATLAB program. The numerical results are validated by comparing them with exact solutions when the bottom topography is flat, and with theoretical results for an uneven topography.


INTRODUCTION
Over the centuries, the curiosity has led the human being to explore and investigate the ocean in several perspectives. Nonetheless, due its immensity and complexity the ocean remains as a strange place for us -what make it an enjoyable subject of studying. Among the topics of research in ocean sciences, we can cite the water waves propagation. A better understanding of the dynamic of water waves is crucial, for instance, to protect coastal communities from ocean hazards such as tsunamis.
Historically, Isaac Newton (1687) was the pioneer in the study of water waves. However, his results lack of mathematical rigorous. Approximately one century later, Leonhard Euler (1757) deduced the equations that model the dynamic of water waves which bears his name [7]. The Euler equations is one of the main models used to investigate the dynamic of surface gravity waves. For instance, these equations describe the motion of a flow of water over rocks, ship wakes, which are waves generated by the passage of a moving loading in the surface of the water [2], and waves generated by storms [13]. One of the ways of solving the Euler equations broadly used in the literature is to construct a conformal mapping that maps the fluid domain onto a strip [8]. This approach transforms the free boundary problem given by the Euler equations into a family of ordinary differential equations which are easier to be solved numerically. Since then, many works have been done on this topic [4,5,8,9,10,12,16]. It is hard to give a comprehensive overview of contributions. For the interested reader, we mention a few articles and references therein for further details. Solitary waves propagating over a flat bottom were computed numerically by Choi and Camassa [5] and later in the presence of a vertically sheared current with constant vorticity by Choi [4]. In this scenario, particle trajectories and the pressure within the bulk of the fluid were computed in Ribeiro-Jr et al [20]. In the presence of an uneven topography Nachbin [16] constructed a conformal mapping to flatten the bottom topography and asymptotically obtained a Boussinesq-type equation with variable coefficients. More recently, particle trajectories for the linear Euler equations were investigated by Flamarion et al. [10] in the presence of a linear sheared current and a variable topography.
Conformal mapping has also been applied in the study of free boundary problems in small scales such as the walking of a fluid droplet (silicon) on the surface of a vibrating bath. A problem which mimics the wave-particle duality present in quantum mechanics [3,6,15]. Using the Schwarz-Christoffel mapping, Nachbin and collaborators [17,18,19] have investigated in details the role of the free surface wave as a spontaneous mechanism of synchronisation of oscillating droplets confined in different types of cavities. Their findings have contributed to a better understand of the hydrodynamic pilot-wave analogy with quantum mechanics and other phenomena as Kuramoto-like synchronisation.
In this work, we study numerically the water wave interactions with a variable topography using the linear Euler equations. The topography models variations and geometric irregularities on the seabed, such variations can be caused by submarine relief, submerged equipments and effects of human activities. The problem is formulated through the use of a conformal mapping which flattens the fluid domain onto a strip. We show that when the average depth of the channel is small compared to the wavelength of the surface wave, the Jacobian of the conformal mapping carries all the underwater obstacle information, then the linear Euler equations can be written in a much simpler system of coordinates. The numerical method is implemented on MATLAB and its validation is performed by comparing numerical solutions with exact ones when the bottom is flat, and with theoretical results for uneven topographies. A similar result on the topic of this article is presented by Artiles and Nachbin [1]. Although we study the same set of equations of these authors there are remarkable differences between the two studies. First, the mathematical formulation considered by Artiles and Nachbin [1] depends on the computation of an explicit formula of the Dirichlet-to-Neumann operator, whereas in this work this requirement is removed. Second, Artiles and Nachbin [1] considered an implicit scheme to solve the initial value problem, while here we use the classical Runge-Kutta forth-order method. In addition, our formulation allows to conclude easily that in the shallow water limit the mapping Jacobian's along the free surface carries the topographic effects.
The study and the numerical code (see Appendix A) presented in this article can be simply extended to investigate other problems in linear water waves. For instance: (i) in the attempt of capture waves that remain stationary over topographic obstacles, which have been known as trapped waves; (ii) to computed the velocity field beneath the free surface and consequently the path of the fluid particles; (iii) to analyse the pattern of waves generated by the passage of a moving disturbance.
This article is organized as follows. In section 1 we present the mathematical formulation of the Euler equations. In section 2 we describe the conformal mapping technique and rewrite the Euler equations in the canonical domain, which is a uniform strip, and present a numerical method to solve them. In section 3 we present the numerical results and the conclusion in section 4.

MATHEMATICAL FORMULATION
We consider a two-dimensional irrotational flow of an inviscid and incompressible fluid with constant density (ρ) in a finite depth channel in the presence of a variable topography (h(x)). On the top of the channel we have a free surface wave ζ (x,t), which is under the action of gravity (g). The governing Euler equations written in terms of the free surface (ζ (x,t)) and the potential velocity (φ (x, y,t)) are [22] To deal mathematically with this problem, it is convenient to consider dimensionless variables.
Using the typical wavelength λ as the horizontal length, h 0 for the vertical length, a for the wave amplitude, agλ / √ gh 0 as the velocity potential scale, λ / √ gh 0 as the time scale as done in [22], we obtain the following dimensionless equations where ε = a/h 0 is the nonlinearity parameter, and µ = h 0 /λ is the shallow-water/long-wave parameter. For computational purposes, it is convenient to consider the new scaling y → µy. In this new scaling, the potential velocity is a harmonic function, thus the Euler equations becomẽ We are interested in investigating solutions of (2.1) in the linear regime (ε = 0) and weakly dispersive regime (µ ≈ 0). Under these assumptions we have the linear system of equations When the bottom is flat, the phase speed of the linear waves are given by [22] c(k) = tanh(µk) We recall that this is for example the speed of the crest of the wave as it propagates. Notice that in the limit µ → 0, all linear waves travel with the same speed c = 1.
In the following section, we present a numerical method to solve (2.2).

CONFORMAL MAPPING AND NUMERICAL METHODS
We recall that one of the main features of conformal mappings is that the Laplace equation is conformally invariant. Therefore, to solve (2.2) we can apply a conformal mapping and solve it in a simpler domain. To this end, we construct a conformal mapping Using the fact that y is a harmonic function and the boundary conditions (3.1) we obtain where Fourier modes are given by and k j = (π/L) j, j ∈ Z. The Cauchy-Riemann equation x ξ = y η yields More details of this conformal mapping can be found in [11].
We choose the vertical strip to be D = µ. The motivation for this choice of D lies on the fact that at η = 0 and in the limit µ → 0 we formally have (3. 2) The mapping's Jacobian is given by when evaluated along the free surface is denoted as Hence, in the shallow-water (µ ≈ 0) limit we obtain from equation (3.2) Now, our goal is to reduce the system (2.2) into a free surface problem with no boundary conditions in the canonical domain. In order to do that, let φ (ξ , η,t) =φ (x(ξ , η), y(ξ , η),t) be the potential velocity in the canonical domain. Notice that the chain rule and the Cauchy-Riemann equations (x ξ = y η and x η = −y ξ ) yields Using the Neumann condition (2.2) 2 and the boundary conditions (3.1) at η = −D in (3.3) we obtain Since φ is a harmonic function and satisfies (3.4), we have the Neumann problem for φ in a strip (3.5) The solution of this problem can be written in terms of Fourier modes and yields On the other hand, at η = 0 from (3.3) we have Substituing in (2.2) we obtain the following set of equation at η = 0 Notice that in the shallow-water limit tanh(k j µ) ≈ k j µ and from (3.6) we have φ η (ξ , 0,t) = −µφ ξ ξ (ξ , 0,t).
Consequently, from (3.7) 1 we obtain that the dynamic of the free surface wave is governed by the classical second order wave equation which implies that the right-going wave solution travels with speed We solve (3.6) and (3.7) in a computational domain [−L, L], with a uniform grid with N points and step ∆ξ = 2L/N. The spatial derivatives are computed using the Fast Fourier Transform (FFT) [21]. In addition, the time advance of the system (3.7) is computed through the Runge-Kutta fourth order method (RK4) with time step ∆t. The numerical method is implemented in MATLAB and the codes can be found in the Appendix A.

NUMERICAL RESULTS
In this section we consider different types of topographies and analyse the interaction of a rightgoing pulse with the topography. To this end, we set the initial wave to be where σ is the effective width of the wave and ξ 0 is its initial location. The initial potential velocity is taken as where c(k) is given by (2.3). This choice of initial data guarantees that the solution travels to the right when the bottom is flat.

Flat bottom
When the bottom is flat, exact solutions of (3.7) can be written in terms of Fourier modes [22] as So, in order to valid our numerical method we compare numerical solutions of (3.7) with formula (4.1). It is worth noting that dispersive waves travel with speed less than c = 1. Figure 1 displays the comparison between the numerical and exact solutions of (3.7) for different values of µ.
Numerically, we see that the dispersion term not only changes the wave speed, but also affects the amplitude of the wave front attenuating it. In addition, when the parameter µ → 0, we see a the Gaussian pulse traveling maintaining its form, no dispersive effects are observed. The relative error in all these simulation is at least order O(10 −7 ).

Bragg Resonance
The Bragg Resonance describes the interaction between a wave and a periodic series of underwater obstacles [14]. This phenomenon is featured by having a reflected wave with wavelength twice as large as the wavelength of the sequence of obstacles underwater. It occurs when variations on the topography are comparable with the free surface wavelength.
In order to verify that the model where λ b is the wavelength of each obstacle underwater. Figure 2 depicts the wave profile of a Gaussian pulse with wavelength λ ≈ 6 after interacting with the obstacle underwater. Each obstacle has wavelength λ b = 2, and the wavelength of a typical reflected waves is roughly λ R ≈ 4. When setting λ b = 1, Figure 3 displays only reflected waves which has wavelength about λ R ≈ 2. These results confirm the ones predicted by the theory. When the wavelength of the initial wave and the obstacle are not comparable, the Bragg Resonance does not occur. To illustrate this, let λ b = 0.5, in this case we observe that only a wave of small amplitude is reflected and its wavelength is around λ R ≈ 4 (see Figure 4). The Gaussian pulse passes over the obstacle and barely feels it.

Waves propagating over a heterogenous medium
When a wave travels with constant speed c 1 and c 2 in two different media, the coefficients of transmission and reflection of these waves are respectively [14] T =   In order to verify that our numerical method captures this phenomenon, we consider the bottom topography where the shallower region is the medium 1 and the deeper region is the medium 2. Once we fix the depth of the channel in the shallow-water limit the speeds c 1 and c 2 can be obtained by the expression (3.8). In this configuration we have c 1 = 1.000000 and c 2 = 0.816496. Thus, the coefficients of transmission and reflection are T = 0.898979 and R = −0.053197 respectively. In order to verify that the numerical method can compute these coefficients correctly, it is enough to compute the numerical speeds of the transmitted and reflected wave. To this end, we compute solutions ζ (ξ , 20) and ζ (ξ , 40) as depicted in Figure 5, and calculate the mean speed of the transmitted and reflected wave using the ξ −position of their local minima and maxima. The numerical speeds of the transmitted and reflected wave are respectively

CONCLUSIONS
In this article, we have presented a numerical method to study linear water waves over an uneven topography using a conformal mapping. We showed that in the shallow-water limit, the Jacobian of the conformal mapping carries all the topography information. This allowed us to approximate the topography by the Jacobian of the conformal mapping in a manner. Consequently, the problem was reformulated in a much simpler system of equations where numerical computations were performed. The results were validated by comparing numerical solutions with exact ones when available and when not, with theoretical results.

A NUMERICAL CODE
In this appendix we present the MATLAB program used to generate Figure 3. The other figures presented in this article can also be reproduced by this code.