Abstract
In this paper, an initial value method for solving a weakly coupled system of two second-order singularly perturbed Convection–diffusion problems exhibiting a boundary layer at one end is proposed. In this approach, the approximate solution for the given problem is obtained by solving, a coupled system of initial value problem (namely, the reduced system), and two decoupled initial value problems (namely, the layer correction problems), which are easily deduced from the given system of equations. Both the reduced system and the layer correction problems are independent of perturbation parameter, . These problems are then solved analytically and/or numerically, and those solutions are combined to give an approximate solution to the problem. Further, error estimates are derived and examples are provided to illustrate the method.
1. Introduction
Singular perturbation problems (SPPs) arise most frequently in diversified fields of applied mathematics. For instance, in fluid mechanics, elasticity, aerodynamics, quantum mechanics, chemical-reactor theory, oceanography, meteorology, modeling of semiconductor devices, and many others in the area. A well known fact is that, the solutions of such problems have a multiscale character, i.e., there are thin transition layer(s) where the solution varies very rapidly, while away from the layer(s) the solution behaves regularly and varies slowly. This leads to boundary and/or interior layer(s) in the solution of the problems. For a detailed discussion on the analytical and numerical treatment of SPPs one can refer the books of Miller et al. [1], O’Malley [2] and Roos et al. [3].
Due to the presence of the layer regions, it has been shown that the classical numerical methods fails to produce good approximations for singular perturbation problems (SPPs). In fact, some numerical techniques employed for solving second-order singularly perturbed boundary value problem (SPBVPs) are based on the idea of replacing this problems by suitable initial value problems (IVPs). The reason for that is, the numerical treatment of a boundary value problem is much more demanding than the treatment of the corresponding IVPs. There are different initial value methods in the literature of SPPs developed for solving SPBVPs, for the detail discussions of such methods one can refer the papers [4–8] and the references there in.
In the past few decades, a considerable amount of works have been reported in the literature of SPPs. However, most of the works connected with the computational aspects are confined to second-order equation. Only few results are reported for higher order and systems of equations. The systems of SPPs have applications in electro analytical chemistry, predator prey population dynamics, modeling of optimal control situations and resistance-capacitor electrical circuits [9]. In recent years, few scholars developed non-classical methods for different classes of systems of singularly perturbed differential equations. A class of systems of singularly perturbed reaction-diffusion equations have been examined by the authors in [10, 11–14]. In the papers [15, 16], a class of strongly coupled systems of singularly perturbed convection–diffusion equations are examined. The scholars in [17–21], considered weakly coupled systems of singularly perturbed convection–diffusion equations with equal or different diffusion parameters. A brief survey of article on the current progress about the numerical treatment of systems of singularly perturbed differential equations is also discussed in [22]. However, most of the methods developed for systems of singularly perturbed problems focus on fitted mesh method, so it is natural to look for an alternative approach for such problems.
In this paper, an initial value method for solving a weakly coupled system of two second-order singularly perturbed convection–diffusion equations exhibiting a boundary layer at one end is proposed. The technique, used in this work, is the careful factorization of original problem into a system of IVPs and two explicit IVPs which are independent of perturbation parameter. First, a system of IVPs is obtained by putting the perturbation parameter to zero, namely the reduced system, which corresponds to the outer solution. Next, using reduction of order together with stretching of variable gives two decoupled IVPs, namely the boundary layer correction problems, which corresponds to the inner solution. And then, the reduced system is solved numerically using fourth-order Runge–Kutta method, whereas, the boundary layer correction problems, are solved analytically. Finally, combining the above two solutions we obtain an approximation for the original problem. In addition, error estimates are derived and examples are provided to illustrate the method.
2. Statement of the Problem
Consider the problem of finding such that
with the boundary conditions
where , and , , , and are given constants, and is the singular perturbation parameter. The coefficient functions are taken to be sufficiently smooth on and satisfying the following conditions:
Under these assumptions the system (1a)–(2b) has a unique solution which exhibits a boundary layer of width on the left side of the underlying interval. The case for , can be put in to (1a) and (1b)by the change of independent variable from to .
The above coupled system of Equations (1a)–(2b) can also be written in vector form as
with the boundary conditions
where
Remark 1. In this paper, we consider only the case where there is one boundary layer at the left end of the interval. The case when the layer occurs at right end, can be analyzed similarly.
Notations. Let , the appropriate norm for studying the convergence of the approximate solution to the exact solution is the maximum norm . In case of vectors , we define
Throughout this paper, (sometimes sub-scripted) denotes generic positive constants independent of the singular perturbation parameter and in the case of discrete problems, also independent of the mesh parameter , these constants may assume different values but remains to be constant.
3. Analytic Results
In this section, a maximum principle, a stability result, and estimates of the derivatives of the system of Equations (1a)–(2b) are presented. First, we consider the following property of the operators and .
Lemma 2 (Maximum principle). Assume that is any sufficiently smooth function such that , and , in , then in .
Proof. Let and be arbitrary points in such that and . Without loss of generality, we assume that and suppose . Clearly and , . Moreover,which is a contradiction. It follows that our assumption is wrong, so that . Similarly, can be dealt. Hence, in .
An immediate consequence of the maximum principle is the following stability result.
Lemma 3 (Stability result). if , then for ,where .
Proof. Define two barrier functions as
where is the unit column vector and .
Clearly, and .
Also
Similarly, it can be proved that . Therefore, by the maximum principle, we obtain for all , which gives the required estimate. □
Now we give bounds on the derivatives of the exact solution for system (1a)–(2b).
Lemma 4. Let be the solution of (1a)–(2b), then and for ,
Proof. By following the approach used in the proof of Theorem 2 of [19] and the technique from [23], we can easily derive this Lemma. □
The solution of the problem (1a)–(2b) can be decomposed into smooth and singular components and respectively, as
where and . Further, the regular component can be written in the form of , where is the solution of the following system
and is the solution of the following system
and is the solution of the following systemwhere and are constants to be chosen such that , for . Thus the regular component is the solution of
and the singular component is the solution of
The following lemma provides the bound on the derivatives of the regular and singular components of the solution .
Lemma 5. For the smooth component and the singular component and their derivatives satisfy the bounds:respectively.
Proof. Using appropriate barrier functions, applying Lemma 2 and adopting the method of proof used in [[5] page 44], the present Lemma can be proved.
4. Description of the Method
In this section, we will obtain the solution of the (1a)–(2b) as a combination of two solutions: outer solution and inner solution.
4.1. Outer Solution
Let , for , be the solution of the reduced problem of (1a)–(2b) given by
Since, the degenerate equation does not satisfy the condition at , therefore, its contribution to the solution of (1a)–(2b) is for those values of which are away from . The problem (17) is therefore termed as outer problem.
For the exact solution of the reduced problem the following theorem gives an error bound.
Theorem 6. Let be the solution of (1a)–(2b) and be its reduced problem solution defined by (17). Thenwhere .
Proof. Consider the following barrier function , where
Its easy to see that further,
for an appropriate choice of , and
next, for the operator we have for an appropriate choice of . Similarly, we can prove that , for all . Therefore, from the maximum principle of Lemma 2, we obtain , and for . Hence the proof of the theorem.
Remark 7. From the above theorem, it is clear that the solution of problem (1a)–(2b) exhibits a strong boundary layer at and further away from the boundary layer region and in particular on , where , for sufficiently small values of , we have.
For the numerical solution of the reduced problem (17) we employ fourth-order Runge–Kutta method for system. Suppose be the numerical solution of the reduced problem obtained from fourth-order Runge–Kutta method, then the maximum error becomes
where is the equal mesh spacing of the domain of the problem.
4.2. Inner Solution
To obtain the inner solution for (1a)–(2b) we will use reduction of order together with streching of variable as follows:
First we rewrite the given problem (1a) for equivalently as
where
From Theorem 6 that the solution satisfies (1a)–(2b) on most part of the interval and away from . Thus by replacing the solution by on the right part of (25), we obtain an asymptotically equivalent approximation as:
where
Integrating both sides of (27) with respect to , gives
where
Using the reduced problem (16) in the above integral, yields
Then, substituting this in to (28), gives us
where and are integration constants. In order to determine ’s, we introduce the condition that the reduced equations of (32) should satisfy the boundary condition at . Thus, we get .
Hence, by substituting in (32), a first-order initial value problem which is asymptotically equivalent to the second-order system of boundary value problems (1a)–(2b) is obtained, and written as:
Next, to compute the solution for the layer part, a new inner variable is introduced by stretching the spatial coordinate , as
Using this stretching transformation in to (33), we obtain
In spite of the simplification, these equations remains first-order differential equation and also regularly perturbed. Thus, for , (24) becomes
These are differential equation for the solution of the layer regions. Moreover, the solutions of (36) are supposed to counter act for the fact that the solutions of the reduced problem do not satisfy the boundary condition at .
Further, using the substitutions, in to (36), we obtain the following boundary layer correction problems
Since these equations are separable linear initial value problems with constant coefficients which can easily be solved analytically, thus and gives
Finally, from standard singular perturbation theory it follows that the solution of the (33) admits the representation in terms of the solutions of the reduced problem (16) and boundary layer correction problem (38), which approximates the solution of the system (1a)–(2b); that is,
The numerical error of the present method has two sources: one from the asymptotic approximation of the modified problem (33) and the other from the numerical approximation of the reduced system (16). We can summarize the results of this section in the form of the following theorem.
Theorem 8. Let be the solution of the problem (1a)–(2b) given by (27) and be its approximate then
where is the equal mesh spacing for the domain of the problem.
Proof. Assume that be the numerical solution of the outer problem determined by making use of the Runge–Kutta method and be the solution of the inner problem whose left boundary conditions are affected by , such that is the approximation for the exact solution of (1a)–(2b) given by . Now using (17) and (27), we obtain,
for an appropriate choice of .
Therefore, we can conclude that , . □
5. Test Problems and Numerical Results
To demonstrate the efficiency and applicability of the proposed method we considered the following two test problems:
Example 5.1. Consider the following boundary value problems for the systems of convention–diffusion equations on :
The exact solution of this problem is
where,
, , and .
Since the exact solution is known, we calculate maximum point wise error by
where is the exact solution and is the numerical solution obtained by using mesh intervals. Tables 1 and 2 display, respectively, the maximum point-wise errors for and for different values of and . The plots of the exact and the approximate solution components for and are shown in Figure 1.
Example 10. Consider the following boundry value problems for the systems of convention–diffusion equations on(0,1): :
Since the exact solution is not known, we calculate maximum point wise error by using the double mesh principle defined by
where and denotes the and components of the numerical solutions obtained by using and meshes points, respectively. Tables 3 and 4 display, respectively, the maximum point-wise errors for and for different values of and . Figure 2 represents the numerical solutions of Example 10, for and .
6. Discussion
In this paper, an initial value method for solving a weakly coupled system of two linear second-order singularly perturbed convection–diffusion equations exhibiting a boundary layer at one end is proposed. The method is some how similar to the asymptotic expansion methods, but differs in detail. The approximate solution of the given problem is obtained by solving a coupled system of initial value problem (namely, the reduced system) and two decoupled initial value problems with constant coefficients (namely, the layer correction problems), which are easily deduced from the the original problem. Both the reduced system and the layer correction problems are independent of perturbation parameter, and therefore, we get easily the numerical solution by solving the reduced system using fourth-order Runge–Kutta method and solving the layer correction problems analytically. The method is simple to apply, very easy to implement on a computer and offers a relatively simple tool for obtaining approximate solution.
We have implemented the method on two test problems to illustrate the theoretical results, and presented the computational results for different values of and in Tables 1–4 and Figures 1 and 2. From the results it is observed that, for very small the present method approximates the exact solution of the problems very well.
Data Availability
No data were used to support this study.
Conflicts of Interest
The authors declare that there are no conflict of interest regarding to the publication of this paper.
Funding
This research work partially supported by Bahir Dar University, College of Science.