Numerical Solution of Quasilinear Singularly Perturbed Problems by the Principle of Equidistribution

Quan Zheng, Fulin Ye

North China University of Technology, Beijing, China.

DOI: 10.4236/jamp.2020.810163                

Abstract

In this paper, the numerical solution and its error analysis of quasilinear singular perturbation two-point boundary value problems based on the principle of equidistribution are given. On the non-uniform grid of the uniformly distributed arc-length monitor function, the solution of the simple upwind scheme is obtained. It is proved that the adaptive simple upwind scheme based on the principle of equidistribution has uniform convergence for small perturbation parameters. Numerical experiments are carried out and the error analysis are confirmed.

Keywords

Quasilinear Singularly Perturbed BVP, Equidistribution, Adaptive Mesh, Uniform Convergence

1. Introduction

We consider a quasilinear singularly perturbed two-point boundary value problem:

$\{\begin{array}{l}Tu\left(x\right):=-\epsilon u^{″}\left(x\right)-b{\left(x,u\left(x\right)\right)}^{\prime }+c\left(x,u\left(x\right)\right)=0,x\in \left(0,1\right),\\ u\left(0\right)=u\left(1\right)=0,\end{array}$ (1)

where, $0<\epsilon \ll 1$ is a sufficiently small positive perturbation coefficient, b and c are sufficiently smooth functions, and for any $x\in [0,1],u\in R$ satisfies $0<\beta \le b_u(x,u)\le {\beta }^{*}$, $0\le b_x(x,u)\le C$ and $0\le c_u(x,u)\le {\gamma }^{*}$, $0\le c_x(x,u)\le C$ respectively where $\beta$, ${\beta }^{*}$, ${\gamma }^{*}$ and C is a generic positive constant, independent of the perturbation parameter $\epsilon$. In this case, the above problem contains a unique solution u, and u has an exponential boundary layer at $x=0$ (see [1] [2] [3]).

Linß considered the problem (1) with $c\left(x,u\left(x\right)\right)=-f\left(x\right)$ in [4]. The problem was discretized using a simple upwind finite difference scheme on adapted meshes using grid equidistribution of monitor functions. Kopteva and Stynes considered the same case in [5]. In [6], Jugal and Srinivasan considered a linear problem and proved the uniform first-order convergence of the numerical solution on adapted meshes by equidistributing the arc-length monitor function. In [7], Linß, Roos and Vulanović considered the problem (1), which was discretized using a nonstandard upwinded first-order difference scheme on generalized Shishkin-type mesh. Linß considered the problem (1) on Bakhvalov-type mesh in [8], on BS-type mesh in [9]. Zheng, Li, and Gao proposed a class of hybrid difference schemes with variable weights on Bakhvalov?Shishkin mesh in [10] and the parameter-uniform second-order convergence of approximating to the solution and the derivative on BS-type mesh and that of nearly second-order on Shishkin mesh were proved.

In this paper, we consider the problem (1) in the general case. The problem is discretized using a simple upwind finite difference scheme on adapted meshes using grid equidistribution of the arc-length monitor function. We will prove the uniform first-order convergence and the numerical experiments verify our analysis.

2. Simple Upwind Scheme and a Posteriori Error Estimate

We set:

$Au(x)=\epsilon u^{\prime }(x)+b(x,u(x))-{\displaystyle {\int }_0^{x}c(x,u(x))dx}$ (2)

So, the problem (1) is equivalent to:

$Tu(x)=-[Au(x){]}^{\prime }$.

For any nodes on interval [0, 1], $0=x_0<x_1<\cdots <x_n=1$, we set:

$h_i=x_i-x_{i-1}$,

$\begin{array}{l}{D}^{-}{u}_i=(u_i-u_{i-1})/h_i,i=1,\dots ,N,\\ D^{+}{u}_i=(u_{i+1}-u_i)/h_{i+1},i=1,\dots ,N.\end{array}$

The simple upwind scheme is:

$T^N{u}_i^N:=-D^{+}\left(A^N{u}_i^N\right)=0,i=1,2,\cdots ,N-1,$ (3)

where

$A^N{u}_i^N:=\epsilon D^{-}{u}_i+b\left(x_i,u_i\right)-{\displaystyle {\sum }_{j=1}^i{h}_{j}c}\left(x_{j-1},u_{j-1}\right).$ (4)

$\left\{u_i^N\right\}$ is the solution computed on the mesh $\left\{x_i\right\}$.

For any grid function $v\left(x\right)$, we define the norm:

${\Vert v\left(x\right)\Vert }_{\infty }={\max }_{x\in \left[0,1\right]}\left|v\left(x\right)\right|$ and ${\Vert v\left(x\right)\Vert }_{*}={\min }_{V:V^{\prime }=v}{\Vert V\left(x\right)\Vert }_{\infty }$.

Easily get: ${\Vert v\Vert }_{*}\le {\Vert v\Vert }_{\infty }$, and ${\Vert D^{-}v\Vert }_{*}\le 2{\Vert v\Vert }_{\infty }$.

Define the linear case of the problem (1):

$Lu\left(x\right):=-{\left[Mu\left(x\right)\right]}^{\prime }=0,x\in \left(0,1\right),u\left(0\right)=u\left(1\right)=0$. (5)

where

$Mu\left(x\right)=\epsilon u^{\prime }\left(x\right)+\bar{b}\left(x\right)u\left(x\right)-{\displaystyle {\int }_0^x\bar{c}\left(x\right)u\left(x\right)dx},$ (6)

$\bar{b}\left(x\right)={\displaystyle {\int }_0^1{b}_u\left(x,w+s\left(v-w\right)\right)ds}$ and $\bar{c}\left(x\right)={\displaystyle {\int }_0^1{c}_u\left(x,w+s\left(v-w\right)\right)ds}$.

where v and w be two arbitrary functions.

To simplify some arguments, we assume that:

$\bar{c}\ge 0,\bar{c}-{\bar{b}}^{\prime }\ge 0$. (7)

The operators L and T are related as follow:

$L\left(v-w\right)=Tv-Tw$. (8)

Lemma 1: Suppose that operator L satisfies inequality (7), then equation (5) has a unique solution $u\in [0,1]$ and ${\Vert u\Vert }_{\infty }\le C{\Vert Lu\Vert }_{*}$.

Lemma 2: Let $v,w\in H^1(0,1)$, and $v(0)=w(0)$, $v(1)=w(1)$, $L(v-w)=0$, then ${\Vert v-w\Vert }_{\infty }\le C{\Vert Tv-Tw\Vert }_{*}$.

Lemma 3: There exists a constant C such that

${\Vert u^N-u\Vert }_{\infty }\le C\underset{i}{\max }\left[{\beta }^{*}\left|u_i^N-u_{i-1}^N\right|+C{h}_i\right]$.

From Lemma 3, we get:

${\Vert u^N-u\Vert }_{\infty }\le C{\max }_i{h}_i\sqrt{1+{\left(D^{-}{u}_i^N\right)}^2}$. (9)

The approach to prove the above lemmas is similar to that in [2], where ${\Vert u^N-u\Vert }_{\infty }\le \underset{i}{\left(2/\beta \right)}\left[{\beta }^{*}{\max }_i\left|u_i^N-u_{i-1}^N\right|+CH\right]$ with $c\left(x,u\left(x\right)\right)=-f\left(x\right)$ was derived.

3. First-Order Uniform Convergence

The adaptive mesh method based on the principle of equal arc-length distribution is applied here.

Define the discrete operator $L^N$ of operator L by

$L^N{u}_i^N:=-D^{+}\left(M^N{u}_i^N\right)=0,i=1,2,\cdots ,N-1,$

where

$M^N{u}_i^N:=\epsilon D^{-}{u}_i+{\bar{b}}_i{u}_i-{\displaystyle \underset{j=1}{\overset{i}{\sum }}{h}_j{\bar{c}}_{i-1}{u}_{i-1}}.$

Let $L^{*,N}$ denote the operator adjoint to $L^N$ with respect to $\left(\cdot ,\cdot \right)$, i.e.,

$\left(L^N{v}_i^N,w_i^N\right)=\left(v_i^N,L^{*,N}{w}_i^N\right)\forall v_i^N,w_i^N\in R^{N-1}.$

Then

$L^{*,N}{u}_i^N=-\epsilon D^{+}{D}^{-}{u}_i+{\bar{b}}_i{D}^{*,-}{u}_i+{\bar{c}}_i{u}_i=0,i=0,1,\cdots ,N,$

where $D^{*,-}{u}_i=\left(u_i-u_{i-1}\right)/h_{i+1}$.

Let $G_{i,j}$ denote the discrete Green’s function associated with the mesh node $x_i$, i.e.,

$L^{*,N}{G}_{i,j}={\delta }_{i,j}^N,i=1,2,\cdots ,N-1,$

and

$G_{0,j}=G_{N,j}=0,j=1,2,\cdots ,N-1,$

where

${\delta }_{i,j}^N=\{\begin{array}{l}\frac{1}{h_{i+1}}i=j,\\ 0otherwise.\end{array}$

We have $0\le G_{i,j}\le \frac{1}{\beta },i=0,\dots ,N;j=1,\dots ,N-1$, which can be proved by discrete comparison principle.

For any mesh function we have the representation

$u_i^N=\underset{j=1}{\overset{N-1}{{\displaystyle \sum }}}{h}_{j+1}{G}_{i,j}{L}^N{u}_i^N.$

It’s easy to prove that $L^{*,N}$ is an M-matrix. It follows that $G_{i,j}$ is an increasing function of i, for $i=1,2,\cdots ,j$. Similarly $G_{i,j}$ is a decreasing function of i, for $i=j,j+1,\cdots ,N$.

Thus, for each $j\in \left\{1,\dots ,N-1\right\}$ we have

$\underset{i=1}{\overset{N}{{\displaystyle \sum }}}\left|G_{i,j}-G_{i-1,j}\right|=2G_{j,j}\le \frac{2}{\beta }.$

We have

${\mathcal{L}}^N={\displaystyle {\sum }_{i=1}^N{h}_i\sqrt{1+{\left(D^{-}{u}_i^N\right)}^2}}\le {\displaystyle {\sum }_{i=1}^N{h}_i\left[1+\left|D^{-}{u}_i^N\right|\right]}=1+{\displaystyle {\sum }_{i=1}^N\left|u_i^N-u_{i-1}^N\right|},$ (10)

where ${\mathcal{L}}^N$ is the arc length of discrete solution $\left\{u_i^N\right\}$.

Set $v,w\in H^1(0,1)$, and $v_0=w_0$, $v_N=w_N$, $L^N\left(v_i-w_i\right)=0$, when

${\bar{b}}_i={\displaystyle {\int }_0^1{b}_u\left(x_i,w_i+s\left(v_i-w_i\right)\right)ds}$, ${\bar{c}}_i={\displaystyle {\int }_0^1{c}_u\left(x_i,w_i+s\left(v_i-w_i\right)\right)ds}$,

we have

$L^N\left(v_i-w_i\right)=T^N{v}_i-T^N{w}_i.$

For $i=1,\dots ,N-1$ one has

$L^N{u}_i^N=L^N\left(u_i^N-\stackrel{\to }{0}\right)=T^N{u}_i^N-T^N\stackrel{\to }{0}=-{\left(T^N\stackrel{\to }{0}\right)}_i,u_0^N=u_N^N=0.$ (11)

Thus, from (10) (11), we have

$\begin{array}{l}{\mathcal{L}}^N\le 1+\underset{i=1}{\overset{N}{{\displaystyle \sum }}}\underset{j=1}{\overset{N-1}{{\displaystyle \sum }}}{h}_{j+1}\left|-{\left(T^N\stackrel{\to }{0}\right)}_j\right|\cdot \left|G_{i,j}-G_{i-1,j}\right|\\ \le 1+{\Vert T^N\stackrel{\to }{0}\Vert }_{\infty }\underset{j=1}{\overset{N-1}{{\displaystyle \sum }}}{h}_{j+1}\underset{i=1}{\overset{N}{{\displaystyle \sum }}}\left|G_{i,j}-G_{i-1,j}\right|\le 1+\frac{2{\Vert T^N\stackrel{\to }{0}\Vert }_{\infty }}{\beta }.\end{array}$

Because the adaptive mesh is equally distributed according to the arc length of the solution, so the arc length from $\left(x_{i-1},u\left(x_{i-1}\right)\right)$ to $\left(x_i,u\left(x_i\right)\right)$ is ${\mathcal{L}}^N/N$.

So: $h_i=x_i-x_{i-1}\le {\mathcal{L}}^N/N\le C_1/N$.

Theorem 4: Set $\left\{u_i^N\right\}$ is the solution of (1) under the grid $\left\{x_i\right\}$ that satisfies arc-length equidistribution, then: ${\Vert u^N-u\Vert }_{\infty }\le C{N}^{-1}$.

Proof: Because $l_i=\sqrt{{\left(x_i-x_{i-1}\right)}^2+{\left(u_i^N-u_{i-1}^N\right)}^2}=h_i\sqrt{1+{\left(D^{-}{u}_i^N\right)}^2}$ and $\left\{u_i^N\right\}$ satisfies arc-length equidistribution, so: $l_i={\mathcal{L}}^N/N=h_i\sqrt{1+{\left(D^{-}{u}_i^N\right)}^2}$.

According to Lemma 3 and (9) (10), we have

${\Vert u^N-u\Vert }_{\infty }\le C\underset{{{}_{1\le i}}_{\le N}}{\max }{h}_i\sqrt{1+{\left(D^{-}{u}_i^N\right)}^2}\le C{N}^{-1}$.

4. Numerical Experiment

Example 1 (see [11]): We study the following quasilinear singularly perturbed two-point boundary value problem:

$\{\begin{array}{l}Tu\left(x\right):=-\epsilon u^{″}\left(x\right)-{\left(e^u\right)}^{\prime }+\frac{\pi }{2}\sin \frac{\pi x}{2}{e}^{2u}=0,forx\in \left(0,1\right),\\ u\left(0\right)=u\left(1\right)=0,\end{array}$

The exact solution of this problem is:

$u\left(x\right)=u_A\left(x\right)+O(\; \epsilon \; )$

where $u_A\left(x\right)=-\ln \left[\left(1+\cos \frac{\pi x}{2}\right)\left(1-\frac{1}{2}{e}^{-x/2\epsilon }\right)\right]$.

For any grid function $\left\{v_i\right\}$ we define the following discrete maximum norm:

$\Vert v\Vert :=\underset{i=1,2,\cdots ,N-1}{\max }\left|v_i\right|$

In order to verify the error estimate, the numerical convergence rate of the solution is calculated as follows:

$r_A^N={\log }_2\left(\frac{\Vert u_A-u^N\Vert }{\Vert u_A-u^{2N}\Vert }\right),$

where, we use the maximum node error $e:=\Vert u_A-u^N\Vert$ approximately replaces $\Vert u-u^N\Vert$, $e_0$ is the maximum nodal error of the numerical solution obtained from the initial mesh. Under the simple upwind difference scheme, the results of the adaptive method with the Bakhvalov-Shishkin grid as the initial grid are shown in Table 1.

From the numerical example, we can see that the adaptive mesh can converge faster and be more accurate when the number of grid points is small in Figure 1, and the convergence result reaches the first-order uniform convergence in Figure 2.

5. Conclusion

In this paper, we consider the general problem (1). A simple upwind difference scheme is used to discretize on adaptive meshs based on grid equidistribution of the arc-length monitor function. The uniform first-order convergence of the algorithm is proved, and numerical experiments verify our analysis.

Figure 1. The maximum node errors of numerical solutions in Example 1 with ε = 10⁻⁵.

Figure 2. Rate of Example 1 with ε = 10⁻⁵.

Table 1. The maximum node errors of numerical solutions in Example 1.

Acknowledgements

The authors thank the support of Natural Science Foundation of China (No. 11471019).

Conflicts of Interest

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

References

[1]	Tollmien, W., Schlichting, H., Grtler, H. and Riegels, F.W. (1961) über Flüssigkeitsbewegung bei sehr kleiner Reibung. Ludwig Prandtl Gesammelte Abhandlungen. Springer, Berlin, Heidelberg. https://doi.org/10.1007/978-3-662-11836-8
[2]	Roos, H.-G., Stynes, M. and Tobiska, L. (2008) Robust Numerical Methods for Singularly Perturbed Differential Equations. Springer-Verlag, Berlin, Heidelberg.
[3]	Kopteva, N.V. (1997) On the Uniform in Small Parameter Convergence of a Weighted Scheme for the One-Dimensional Time-Dependent Convection-Diffusion Equation. Computational Mathematics and Mathematical Physics, 37, 1173-1180.
[4]	Linβ, T. (2001) Uniform Pointwise Convergence of Finite Difference Schemes Using Grid Equidistribution. Computing, 66, 27-39. https://doi.org/10.1007/s006070170037
[5]	Kopteva, N. and Stynes, M. (2001) A Robust Adaptive Method for a Quasi-Linear One-Dimensional Convection-Diffusion Problem. SIAM Journal on Numerical Analysis, 39, 1446-1447. https://doi.org/10.1137/S003614290138471X
[6]	Jugal, M. and Srinivasan, N. (2010) Parameter-Uniform Numerical Method for Global Solution and Global Normalized Flux of Singularly Perturbed Boundary Value Problems Using Grid Equi-distribution. Computers & Mathematics with Applications, 60, 1924-1939. https://doi.org/10.1016/j.camwa.2010.07.026
[7]	Linβ, T., Roos, H.G. and Vulanovic, R. (2000) Uniform Pointwise Convergence on Shishkin-Type Meshes for Quasi-Linear Convection-Diffusion Problems. SIAM Journal on Numerical Analysis, 38, 897-912. https://doi.org/10.1137/S0036142999355957
[8]	Linβ, T. (2001) Sufficient Conditions for uniform Convergence on Layer-Adapted Grids. Applied Numerical Mathematics, 37, 241-255. https://doi.org/10.1016/S0168-9274(00)00043-X
[9]	Linβ, T. (2001) Uniform Second-Order Pointwise Convergence of a Finite Difference Discretization for a Quasilinear Problem. Computational Mathematics & Mathematical Physics, 41, 947-958.
[10]	Zheng, Q., Li, X. and Gao, Y. (2015) Uniformly Convergent Hybrid Schemes for Solutions and Derivatives in Quasilinear Singularly Perturbed BVPs. Applied Numerical Mathematics, 91, 46-59. https://doi.org/10.1016/j.apnum.2014.12.010
[11]	Vulanovic, R. and Hovhannisyan, G. (2006) A Posteriori Error Estimates for One- Dimensional Convection-Diffussion Problems. Computers & Mathematics with Applications, 51, 915-926. https://doi.org/10.1016/j.camwa.2005.11.027