Solving Systems of Singularly Perturbed Convection Diffusion Problems via Initial Value Method

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, . ese 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.


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.
is 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). e reason for that is, the numerical treatment of a boundary value problem is much more demanding than the treatment of the corresponding IVPs. ere 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][5][6][7][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. e 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][12][13][14].
In the papers [15,16], a class of strongly coupled systems of singularly perturbed convection-diffusion equations are examined. e scholars in [17][18][19][20][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. e 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.

Statement of the Problem
Consider the problem of finding 푦 1 , 푦 2 ∈ 퐶 0 Ω ∩ 퐶 2 (Ω) such that with the boundary conditions where Ω = (0, 1), Ω = [0, 1] and 푦 = 푦 1 , 푦 2 푇 , 1 , 2 , 1 and 2 are given constants, and 0 < 휀 ≪ 1 is the singular perturbation parameter. e 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 le side (푥 = 0) of the underlying interval. e case 푎 (푥) ≤ 훼 < 0 for 푖 = 1, 2, can be put in to (1a) and (1b) by the change of independent variable from to 1 − 푥. e 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 le end of the interval. e case when the layer occurs at right end, can be analyzed similarly.
Notations. Let 푦 : 퐷 → R, the appropriate norm for studying the convergence of the approximate solution to the exact solution is the maximum norm , we define roughout 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.

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 1 and 2 .
Lemma 2 (Maximum principle). Assume that 휋(푥) is any sufficiently smooth function such that 휋(0) ≥ 0, 휋(1) ≥ 0 and Proof. Let * and * be arbitrary points in (0, 1) such that and 휋 2 푦 * = min 푥∈Ω 휋 2 (푥) . Without loss of generality, we assume that which is a contradiction. It follows that our assumption An immediate consequence of the maximum principle is the following stability result.
and 푢 1 = 푢 11 , 푢 12 푇 is the solution of the following system and 푢 2 = 푢 21 , 푢 22 푇 is the solution of the following system where 1 and 2 are constants to be chosen such that  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.

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.
be the solution of the reduced problem of (1a)-(2b) given by Since, the degenerate equation does not satisfy the condition at 푥 = 0, therefore, its contribution to the solution of (1a)-(2b) is for those values of which are away from 푥 = 0. e problem (17) is therefore termed as outer problem.
For the exact solution of the reduced problem the following theorem gives an error bound.
Its easy to see that 휑 ± 푖 ∈ 퐶 1 (Ω) ∩ 퐶 0 Ω further, for an appropriate choice of 1 , and next, for the operator 1 we have Theorem 6. Let 푦(푥) be the solution of (1a)-(2b) and 푢 0 (푥) be its reduced problem solution defined by (17). en Proof. Consider the following barrier function for an appropriate choice of 2 . Similarly, we can prove that 퐿 2 휑 ± (푥) ≥ 0, for all 푥 ∈ Ω. erefore, from the maximum principle of Lemma 2, we obtain 휑 ± 푖 (푥) ≥ 0, ∀푥 ∈ Ω and for 푖 = 1, 2. Hence the proof of the theorem. For the numerical solution of the reduced problem (17) we employ fourth-order Runge-Kutta method for system. Suppose 푈 0 = 푈 01 , 푈 02 be the numerical solution of the reduced problem obtained from fourth-order Runge-Kutta method, then the maximum error becomes where ℎ = 1/푁 is the equal mesh spacing of the domain of the problem.

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 푖 = 1, 2 equivalently as ᐈ ᐈ ᐈ ᐈ ≤ 퐶휀, ∀푥 ∈ Ω and for 푖 = 1, 2. (24) From eorem 6 that the solution 푢 0 (푥) satisfies (1a)-(2b) on most part of the interval [0, 1] and away from 푥 = 0. us by replacing the solution 푦(푥) by 푢 0 (푥) 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 en, substituting this in to (28), gives us Proof. Assume that 0 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 le boundary conditions are affected by 0 , such that 푌 = 푈 0 + 푊 is the approximation for the exact solution of (1a)-(2b) given by 푦 = 푢 0 + 푤 + 퐶 1 휀. Now using (17) and (27), we obtain, for an appropriate choice of . erefore, we can conclude that

Test Problems and Numerical Results
To demonstrate the efficiency and applicability of the proposed method we considered the following two test problems: Example 9. Consider the following boundary value problems for the systems of convention-diffusion equations on (0, 1): e exact solution of this problem is where 1 and 2 are integration constants. In order to determine 's, we introduce the condition that the reduced equations of (32) should satisfy the boundary condition at 푥 = 1. us, we get 푘 1 = 푘 2 = 0. Hence, by substituting 푘 1 = 푘 2 = 0 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. us, for 휀 = 0, (24) becomes ese 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 푥 = 0.
Further, using the substitutions, 푊 푖 (푡) = 푤 푖 (푡) − 푢 0푖 (0) 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, e 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.
Journal of Applied Mathematics 6 where 푁 푖,푗 and 2푁 푖,푗 denotes the 푡ℎ and 2푖 푡ℎ components of the numerical solutions obtained by using and 2푁 meshes points, respectively. Tables 3 and 4 display, respectively, the maximum point-wise errors for 1 and 2 for different values of and . Figure 2 represents the numerical solutions of Example 10, for 푁 = 1024 and 휀 = 0.05.

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. e method is some how similar to the asymptotic expansion methods, but differs in detail. e 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. e 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.