High Resolution Compact Finite Difference Schemes for Convection Dominated Problems

Abdullah Shah, Saher Akmal Khan, Najib Ullah

Department of Mathematics, COMSATS University Islamabad, Islamabad, Pakistan.

DOI: 10.4236/wjm.2019.912017                

Abstract

In this short article, the upwind and central compact finite difference schemes for spatial discretization of the first-order derivative are analyzed. Comparison of the schemes is provided and the best discretization scheme for convection dominated problems is suggested.

Keywords

Upwind Scheme, Upwind Compact Scheme, Flux-Difference Splitting, Convection Dominated Problems

1. Introduction

With the ever-increasing interest in numerical calculations demanding high accuracy for a wide range of length scales, such as large-eddy simulation and direct numerical simulation of turbulence, high-order numerical methods are desired. Particularly, high-order finite difference, finite volume, and finite element methods have received more attention in handling complex problems. These high-order methods try to achieve high accuracy and avoid spurious oscillations and are usually characterized by their self-adaptive nature. The use of high-order methods is particularly warranted by the need to simulate flows containing discontinuous phenomena, such as fluid interfaces and steep shear layers. The compact high-order finite difference schemes provide an effective way of combining the robustness of finite difference schemes and the accuracy of spectral methods [1] [2] [3]. Generally, the computation of derivatives in compact finite differences is implicit in the sense that the derivative values at a particular node are computed not only from the function values but also from the values of the derivative at the neighboring nodes [4]. Compared to non-compact counterparts of the same order of accuracy, compact schemes utilize a smaller stencil, have smaller truncating errors, and give better resolution especially at higher wave numbers [5] [6]. Compact finite difference schemes can generally be classified into two broad categories: upwind and central. The upwind compact schemes inherently possess the needed dissipation to control the numerical instabilities. Fu and Ma [7] have developed some upwind compact schemes which are successfully implemented by Shah et al. [8] [9] [10] [11] for solving fluid flow problems. As these schemes possess appropriate dissipation to prevent non-physical oscillations, they seem to be suitable for solving the convection dominated problems. N.B. Ali et al. [12] used implicit and explicit third and fifth-order upwind compact schemes for solving the level set equation. De V. E. and Eswaran, V. [13] have studied some optimized upwind and upwind compact schemes for the solution of acoustic wave problem. Central compact schemes have the advantage of achieving high-order accuracy with fewer grid points in the stencil, but they are non-dissipative, and using central compact schemes on non-staggered meshes for convection terms might cause numerical oscillations even for flows without discontinuities. Reducing or removing such oscillations requires the introduction of dissipation terms or the use of filtering approach [14]. Resolution characteristics imply how compact finite difference approximation represents the exact result over the full range of length scales that can be realized for a given mesh [15]. This work aims to study different compact schemes to find the scheme more suitable for solving convection dominated problems.

2. Model Problem

In order to examine approximating behaviors of various numerical schemes, the following linear convection equation (also known as one-way wave equation) is considered.

$\frac{\partial u}{\partial t}+c\frac{\partial u}{\partial x}=0,\mathrm{<mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>\; <mi>\; </mi>}c>0.$ (1)

The semi discrete form Equation (1) is

$\frac{\partial u_j}{\partial t}+c\frac{\partial u_j}{\partial x_j}=0.$ (2)

The solution of Equation (1) represented by $u\left(x\mathrm{,}t\right)$ by a typical Fourier mode is given by:

$u\left(x,t\right)={\hat{u}}_k\left(t\right){\text{e}}^{ikx},$ (3)

${\hat{u}}_k$ is the Fourier mode of the wave number k and $i=\sqrt{-1}$, the exact spatial differentiation of Equation (3) is represented by;

$u^{\prime }=i\left(kh\right)\frac{{\hat{u}}_k}{h}{\text{e}}^{ikx},$ (4)

where the wave number is scaled by the grid size $h=\frac{l}{n}$, where l is the length of

domain and n is the number of grids. By analogy the numerical approximation of the derivative is written as [13]

$u^{\prime }=k_{eq}\left(kh\right)\frac{{\hat{u}}_k}{h}{\text{e}}^{ikx}=\left(k_r\left(kh\right)+i\left(k_i\right)kh\right)\frac{{\hat{u}}_k}{h}{\text{e}}^{ikx}.$ (5)

The exact solution of Equation (1) is $u\left(x\mathrm{,}t\right)={\text{e}}^{ik\left(x-ct\right)}$, and the exact solution

of Equation (2) can be written as $u\left(x_j,t\right)={\text{e}}^{-k_r\frac{ct}{\Delta x}}{\text{e}}^{ik\left(x_j-\frac{k_i}{k\Delta x}ct\right)}$, where the modified

wave number $k_{eq}=k_r+i{k}_i$. $k_i$ is related to the phase speed in the numerical solution, and $k_r$ is related to the numerical damping of a difference scheme. Fourier analysis of different discretization schemes allows us to choose the best scheme.

2.1. Upwind Compact Scheme

In this subsection, third and fifth-order upwind compact and upwind explicit schemes are analyzed. For the third-order upwind compact scheme [16], we have

$\frac{2}{3}{u}_j+\frac{1}{3}{u}_{j-1}=\frac{u_{j+1}+4u_j-5u_{j-1}}{6},$ (6)

that satisfy the relation

$\begin{array}{c}{u}_j=\frac{\frac{1}{6}{\text{e}}^{i\alpha }+\frac{2}{3}-\frac{5}{6}{\text{e}}^{-i\alpha }}{\frac{2}{3}+\frac{1}{3}{\text{e}}^{-i\alpha }}\hat{u}\left(t\right){\text{e}}^{ik{x}_j}\\ =\frac{\cos \alpha +i\sin \alpha +4-5\left(\cos \alpha -i\sin \alpha \right)}{4+2\left(\cos \alpha -i\sin \alpha \right)}\hat{u}\left(t\right){\text{e}}^{ik{x}_j}\\ =\frac{{\left(1-\cos \alpha \right)}^2+i\sin \alpha \left(8+\cos \alpha \right)}{5+4\cos \alpha }\hat{u}\left(t\right){\text{e}}^{ik{x}_j}\end{array}$ (7)

where

$k_r=\frac{{\left(1-\cos \alpha \right)}^2}{5+4\cos \alpha },k_i=\frac{\sin \alpha \left(8+\cos \alpha \right)}{5+4\cos \alpha }.$ (8)

Similarly, for the fifth-order upwind compact scheme [7], we have

$\frac{3}{5}{u}_j+\frac{2}{5}{u}_{j-1}=\frac{-u_{j+2}+12u_{j+1}+36u_j-44u_{j-1}-3u_{j-2}}{60},$ (9)

with

$\begin{array}{c}{u}_j=\frac{-\frac{1}{60}{\text{e}}^{2i\alpha }+\frac{12}{60}{\text{e}}^{i\alpha }+\frac{36}{60}-\frac{44}{60}{\text{e}}^{-i\alpha }-\frac{3}{60}{\text{e}}^{-2i\alpha }}{\frac{36}{60}+\frac{24}{60}{\text{e}}^{-i\alpha }}\hat{u}\left(t\right){\text{e}}^{ik{x}_j}\\ =\frac{1}{117+108\cos \alpha }(24-6{\cos }^3\alpha -18{\sin }^2\alpha -18\sin \alpha \\ +i\left(6{\sin }^3\alpha +45\sin \alpha \cos \alpha +180\sin \alpha \right))\hat{u}\left(t\right){\text{e}}^{ik{x}_j}\end{array}$ (10)

where

$k_r=\frac{2\left(1-3\cos \alpha +3{\cos }^2\alpha -{{\displaystyle \cos }}^3\alpha \right)}{3\left(13+12\cos \alpha \right)},k_i=\frac{\sin \alpha \left(2{\sin }^2\alpha +15\cos \alpha +60\right)}{3\left(13+12\cos \alpha \right)}.$ (11)

For the explicit third-order upwind scheme [7],

$k_r=\frac{1}{6}\left(3-4\cos \alpha +\cos 2\alpha \right),k_i=\frac{1}{6}\left(8\sin \alpha -\sin 2\alpha \right),$ (12)

and for the explicit fifth-order upwind scheme, we have,

$k_r=\frac{10-15\cos \alpha +6\cos 2\alpha -\cos 3\alpha }{30},k_i=\frac{45\sin \alpha -9\sin 2\alpha +\sin 3\alpha }{30}.$ (13)

Figure 1 shows variations of $k_r$ and $k_i$ with the reduced wave number $\alpha$ for the above four schemes. We can see the fifth-order schemes can approximate the exact damping ( $k_r^E=0$ ) to higher waver numbers than the third-order schemes, and the compact schemes can approximate the exact dispersion relation ( $k_i^E=\alpha$ ) better than the non-compact schemes.

Table 1 gives the upper limit of the reduced wave number, which corresponds to a point in Figure 1 where $k_r$ or $k_i$ begins to reach 2% errors relative to their exact solutions respectively. Larger upper limit implies fewer grid points are needed to resolve a given physical structure. For example, to approximate the exact wave speed within 2% error, the ratio of grid points needed by the 5th-order upwind compact scheme to those needed by the 5th-order upwind biased scheme is $1.25/1.71=0.73$ in one dimensional case. In three-dimensional case, this ratio becomes ${\left(1.25/1.71\right)}^3=0.39$, resulting in significant saving in computer resources.

2.2. Central Compact Schemes

In this section, various compact finite difference schemes are studied. The family of cell centered central compact schemes given by Lele et al. [3] is given by:

$\nu {u^{\prime }}_{i-2}+\mu {u^{\prime }}_{i-1}+{u^{\prime }}_i+\mu {u^{\prime }}_{i+1}+\nu {u^{\prime }}_{i+2}=c\frac{u_{i+3}-u_{i-3}}{6h}+c\frac{u_{i+2}-u_{i-2}}{4h}+c\frac{u_{i+1}-u_{i-1}}{2h}$ (14)

The order of these schemes can be based parameters values as shown in Table 2.

Taking Fourier transform of Equation (14), we have

Figure 1. Variations of $k_r$ and $k_i$ vs. $\alpha$ for the compact and non-compact schemes.

Table 1. Upper limits of the reduced wave number when $k_r$ and $k_i$ of the difference schemes first exceed 2% errors relative to exact solutions.

Table 2. Values of parameters involved for the central compact scheme.

$k_i=\frac{a\sin \alpha +\frac{b}{2}\sin 2\alpha +\frac{c}{3}\sin 3\alpha }{1+2\mu \cos \alpha +2\nu \cos 2\alpha };k_r=0.$ (15)

The different values of $k_i$ are given in Table 3.

The difference between modified wave number and exact wave number is very small, therefore these schemes have spectral like resolution. The comparison of various central compact schemes is presented in Figure 2. The eighth-order central compact scheme seems to follow the exact wave number more closely than all other central compact schemes, though it has a broader stencil width.

2.3. Comparison of Upwind and Central Compact Scheme

In this subsection, the upwind and central compact schemes are compared based upon the resolution characteristics $k_i$ vs $\alpha$. For this purpose, two upwind compact schemes and two central compact schemes are selected from the previous sections.

The comparison plot for $k_i$ vs $\alpha$ is shown in Figure 3.

The comparison of the scheme enables us to find the scheme best suitable from the chosen schemes. Figure 3 shows that the upwind compact schemes give the better resolution amongst all the schemes while central compact schemes have poor resolution. So in order to improve the resolution of central schemes, filtering is required.

3. Conclusion

We have analyzed upwind, upwind compact and central compact schemes of

Table 3. Values of parameters involved for the central compact scheme.

Figure 2. Comparison of various central compact schemes for $k_i$ vs $\alpha$.

Figure 3. Comparison of upwind and central compact schemes $k_i$ versus $\alpha$.

different order accuracy for numerical investigation of convection equation. It is observed that the use of the upwind compact scheme makes the numerical solution more stable as compared with the central scheme and can be used for convection dominated problems. A comparison is also given with non-compact schemes of the same order of accuracy with almost the same computational cost.

Acknowledgments

The work of Abdullah Shah is supported by HEC under NRPU No. 7781 and PSF No. 5651.

Conflicts of Interest

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

References

[1]	Hirsh, R.S. (1975) High Order Accurate Difference Solutions of Fluid Mechanics Problems by a Compact Differencing Technique. Journal of Computational Physics, 19, 90-109. https://doi.org/10.1016/0021-9991(75)90118-7
[2]	Rubin, S.G. and Khosla, P.K. (1977) Polynomial Interpolation Methods for Viscous Flow Calculations. Journal of Computational Physics, 2, 217-244. https://doi.org/10.1016/0021-9991(77)90036-5
[3]	Lele, S.K. (1992) Compact Finite Difference Schemes with Spectral-Like Resolution. Journal of Computational Physics, 103, 16-42. https://doi.org/10.1016/0021-9991(92)90324-R
[4]	Sheng, T.Y. (1991) Runge-Kutta Methods Combined with Compact Difference Schemes for the Unsteady Euler Equations. Center for Modeling of Turbulence and Transition-Research Briefs, 93, 15802.
[5]	Zhou, Q., Yao, Z.F. and Shen, M.Y. (2007) A New Family of High-Order Compact Upwind Difference Schemes with Good Spectral Resolution. Journal of Computational Physics, 227, 1306-13391. https://doi.org/10.1016/j.jcp.2007.09.008
[6]	Shah, A., Yuan, L. and Islam, S. (2012) Numerical Solution of Unsteady Navier-Stokes Equations on Curvilinear Meshes. Journal of Computer and Mathematics with applications, 63, 1548-1556. https://doi.org/10.1016/j.camwa.2012.03.047
[7]	Fu, D.X. and Ma, Y.W. (1997) A High-Order Accurate Difference Scheme for Complex Flowfields. Journal of Computational Physics, 134, 1-5. https://doi.org/10.1006/jcph.1996.5492
[8]	Shah, A. and Yuan, L. (2010) Upwind Compact Finite Difference Scheme for Time-Dependent Incompressible Navier-Stokes Equations. Applied Mathematics and Computation, 215, 3201-3213. https://doi.org/10.1016/j.amc.2009.10.001
[9]	Rizwan, M., Shah, A. and Yuan, L. (2016) A Central Compact Scheme for Numerical Solution of Two-Phase Incompressible Flow Using Allen-Cahn Phase Field Model. Journal of the Brazilian Society of Mechanical Sciences and Engineering, 38, 433-441. https://doi.org/10.1007/s40430-015-0342-4
[10]	Khan, S. and Shah, A. (2019) Simulation of the Two-Dimensional Rayleigh-Taylor Instability Problem by Using Diffuse-Interface Model. AIP Advances, 9, Article ID: 085312. https://doi.org/10.1063/1.5100791
[11]	Saeed, S., Shah, A. and Khan, S. (2018) Numerical Investigation of Bubbles Coalescence in a Shear Flow with Diffuse-Interface Model. Heliyon, 4, e01024. https://doi.org/10.1016/j.heliyon.2018.e01024
[12]	Borujerdi, A.N. and Kebriaee, A. (2003) Upwind Compact Implicit and Explicit High Order Finite Difference Scheme for Level Set Techniques. International Journal for Computational Methods in Engineering Science and Mechanics, 3, 308-318.
[13]	De, A.K. and Eswaran, V. (2006) Analysis of a New High Resolution Upwind Compact Scheme. Journal of Computational Physics, 218, 398-416. https://doi.org/10.1016/j.jcp.2006.02.020
[14]	Gaitonde, D.V. and Visbal, M.R. (2000) Pade-Type Higher-Order Boundary Filters for the Navier Stokes Equations. AIAA Journal, 38, 2103-2112.
[15]	Liang, X., Zhang, S., Zhang, H. and Shu, C. (2013) A New Class of Central Compact Scheme with Spectral Like Resolution I. Linear Schemes. Journal of Computational Physics, 248, 235-256. https://doi.org/10.1016/j.jcp.2013.04.014
[16]	Shah, A. and Yuan, L. (2009) Flux-Difference Splitting Based Upwind Compact Scheme for the Incompressible Navier-Stokes Equations. International Journal for Computational Methods in Engineering Science and Mechanics, 61, 552-568. https://doi.org/10.1002/fld.1965