Advances in Mathematical Physics
Volume 2018 (2018), Article ID 3469534, 10 pages
https://doi.org/10.1155/2018/3469534

Research Article

CIP Method of Characteristics for the Solution of Tide Wave Equations

Yafei Nie,¹ Kai Fu,² and Xianqing Lv¹

¹Physical Oceanography Laboratory, Ocean University of China, Qingdao 266100, China
²School of Mathematical Sciences, Ocean University of China, Qingdao 266100, China

Correspondence should be addressed to Kai Fu; kfu@ouc.edu.cn

Received 2 February 2018; Accepted 26 May 2018; Published 2 July 2018

Academic Editor: Zhi-Yuan Sun

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

Abstract.

The CIP-MOC (Constrained Interpolation Profile/Method of Characteristics) is proposed to solve the tide wave equations with large time step size. The bottom topography and bottom friction, which are very important factor for the tidal wave model, are included to the equation of Riemann invariants as the source term. Numerical experiments demonstrate the good performance of the scheme. Compared to traditional semi-implicit (SI) finite difference scheme which is widely used in tidal wave simulation, CIP-MOC has better stability in simulating large gradient water surface change and has the ability to use much longer time step size under the premise of maintaining accuracy. Besides, numerical tests with reflective boundary conditions are carried out by CIP-MOC with large time step size and good results are obtained.

1. Introduction

Tide wave is the water movement caused by the gravitational force of the celestial body. The control equations of tidal wave model are the original shallow water equations including source terms (e.g., bottom friction, bottom topography and horizontal eddy viscosity term) [1]. Therefore, getting the solution of shallow water equations with source terms is of great significance.

In the simulation of ocean tidal waves, Eulerian schemes are widely used, for example, Backhaus [2] and Casulli [3] used semi-implicit scheme (hereafter SI) for the solution of shallow water equations; Lv and Zhang [4] used the semi-implicit scheme to solve tide wave equations and their computational format was used to study bottom friction coefficients [5] and tidal open boundary conditions [6]; the finite volume method is also frequently used [7, 8]. However, since the phase speed of the gravity waves is so fast, time step is strictly limited by CFL condition when using Eulerian schemes [9].

There are growing interests in finding methods that enable using large time steps to solve shallow water equations. Erbes [9] proposed a semi-Lagrangian method of characteristics (MOC) with quadratic interpolation to solve the shallow water equations with a large CFL number, but there was dispersion error produced. Ogata and Yabe [10] applied the CIP method to shallow waters combined with the MOC and showed low dispersion error and low numerical damping even with large CFL number. Toda et al. [11] proposed a new scheme by adopting the Conservative Semi-Lagrangian (CSL) based on the CIP method which shows good conservation. However, applying this method in realistic simulation of ocean tidal wave will meet a difficulty treating source terms including bottom topography, bottom friction, and boundary conditions.

This paper is designated to apply the CIP-MOC [9] to the solution of tide wave equations which are based on shallow water equations and include source terms of bottom topography and friction. The ideal numerical experiment shows that the CIP-MOC method outperforms SI scheme in stability and dispersion errors when using large time step. Problem with reflection boundary condition is also discussed.

This paper is organized as follows. Section 2 presents the CIP-MOC method. In Section 3, numerical experiments are carried out and comparison is made between the CIP and SI. Conclusion is given in Section 4.

2. Model and Method

2.1. CIP Method

CIP method is a compact, robust, less diffusive, and high-order scheme in computational fluid dynamics [12]. Compared with other traditional interpolation methods, the internal information of grid is better described by CIP which uses spatial derivative in the description of function distribution in grid cells. We consider the following advection problem:

the solution is [13]

If and are known at each grid point, the profile between two adjacent points can be approximated by cubic polynomial

where . Consider four constraints

then coefficients can be given as

where , , and .

Thus, the profile of and at next time step can be obtained as

where .

2.2. One-Dimensional Tide Wave Equations

Let bottom friction term be the linear form, then the one-dimensional tide wave equations in a primitive form are written as [14]

where is time, is Cartesian coordinate, is undisturbed water depth, is sea surface elevation above the undisturbed sea level, is the flow velocity, is the bottom friction coefficient, and is the gravitational acceleration.

To solve (9) with the method of characteristics and preserve the balance between the depth gradient and the bottom effect, the surface gradient method [15] is used here. Replace the variable with the water level , where represents the bottom topographic height (see Figure 1), we get the equivalent of (9) which can be written as

Figure 1: Shallow water variables.

In a vector-matrix form, (10) can be rewritten as

The eigenvalue matrix and eigenvector matrix of coefficient matrix are

where and . Thus, matrix can be diagonalized as , where

is the inverse of . Then (11) can be written as

where represents the right hand side of (11).

Two characteristics about Riemann invariants can be obtained from (14) as

By using CIP method, we have

where and are calculated by

the subscript CIP in variables means that they are interpolated with the CIP scheme, and the bottom friction term is discretized with the Crank-Nicholson method which is implicit in time. and at time step can be obtained by solving linear equations (16) (see Figure 2):

Figure 2: Space-time diagram that shows the characteristics . The Riemann invariants consist of and , are characteristic speeds.

Time integration of the bottom topographic term in (18) and (19) can be approximated using as follows [16]:

When applying the CIP method, spatial gradient at the next time step is also needed. Take spatial derivative of (15)

includes the terms related to spatial derivatives of , , and . Considering that the influence of is not significant in this problem, we neglect here. Then (21) can be solved as

then values and derivatives of each value at time step are obtained.

2.3. Two-Dimensional Tide Wave Equations

Assuming that pressure is hydrostatic and density is constant, the depth averaged two-dimensional tidal model without horizontal eddy viscosity is given as follows [17]:

where and are and components of flow velocity and represents the Coriolis parameter.

In the same way as one dimension, replace variable with water level ; we get the equivalent of (24) which can be written in vector-matrix form as

where

Since matrices and are not commutative; thus, there is no matrix to diagonalize them at the same time. Split (25) into two sequential phases; the source term is divided into two terms added to each directional phase:

where

In order to get high accuracy in time, we calculate alternately in the x- and y-directions as follows:

where and represent operation of (27) and (28), respectively.

Here we give the solution process in x-direction. Similar to the one-dimensional case, left multiply (27) by the eigenvector matrix :

where

Equation (31) leads to the following three equations:

where are the Riemann invariants. These three equations all have the form of convective equations that can be calculated by CIP method. The Coriolis force, bottom topography, and bottom friction term are added to each Riemann invariant along characteristic line. Thus, the Riemann invariants and can be discretized as follows:

So we have

where , , ,

Similar to the derivation in one dimension, neglecting the influence of , the derivatives of each value at time step are obtained

The operation in the x-direction is completed.

The operation in y-direction is essentially the same as the x-direction and the y derivative is computed by linear interpolation

where .

3. Numerical Results

3.1. A Small Perturbation of One-Dimensional Steady State Water

To demonstrate that the CIP method can be applied to the computation of ocean tidal waves with bottom topography and bottom friction, we first consider the quasi-stationary test case given in [18]. The difference is that we add the effect of the bottom friction here.

The bottom topography is given by

The initial conditions are

where is constant amplitude of the perturbation. According to the propagation mechanism of wave, the small disturbance will be split into two smaller waves propagating to both sides. However, because of the Gibbs phenomenon, the small pulse of perturbation is a challenge for many numerical schemes [19]. Numerical experiments are carried out with both the CIP and SI schemes; the surface level with cells at the final time is shown in Figure 3. We keep as the SI scheme is strictly restricted by the CFL condition, and when using the CIP.

(a)

(b)

(c)

(d)

Figure 3: A small perturbation of one-dimensional steady state water with a pulse . (a) shows SI scheme with CFL=0.63; (b), (c), and (d) are the results of the CIP method with CFL=0.63, 1.57, and 2.61.

Since flow velocity and are far smaller than wave velocity , the CFL number is defined as

The experimental results are similar to those in [19, 20]. It can be seen that nonphysical oscillations occur when simulating marching wave using the SI scheme, waves are broken in the process of travel, and the amplitude of nonphysical oscillations increases gradually in time, which also means that the scheme does not have the ability to simulate large gradient surface change. However, the shape of the simulated wave is smooth and has no oscillations when using CIP until . Figure 4 gives the results computed by CIP with CFL=2.61 for k=0, 1, and 2. It can be seen that the separated wave height is less than half of the initial perturbation because of the bottom friction. Numerical experiments show that the CIP method can successfully simulate large gradient water surface with both bottom topography and bottom friction using large time step sizes.

Figure 4: Effect of bottom friction. Results computed by CIP with CFL=2.61 for k=0, 1, and 2.

3.2. Smooth Surface Wave Propagation Problem

Since boundary condition plays a significant role in the shallow water flow of an actual marine area, we consider a smooth surface wave propagation problem here with reflective boundary conditions. The conditions for the x-direction are

Since the spatial gradient is also required when using the CIP method, the spatial gradient condition at the boundary is obtained according to (48):

When the traced upstream departure point on the shore, the corresponding actual upstream departure point in the waters is approximated using reflective boundary conditions (see Figure 5).

Figure 5: Space-time diagram that shows the reflecting boundary conditions. The filled circle represents the approximated upstream departure point, while the hollow circle represents the virtual upstream point on the shore.

The initial conditions are given as follows:

where . Let the bottom friction coefficient and the constant Coriolis parameter be .

The computation domain is , and and are set to 1.0. The contour map of the water surface at t=20 and t=40 is shown in Figures 6 and 7, respectively. It can be seen that at t=20 the contour line is not a complete circle because of the impact of the bottom topography (similar to the one-dimensional case). At , the reflected waves are superimposed on each other and eventually form what is shown in Figure 7. The results show that when time step increases to three times to the SI method, CIP-MOC scheme still simulate the change process of water surface very well, even considering the bottom topography (see Figure 8), boundary, and Coriolis effect.

(a)

(b)

(c)

(d)

Figure 6: Smooth surface wave propagation problem. The contour of the water level H at t=20. (a) is the result of SI with CFL=0.63; (b), (c), and (d) are the results of CIP-MOC with CFL=0.63, 1.42, and 1.90.

(a)

(b)

(c)

(d)

Figure 7: Smooth surface wave propagation problem. The contour of the water level H at t=40. (a) is the result of SI with CFL=0.63; (b), (c), and (d) are the results of CIP-MOC with CFL=0.63, 1.42, and 1.90.

(a)

(b)

Figure 8: The comparison of cross-sections between topographic (a) and nontopographic conditions (b) with different CFL numbers.

4. Conclusion

In this paper, we apply the CIP-MOC to the solution of tide wave equations using large time step size. Bottom topography, bottom friction, and Coriolis term are included to the equation of Riemann invariants as the source term. One-dimensional experiment shows that the CIP-MOC scheme can keep the shape very well and has no oscillations even when the time step is increased to four times of the SI scheme in the simulation of large gradient water surface. Bottom friction effect is also tested in the experiment. In two-dimensional case, the problem with reflected boundary conditions is considered. Numerical experiments show that the CIP-MOC scheme is succeeded in simulating the change process of water surface using time step that is three times of the SI method. All the results validated the good stability and low dispersion error of the CIP-MOC method.

Data Availability

All results presented in the article were produced from model simulations. Therefore, there is no data to be made available. Researchers who wish to replicate the study will use the equations and parameters described in the article. With such equations and parameters, researchers can use modeling simulations to replicate the figures presented in the article.

Conflicts of Interest

The authors declare that there are no conflicts of interest regarding the publication of this paper.

Acknowledgments

Partial support for this research was provided by the National Key Research and Development Plan through Grant 2016YFC1402304, the Key Research and Development Plan of Shandong Province through Grants 2016ZDJS09A02 and ZR2016AB16, and the National Natural Science Foundation of China through Grant 41606006.

References

J. Pedlosky, Geophysical fluid dynamics, Springer Science & Business Media, 2013.

J. O. Backhaus, “A semi-implicit scheme for the shallow water equations for application to shelf sea modelling,” Continental Shelf Research, vol. 2, no. 4, pp. 243–254, 1982.

V. Casulli, “Semi-implicit finite difference methods for the two-dimensional shallow water equations,” Journal of Computational Physics, vol. 86, no. 1, pp. 56–74, 1990.

X. Lu and J. Zhang, “Numerical study on spatially varying bottom friction coefficient of a 2D tidal model with adjoint method,” Continental Shelf Research, vol. 26, no. 16, pp. 1905–1923, 2006.

H. Pan, Z. Guo, and X. Lv, “Inversion of tidal open boundary conditions of the M2 constituent in the bohai and yellow seas,” Journal of Atmospheric and Oceanic Technology, vol. 34, no. 8, pp. 1661–1672, 2017.

Z. Guo, H. Pan, W. Fan, and X. Lv, “Application of surface spline interpolation in inversion of bottom friction coefficients,” Journal of Atmospheric and Oceanic Technology, vol. 34, no. 9, pp. 2021–2028, 2017.

R. H. Weisberg and L. Zheng, “A finite volume coastal ocean model simulation of the tide, buoyancy, and wind-driven circulation of Tampa Bay,” Journal of Geophysical Research, vol. 111, article C01005, 2006.

C. Chen, H. Huang, R. C. Beardsley, H. Liu, Q. Xu, and G. Cowles, “A finite volume numerical approach for coastal ocean circulation studies: comparisons with finite difference models,” Journal of Geophysical Research: Oceans, vol. 112, no. 3, Article ID C03018, 2007.

G. Erbes, “A Semi-Lagrangian Method of Characteristics for the Shallow-Water Equations,” Monthly Weather Review, vol. 121, no. 12, pp. 3443–3452, 1993.

Y. Ogata and T. Yabe, “Multi-dimensional semi-lagrangian characteristic approach to the shallow water equations by the CIP method,” International Journal of Computational Engineering Science, vol. 5, no. 3, pp. 699–730, 2004.

K. Toda, Y. Ogata, and T. Yabe, “Multi-dimensional conservative semi-Lagrangian method of characteristics CIP for the shallow water equations,” Journal of Computational Physics, vol. 228, no. 13, pp. 4917–4944, 2009.

T. Yabe and T. Aoki, “A universal solver for hyperbolic equations by cubic-polynomial interpolation. I. One-dimensional solver,” Computer Physics Communications, vol. 66, no. 2-3, pp. 219–232, 1991.

T. Yabe and E. Takei, “A new higher-order Godunov method for general hyperbolic equations,” Journal of the Physical Society of Japan, vol. 57, no. 8, pp. 2598–2601, 1988.

Y. Anle and L. Fengqi, Physical Oceanography, China Ocean University Press, 1992.

J. G. Zhou, D. M. Causon, C. G. Mingham, and D. M. Ingram, “The surface gradient method for the treatment of source terms in the shallow-water equations,” Journal of Computational Physics, vol. 168, no. 1, pp. 1–25, 2001.

R. Akoh, S. Ii, and F. Xiao, “A CIP/multi-moment finite volume method for shallow water equations with source terms,” International Journal for Numerical Methods in Fluids, vol. 56, no. 12, pp. 2245–2270, 2008.

R. W. Lardner, “Optimal assimilation of current and surface elevation data in a two-dimensional numerical tidal model,” Applied Mathematical Modelling, vol. 19, no. 1, pp. 30–38, 1995.

R. J. LeVeque, “Balancing source terms and flux gradients in high-resolution Godunov methods: the quasi-steady wave-propagation algorithm,” Journal of Computational Physics, vol. 146, no. 1, pp. 346–365, 1998.

Q. Zhu, Z. Gao, W. S. Don, and X. Lv, “Well-balanced hybrid compact-WENO scheme for shallow water equations,” Applied Numerical Mathematics, vol. 112, pp. 65–78, 2017.

Y. Xing and C.-W. Shu, “High order finite difference WENO schemes with the exact conservation property for the shallow water equations,” Journal of Computational Physics, vol. 208, no. 1, pp. 206–227, 2005.