An adaptive grid method for a singularly perturbed convection-diffusion equation with a discontinuous convection coefficient

1.   Introduction

Singularly perturbed problems have been widely used to describe various models of physics and engineering ^[1,2]. Typical examples include the Navier-Stokes equation with large Reynolds number in fluid dynamics, the convective heat transport with a large Péclet number, etc. A notable feature of these problems is that the high-order derivative term is multiplied by a small positive parameter $\varepsilon$. In general, when $\varepsilon\rightarrow 0$, the solutions of these problems exhibit boundary layers or inner layers, which are basically thin regions in the neighbourhood of the boundary or interior of the domain. Thus, it is a challenge to obtain a reliable numerical solution due to the existence of boundary layer and/or inner layers in the continuous solution. To overcome this difficulty, one effective method is to use some nonuniform meshes that are fine where layers appear in the solution. To the best of our knowledge, such meshes can be divided into two classes: layer-adapted meshes (Shishkin mesh and Bakhvalov mesh) and adaptive grids, see the monographs ^[2] and references therein for dedicated discussions on such meshes and relevant numerical methods. In the past two decades, layer-adapted mesh approaches have attracted considerable attention in the numerical methods of singularly perturbed problems community, especially for the singularly perturbed ordinary differential equations with a discontinuous convection coefficient. The authors in ^[3] considered a singularly perturbed convection-diffusion problem with a discontinuous convection coefficient and presented an almost first-order uniformly convergent numerical scheme on a Shishkin-type mesh. Cen ^[4] considered a singularly perturbed convection-diffusion equations with discontinuous convection coefficient and constructed a second-order hybrid difference scheme on a Shishkin mesh. Shanthi et al. ^[5] developed a classical upwind finite difference scheme on a Shishkin mesh for a singularly perturbed second-order ordinary differential equation with two parameters and a discontinuous source term and proved that their proposed method was almost first-order uniformly convergent. In ^[6,7,8], the authors constructed some parameter-uniform hybrid finite difference schemes on a Shishkin-type mesh for singularly perturbed convection-diffusion problems with a discontinuous source term and a discontinuous convection coefficient. Tamilselvan and Ramanujam ^[9] constructed a finite difference scheme on a Shishkin mesh to solve a weakly coupled system of two singularly perturbed convection-diffusion equations with discontinuous convection coefficients and gave the rigorous proof of the parameter uniform convergence in the global maximum norm. Pathan and Vembu ^[10] proposed a parameter-uniform hybrid numerical method on a Shishkin mesh to solve a weakly coupled system of two singularly perturbed convection-diffusion equations with discontinuous convection coefficients and source terms. Aarthika et al. ^[11] considered a two-dimensional singularly perturbed reaction-diffusion equation with a discontinuous source term and constructed a hybrid finite difference method on a piecewise uniform Shishkin mesh. Further, they proved that the hybrid finite difference method was almost second-order uniformly convergent with respect to the perturbation parameter.

While the $\varepsilon$-uniformly convergent layer-adapted grid methods have been applied successfully to singularly perturbed ordinary differential equations with a discontinuous coefficient, a lot of researchers pay attention to discuss the layer-adapted mesh approaches for singularly perturbed time-dependent/independent problems with a discontinuous coefficient (or discontinuous source term), see, e.g., ^{[12,13,14,15,16]}. On this basis, the layer-adapted grid method for a system of singularly perturbed parabolic problems with a discontinuous coefficient (or discontinuous source term) has been discussed in the literature. For example, Rao and Chawla ^[17] considered a parameter-uniform numerical method for a parabolic system with an arbitrary number of linear singularly perturbed equations of reaction-diffusion type coupled in the reaction terms with a discontinuous source term and proved that their numerical method was uniformly convergent of first order in time and almost second order in the spatial variable. The authors in ^[18] proposed a finite difference scheme on a Shishkin mesh to the solution of a two-parameter singularly perturbed convection-diffusion-reaction system of partial differential equations with discontinuous coefficients. Rao and Chaturvedi ^[19] analyzed a numerical method for a coupled system of two singularly perturbed parabolic semilinear reaction-diffusion equations having discontinuous source terms and proved that the proposed method was parameters-uniformly convergent of first-order in time and almost second-order in space.

It should be pointed out that there have been extensive studies on layer-adapted grid approach of singularly perturbed problems with a discontinuous coefficient (or discontinuous source term). However, to the best of our knowledge, limited work has been done in the adaptive grid algorithm based on the a posteriori error estimation for these problems. Hence, the main body of this text is to develop an adaptive grid method for the following singularly perturbed convection-diffusion equation with a discontinuous convection coefficient

where $0 < \varepsilon\ll1$ is a small positive parameter, $g_l$, $g_r$ are two given constants, $d$ is a jump point in any function with $[\omega(d)] = \omega(d^+)-\omega(d^-)$ and $a, f\in C^2(\Omega^d)$. Furthermore, there exist three positive constants $\alpha_i(i = 1, 2)$ and $C$, such that

These hypotheses guarantee that this problem (1.1) has a solution $u\in C^1(\Omega)\bigcap C^2(\Omega^d)$ (see [3, Theorem 1]. Moreover, there exist an interior layer in the vicinity of the point of discontinuity $x = d$.

The constructive organization of this article is as follows: some facts about the exact solution $u$ and the corresponding discretization scheme of problem (1.1) are listed in Section 2. Moreover, the stability bound for the calculated solution $u_i^N$ (on an arbitrary grid) is given. Then, in Section 3, a maximum norm a posteriori error estimate is derived, see Theorem 3.1, which is the most fundamental result of our paper. In Section 4, we construct an adaptive grid generation algorithm by using our presented a posteriori estimate to monitor-function equidistribution. Numerical results are presented in Section 5 that sustain our theoretical estimate. Finally, Section 6 is a summary of our conclusions.

Notation 1.1. Throughout the paper, $C$ will denote a generic positive constant that is independent of $\varepsilon$ and of the mesh parameter $N$. It may take different values in different places. For a given continuous function $v(x)$ on $\bar{\Omega} = [0, 1]$, the $L_{\infty}$ norm is defined by $\left\|v(\cdot)\right\|_{\infty, \bar{\Omega}} = \sup \left|v(x)\right|$. Meanwhile, for a given mesh function $v^N = \left(v_0^N, v_1^N, \cdots, v_N^N\right)$, we define the discretization maximum norm $\left\|v^N\right\|_{\infty} = \max\limits_{0\leq i\leq N}\left|v_i^N\right|$.

2.   Preliminary results

In view of ^[3], the differential operator $\mathcal{L}$ fulfillment maximum principle. Thus, Eq (1.1) has a unique solution $u$, which has the following bound.

Lemma 2.1. ^[3] Let $u\in C(\bar{\Omega})\cap C^2(\Omega^d)$ be the exact solution of problem (1.1). Then

where $\gamma = \min\left\{\alpha_{1}/d, \alpha_{2}/(1-d)\right\}$.

Corollary 2.1. For any two arbitrary functions $v(x)$ and $w(x)$, satisfy $v(0) = w(0), v(1) = w(1)$ and

where $F(x)$ is a piecewise continuous function, then

Proof. For any two arbitrary functions $v(x)$ and $w(x)$, let $\mu = v-w$. Obviously $\mu(0) = \mu(1) = 0$. Then the desire result can be followed from Lemma2.1.

□

For our numerical method, we construct an arbitrary mesh $\bar{\Omega}^N_d = \left\{x_i\right\}_{i = 0}^N$ with $x_s = d$, $x_0 = 0$ and $x_N = 1$, where $1 < s < N$ is a index. Let $h_i = x_i-x_{i-1}, \; i = 1, \cdots, N$ be the local mesh step size. Then for a given grid function $\left\{z_i\right\}_{i = 0}^N$, we define some difference operators as follows:

where $\hbar_i = \frac{h_{i+1}+h_i}{2}, \; i = 1, \cdots, N-1.$ Furthermore, we discretize the problem (1.1) on the above mesh $\bar{\Omega}^N_d$ using the following finite difference scheme:

where $u_i^N$ is the approximation solution of $u(x)$ at point $x = x_i$, $a_{i-\frac{1}{2}} = a\left((x_{i-1}+x_i)/2\right), \; \; a_{i+\frac{1}{2}} = a\left((x_{i}+x_{i+1})/2\right)$, $f_{i-\frac{1}{2}}, \; f_{i+\frac{1}{2}}$ are similar to $a_{i-\frac{1}{2}}, \; a_{i+\frac{1}{2}}$.

Lemma 2.2. ^[3] Let $u^N$ be the solution of the discrete scheme (2.2). Then

where $\gamma = \min\left\{\alpha_{1}/d, \alpha_{2}/(1-d)\right\}$.

3.   A posteriori error estimation

In order to obtain a posteriori error estimation of the solution $u^N$ of Eq (2.2), we first define a piecewise quadratic function $\tilde{u}^N(x)$ on the interval $J_i = [x_{i-1}, x_{i}], \; i = 1, 2, \cdots, N$ as follows:

Then

Theorem 3.1. Let $u(x)$ be exact solution of the problem (1.1), $u^N_i$ be the discrete solution of the discrete scheme (2.2) and $\tilde{u}^N(x)$ be the piecewise quadratic function defined by (3.1). Then we have

where

Proof.

First, for $\forall x\in(x_{i-1}, x_i), \; 1\leq i\leq s$, it follows from Eqs (1.1), (2.2) and (3.3) that

Then applying Taylor formula to $a(x)$ and $f(x)$ at $x = x_{i-\frac{1}{2}}$, we obtain

where $\xi_1, \xi_2\in(x_{i-\frac{1}{2}}, x)$. For $x\in(x_{i-1}, x_i), \; s\leq i < N$, similar to Eqs (3.4) and (3.5), it is easy to get

Furthermore, for $x\in(x_{N-1}, x_N)$, one has

Thus, it follows from (3.4)-(3.7) that

Finally, using Corollary 2.1 and Eq (3.8), yields,

This completes the proof.

□

4.   An adaptive grid algorithm

To this day, there are many researchers have studied adaptive grid methods for singularly perturbed problems with a continuous convection coefficient and a source term, and have made many significant research results (see ^{[20,21,22,23,24,25,26]} for example). Here, the main contribution of this text is to design an adaptive grid method to solve the singularly perturbed convection-diffusion equation (1.1) for the first time.

As far as we known, for a given positive monitor function $M_i\left(\cdot\right)$, the hinge technique of adaptive grid method is to discover a grid $\left\{x_i\right\}_{i = 0}^N$ such that

Here, Eq (4.1) is called the discrete mesh equidistribution principle. Here, according to the output of Theorem 3.1, we select the monitor function as below:

In an effort to the equidistributed mesh $\bar{\Omega}^N$ and the corresponding numerical solution $u_i^N$, we design the grid generation algorithm as follows:

Step 1. For a given positive integer $N$, choose an initial uniform mesh $\bar{\Omega}^{N, (0)}_d = \left\{0 = x_{0}^{(0)} < x_1^{(0)} < \cdots < x_N^{(0)} = 1\right\}$.

Step 2. For $k = 0, 1, \cdots$ and the grid $\bar{\Omega}^{N, (k)}_d$, there exists a index $J$ such that the grid point $x_J^{(k)}$ satisfying

Then let $x_J^{(k)} = d$.

Step 3. Let $\left\{u_i^{N, (k)}\right\}_{i = 0}^N$ be the solution of discretization scheme (2.2) on $\bar{\Omega}^{N, (k)}_d$. Set

where $M^{(k)}_i = 1+\left|D^-u^{N, (k)}_i\right|+\sqrt{\left|D^2u^{N, (k)}_i\right|}$ and $M^{(k)}_0 = M^{(k)}_1, \; \; M^{(k)}_N = M^{(k)}_{N-1}$.

Step 4. Let $h_{i}^{(k)} = x_{i}^{(k)}-x_{i-1}^{(k)}, \; i = 1, \cdots, N$ and set $\Phi_{0}^{(k)} = 0$ and

Then let $\phi^{(k)}(s)$ be the piecewise linear interpolant function through $\left(\Phi_{i}^{(k)}, x_{i}^{(k)}\right)$ and generate the new grid $\bar{\Omega}^{N, (k+1)}$ by using $x_{i}^{(k+1)} = \phi^{(k)}(Y_{i}^{(k)})$.

Step 5. Choose a positive constant $\epsilon$, if the stopping criterion

holds true, go to Step 6, otherwise return to Step 2.

Step 6. Set $\bar{\Omega}^{*}_d = \bar{\Omega}^{N, (k+1)}_d$ and $\left\{u_i^{*}\right\}_{i = 0}^N = \left\{u_i^{N, (k+1)}\right\}_{i = 0}^N$ then stop.

5.   Numerical results and discussion

In order to verify our theoretical result, we took into account the following test question

where

Since the exact solution of this problem (5.1) is not known, we use the following double mesh principle to calculate the maximum point-errors and the corresponding convergence rates: Let $u^N_i$ and $u^{2N}_i$ be the numerical solutions of the discrete scheme (2.2) on mesh $\bar{\Omega}^N$ and $\bar{\Omega}^{2N}$, respectively, where the grid $\bar{\Omega}^{2N}$ is obtained by bisecting the original mesh $\bar{\Omega}^N$. Then the errors and rates of convergence are computed in the usual way:

Next, for $\varepsilon = 2^{-2j}, j = 2, \cdots, 8$ and $N = 2^k, \; k = 7, \cdots, 11$, we use our presented adaptive grid method to solve this test problem. The error and rates of convergence for the numerical solutions are displayed in Table 1. It is shown from Table 1 that the accuracy of our adaptive grid method is first-order, which is confirmed our theoretical result given in Theorem 3.1. Moreover, to illustrate the advantageous of our adaptive grid method, for $\varepsilon = 2^{-3}, 2^{-15}$ and the same $N$, Table 2 gives the numerical results obtained by using our adaptive grid and the Shishkin mesh, respectively, which are evidences that the maximum point-wise errors of our presented adaptive grid method are much better than that obtained by the Shishkin mesh.

To verify our adaptive grid generation algorithm given in Section 4, Figure 1 represents the grid iteration process for $N = 64$ and $\varepsilon = 2^{-8}$. It is shown that the solution of this test problem has a interior layer at $x = 0.5$.

6.   Conclusions

This paper mainly discussed an a posteriori error estimation in maximum norm for a finite difference scheme to the singularly perturbed convection-diffusion equation with a discontinuous convection coefficient. To deal with the jump point $x = d$, we designed an adaptive grid generation algorithm based on the presented a posteriori error estimation and the mesh equidistribution principle. It should be pointed out that the proposed adaptive grid algorithm used in this work can be extended to structure an adaptive grid approach that applies to the other singularly perturbed problems with a discontinuous coefficient and a source term.

Use of AI tools declaration

The authors declare they have not used Artificial Intelligence (AI) tools in the creation of this article.

Acknowledgments

This work is supported by the Natural Science Foundation of Guangxi Province (2020GXNSFAA159010) and the Guangxi University Young and Middle-aged Teachers' Basic Research Capability Enhancement Project (2023KY0394), and Innovation Project of Guangxi Graduate Education (YCSW2023438).

Conflict of interest

The authors declare that there is no conflict of interest.