_{1}

^{*}

Finite difference techniques are widely used for the numerical simulation of time-dependent partial differential equations. In order to get better accuracy at low computational cost, researchers have attempted to develop higher order methods by improving other lower order methods. However, these types of methods usually suffer from a high degree of numerical dispersion. In this paper, we review three higher order finite difference methods, higher order compact (HOC), compact Padé based (CPD) and non-compact Padé based (NCPD) schemes for the acoustic wave equation. We present the stability analysis of the three schemes and derive dispersion characteristics for these schemes. The effects of Courant Friedrichs Lewy (CFL) number, propagation angle and number of cells per wavelength on dispersion are studied.

Real life situations in many disciplines including engineering, physics, economics, biosciences, etc, can be described through mathematical models. Differential equations play an important role in modelling interactions in physics, engineering, and biological processes, from ground motion (earthquake) to interactions between neurons. Differential equations can either be categorized as ordinary differential equations (ODEs), which are for example, used to model dynamical systems and their applications to science and engineering or partial differential equations (PDEs), which most importantly used to describe a wide variety of phenomena in physics; such as sound, heat, fluid dynamics, elasticity and in most field of engineering in general. However, most differential equations, such used to solve real-life problems, have no analytical solutions or are unrealistic to solve analytically. Therefore, solutions for these problems are attempted by means of numerical approximations. To this end various numerical methods are developed to provide better approximate solutions to these equations.

Numerical techniques are a powerful tool for handling large systems of equations involving complex geometries that are otherwise impossible to handle analytically [

The aim of this research is to explore a higher order compact finite difference scheme with comparably low dispersive characteristics. Therefore, in this paper, firstly, we discussed a conventional higher order finite difference schemes. Then we discussed another higher order finite difference scheme developed by Das [

The two dimensional version of the wave equation with velocity v and acoustic pressure u in homogeneous media can be written as

where

and the boundary condition

To discretize the above PDE we consider a uniform grid

The central finite difference operators for second derivatives are written as [

A higher order finite difference approximation for the second derivative is given by Liu [

where

Though this approximation is accurate up to order 2M, its not efficient because when M becomes large, the stencil that need to be computed becomes large.

The standard second-order finite difference approximation for the acoustic wave equation in (1) using the above operators is

where

A fourth order accurate implicit finite difference scheme for one dimensional wave equation is presented by Smith [

Consider two dimensional wave equation, using Taylor’s series expansion of

If u is a solution of (1), then we have the following

and

Substituting equations (13) and (14) in to equation (12) and using the notations

and

M | |||||
---|---|---|---|---|---|

1 | −2 | 1 | - | - | - |

2 | - | - | |||

3 | - | ||||

4 |

we obtain

Following Strikwerda [

and

We also have

an

Substituting (16)-(19) into (15) gives a scheme which is fourth order accurate

where

Applying the operator

to both sides of (20), we obtain

Equation (21) can be written as an ADI-scheme and can be referred to as HOC-ADI, which can be presented more precisely as,

We have introduced the higher order finite difference alternating direction implicit scheme for two-dimension- al wave equation. In the next section we will derive another higher order scheme based on Padé approximations.

Therefore,

Similarly the discretization in the t, y and z directions are

substituting (27)-(29) into (1) we have

The scheme presented in (30) is a 4^{th}-order accurate both in time and space for the 2-dimensional acoustic wave equation based on Padé approximation.

A compact finite difference scheme comprises of adjacent point stencils of which differences are taken at the middle node, therefore typically 3, 9 and 27 nodes are used for compact finite difference descretization in one, two and three dimensions, respectively. Whereas a non-compact stencil comprises of any number of nodes which are near the vicinity of the middle node. Usually, when the order of accuracy of a non compact scheme is greater than two, then the scheme is said to be higher order.

The standard alternating direction implicit scheme which is fourth order both in time and space is given in equa- tion (22)-(23). Recently, Das [

coustic wave equation using Padé approximation, we will discuss this method in detail.

Multiplying both sides of (30) by the operator

we get

which upon simplification leads to

Further simplification give

Neglecting the term

scheme, which is fourth order in time and space,

This scheme needs to be solved in two directions simultaneously, however it is easier to solve it in one direction at a time, therefore we can write it as

The above scheme is based on Pad’e approximation and which can be referred to as CPD-ADI scheme.

Substituting (16) into (36) we obtain

Equation (37) is a tridiagonal system and therefore it can be solved easily by using any tri-diagonal solver.

In the next section we will derive a numerical scheme which is a combination of the compact and non-com- pact operators that we introduced in chapter two.

The higher order compact scheme (35)-(36) derived using Padé approximation is efficient and accurate up to order 4. However, the scheme suffers from a moderate numerical dispersion [

Therefore, for the two-dimensional wave equation we have

Following the same procedure as we did for (32)-(34) we obtain

Therefore, for implementation, this can further be written as an ADI scheme

Equations (40)-(41) is an ADI scheme obtained by mixing compact and non-compact difference operators and therefore it can be referred to as non-compact padé based scheme (NCPD-ADI).

For the implementation, we will apply Equation (40)-(41) for the interior nodes, whereas at the nodes near the horizontal boundary we will use (35)-(36). Therefore, the combination of the two schemes can be called compact Padé based interlinked alternating direction implicit scheme (IPD-ADI).

At this stage, we would like to mention that a typical numerical method may produce spurious solution, a phenomena which is greatly attributed to dispersion. Thus, in the next section we will study this phenomena in detail.

Finite difference schemes are widely used in solving time dependent partial differential equations numerically, however in approximation of derivatives using Taylor’s series, truncation of terms produces inaccuracies in finding solutions for such differential equations. Due to these inaccuracies undesirable ripples are produced in the process of discretization. Eventhough such equations are not frequency dependent, due to truncation of higher order terms frequency dependent numerical errors are produced. This numerical phenomena is called dispersion. Thus, even for non-dispersive partial differential equations their discrete approximations are dispersive. Moreover in descritizing a differential equation, depending on our interest in accuracy, we ignore higher order terms, due to these inaccuracies errors are introduced which results in instability of the finite difference scheme.

Dispersion RelationWe will consider a scalar wave Equation (1) which has solution of the form

Then its discrete form is given by

where

discrete directional wave numbers are given by

and

The dispersion relation is given by

where

In numerical analysis the relation between the numerical angular frequency ω and k, v and numerical parameters τ and h is called numerical dispersion relation.

Dispersion error therefore can be considered as the difference between the numerical and analytical frequencies or the ratio of these frequencies. Consequently, it can be measured as the difference between, or the ratio of numerical and analytical velocities. Since the latter method is the most convenient to work with, we will use this method for our dispersion measurement. Thus, mathematically,

where v is numerical velocity and δ is the normalized phase error.

Therefore

Applying finite difference operators to (43) gives as

For two-dimensional HOC-ADI scheme the phase relation is derived as follows. Substituting (43) into (21) we obtain

Similarly the normalized phase relation for CPD and NCPD schemes in two dimensions respectively are given by

and

where

The plots for these normalized dispersion relations are given in the subsequent sections.

A numerical scheme is said to be stable if the errors produced, due to descritization, at one time-step of the calculation do not propagate as computation continues. For a stable scheme, the errors decay and eventually damp out, whereas for the unstable scheme the errors grow with time.

In this section we will analyse stability for HOC, CPD and NCPD schemes using Von Neumann stability analysis.

Therefore, the CFL condition for the HOC scheme is derived as follows

Equation (43) can be written as

Substituting (56) in to (21) we obtain

We consider

where

and

For stability the root condition for (58) should be satisfied. Therefore,

Similarly for the CPD scheme we find the stability condition as

where

and

And for the NCPD scheme the stability condition is

where

and

Hence from (61), (62) and (65) we can observe that the three schemes (HOC, CPD and NCPD) are conditionally stable. Note that most of the other finite difference schemes for the acoustic wave equation are conditionally stable, see e.g., [

In the previous section, we have discussed the dispersion relation for three finite difference schemes which are HOC, CPD and NCPD schemes. In the next two sections, we will provide a comparison of the dispersion relations of these schemes. We present plots of the dispersion relations as a function of propagation angle θ and number of cells per wavelength (kh) for these schemes. Therefore, we discussed the normalized phase error with respect to different parameters and we approached the comparison of the dispersion errors by plotting the norm of the errors and relative phase errors.

Normalized Phase ErrorWe can observe that the numerical dispersion of the above mentioned methods are dependent on three variables. Therefore the normalized phase errors for these methods are studied based on these three factors.

• CFL numbers (r)

• Propagation angle (θ)

• Number of cells per wavelength (kh)

• Effect of different CFL number

For the analysis given here, we keep the grid size equal in all spatial directions, i.e.,

In

In

1) Effect of different propagation angles (θ)

In this subsection, we present plots of the three schemes to study the effect of propagation angle on numerical

dispersion. We choose

In Figures 5-7, we provide the normalized phase error with different propagation angles of the wave. We plotted the phase error versus kh for different values of θ varying from

2) Effect of different number of cells per wavelength (kh)

To examine the effect of number of cells per wavelength, we plot the dispersion errors versus θ for different values of number of cells per wavelength for the three schemes.

In Figures 8-10 we presented the normalized phase error versus propagation angle θ for different values of kh. It is observed that the phase error increases as the propagation angle increases for all the schemes. Though all

the schemes show lower dispersion error for

In this section we presented plots of the

Presented in

The relative phase error is calculated with the following formula Relative error =

the reference value and

In

We have discussed the relative phase errors for different schemes in the individual plots above. As we have discussed in the earlier section, phase error is dependent of three factors. Therefore, it is better to see the relative phase error for the combined effect of these factors. The plots in Figures 13-15 show two-dimensional relative

phase errors, where kh represents polar radius, θ is the polar angle and the plot shows for four values r.

From Figures 13-15 we clearly observe that smaller r values show lower dispersion. It is also shown that waves whose propagation direction makes an angle of 45 degree with the axis are the most accurate for CPD and NCPD schemes, whereas for HOC scheme waves propagating in the axial direction shows relatively low error

for some values of r. Moreover, as clearly indicated, CPD and NCPD schemes show large low dispersion region than HOC scheme.

In this thesis, three ADI schemes are presented for 2D and 3D acoustic wave equation with constant velocity. We derived HOC-ADI, CPD-ADI and NCPD-ADI schemes. The former scheme is based on Taylor approximation and the latter schemes are based on Padé approximation. We perform stability analysis for the mentioned

schemes. Dispersion relations for these schemes have been also obtained followed by a thorough discussion of these relations by plotting them for different parameter sets with respect to different values of kh and θ. From this analysis, it is observed that an increase in these parameters increases the dispersion error and vice-versa. It is observed that numerical dispersion is dependent of these parameters in general for the mentioned schemes.

Comparison of phase errors for these schemes is also carried out with the given parameters. Phase error comparison of the three schemes showed that NCPD scheme has relatively lower dispersion error than CPD and HOC schemes. However, in this work we could only compare dispersion properties for three higher order numerical schemes, therefore in the future, we will analyse dispersion properties of a large group of schemes so that we will be able to explore higher order accurate low dispersive schemes.

Yahya AliAbdulkadir, (2015) Comparison of Finite Difference Schemes for the Wave Equation Based on Dispersion. Journal of Applied Mathematics and Physics,03,1544-1562. doi: 10.4236/jamp.2015.311179