A subdivision-based approach for singularly perturbed boundary value problem

A numerical approach for solving second order singularly perturbed boundary value problems (SPBVPs) is introduced in this paper. This approach is based on the basis function of a 6-point interpolatory subdivision scheme. The numerical results along with the convergence, comparison and error estimation of the proposed approach are also presented.


Introduction
The singularly perturbed boundary value problems (SPBVPs) frequently occur in the different areas of physical phenomena. Specifically, these occur in the fields of fluid dynamic, elasticity, neurobiology, quantum mechanics, oceanography, and reactor diffusion process. These problems often have sharp boundary layers. These boundary layers usually appear as a multiple of the highest derivative. Their small values cause trouble in different numerical schemes for the solution of SPBVPs. Therefore, it is important to find the numerical and analytic solutions of these types of problems. The different second order SPBVPs have different expressions but we deal with the following: where 0 < ε < 1, while p(u), q(u), g(u) are bounded and real valued functions. g(u), α 0 , α 1 depends on ε. We may refer to Ascher et al. [1] for more details as regards such a type of SPBVPs.
Here, we first present short review of different methods for the solution of second order SPBVPs; then we discuss a subdivision-based solution of SPBVPs.
The second order SPBVPs were solved based on cubic spline scheme by Aziz and Khan [2,6] in 2005. These problems were also solved by Bawa and Natesan [3] in the same year.
They have used quintic spline based approximating schemes. Kadalbajoo and Aggarwal [5] and Tirmizi et al. [17] solved self-adjoint SPBVPs by using B-spline collocation and nonpolynomial spline function schemes, respectively. Kumar and Mehra [7] and Pandit and Kumar [14] solved SPBVPs by a wavelet optimized difference and uniform Haar wavelet methods, respectively.
The second order SPBVPs were also solved by [9,12,13]. They have used finite difference scheme for the solution.
The linear SPBVPs were solved by [4,11,15]. They have used interpolating subdivision schemes for this purpose. The solution of second order SPBVPs by subdivision techniques did not reported yet. We develop an algorithm by using a 6-point interpolating subdivision scheme (6PISS) [8]. We have where the scheme is C 2 -continuous for 0 < μ < 0.042. It has support width (-5, 5). It has fourth order of approximation. It satisfies the 2-scale relation where Here is the layout of the rest of the work. In Sect. 2, we first find the derivatives of ρ(u) then by using them we develop the collocation algorithm. The convergence of the method is discussed in Sect. 3. In Sect. 4, we present the numerical solutions of different problems. The comparison of the solutions obtained by different methods is also offered. Section 5 deals with our conclusion.

The numerical algorithm
In this section, we develop an algorithm to deal with second order SPBVPs. First we discuss the derivatives of 2-scale relations known as basis functions of the subdivision scheme.

Derivatives of 2-scale relations
The 6PISS is C 2 -continuous by [8], so its 2-scale relations ρ(u) are also C 2 -continuous. First we find the eigenvectors (both left and right) of the subdivision matrix of 6PISS then we find the derivatives of ρ(u). For simplicity, we choose μ = 0.04 to find eigenvectors. We use a similar approach to [4,11] to find the derivatives. The first two derivatives of the 2-scale relation are given in Table 1.

The 6PISS based algorithm
Let m be the indexing parameter which might be equal to or greater than the last right end integral value of the right eigenvector corresponding to the eigenvalue 1 2 of the subdivision matrix for (3). Some useful notations depending on the indexing parameter m are defined as h = 1/m and υ κ 1 = κ 1 /m = ih with κ 1 = 0, 1, 2, . . . , m. Finally, we suppose that the approximate solution of (1) is where {d κ 1 } are the unknowns to be determined; then with given boundary conditions at both ends of the interval From (6), we have We get the system of equations by using (6) and (8) in (7), where A κ = p(υ κ ), q κ = q(υ κ ) and g κ = g(υ κ ). This implies This further implies m+4 where κ = 0, 1, 2, . . . , m and υ κ 1 = ih or υ κ = jh. By using the notation "ρ(κ 1 ) = ρ κ 1 ", (9) can be written as m+4 As we observe from Table 1, The above system of equations is summarized in the following proposition.
The column matrices D and G 1 are given by and The system (14) in its present form does not have a unique solution. We need eight extra equations to get its unique solution. Luckily, two equations can be obtained from (2), i.e., D(0) and D(1) and for the remaining six equations, we move to the next section.

End point constraints
If the data points are given then the 6PISS is suitable to fit the data with a fourth order of approximation. So we use the fourth order polynomials to define the constraints at the end points. Here we suggest two types of polynomials i.e. the simple cubic polynomial (i.e. a polynomial of order 4) and cardinal basis function-based cubic polynomial, to get the constraints.

C-1: Constraints by polynomial of degree three:
We use the fourth order polynomial Since by (6), D(υ κ 1 ) = d κ 1 for κ 1 = 1, 2, 3 then, by replacing υ κ 1 by -υ κ 1 , we have Hence, we get the following three constraints defined at the left end points: A similar procedure is adopted for the right end i.e. we can define d κ 1 = C 1 (υ κ 1 ), κ 1 = m + 1, m + 2, m + 3 and So the following three constraints are defined at the right end: C-2: Constraints by cardinal basis functions: The following fourth order polynomial C 2 (υ) can be used to find the left points d -3 , d -2 , d -1 : while the basis functions are given by and for t = 0, 1 A similar procedure is adopted for the right end points and and for t = m, m + 1
By C-1 constrants: If we use (2), (12), (18) and (19) then the system can expressed as where S 1 = (S T L 1 , S T , S T R 1 ) T , S is defined by (15). The matrix [S L 1 ] 4×(m+9) is defined as its first three rows and the fourth row are obtained from (18) and (2), respectively, its first row and the last three rows are obtained from (2) and (19), respectively, while the matrices D and G 1 are defined in (16) and (17), respectively.

Existence of the solution
The matrices S 1 and S 2 involved in the systems (22) and (24) are non-singular. Their nonsingularity can be checked by finding their eigenvalues. We notice that for m ≤ 500 the eigenvalues are nonzero. By [16], these are non-singular. Their singularity is not guaranteed for m > 500.

Error estimation of the algorithm
This section discussed the mathematical results as regards the convergence of the proposed method.

Solutions of second order SPBVPs and discussions
In this section, we consider second order SPBVPs and find their numerical solutions by using different algorithms. Since we have developed two linear systems i.e. (22) and (24) for approximate solutions of the SPBVPs, both systems have been used for solutions. We also give a comparison of solutions by computing the maximum absolute errors of the analytic and approximate solutions.
Example 4.1 This type of problem has also solved by [2,3,5,6,14], where the boundary conditions of the above problem are its analytic solution is cos 2 (πυ).        where the boundary conditions of the above problem are its analytic solution is Example 4.3 Take the problem already solved by [7,9], here 0 ≤ υ ≤ 1 and its analytic solution is

Discussion and comparison
We solve SPBVPs by our algorithm and summarized the results in the following form.  Tables 7 and 10 show the MAE while Tables 9 and 10 present the comparison of MAE with [10,12,13]. This shows our results are better. Figure 4 shows the solutions while  Tables 11 and 12. We compare our result with the results of [7,9] and found them to be better. The graphical representation is given in Figs. 8 and 9 • From these results we conclude that the condition C-2 gives better results than the condition C-1. • If we keep m fixed, then MAE increases with the increase of ε. It is also observed that if we keep ε fixed, then MAE decreases with the increase of m.

Conclusions
In this paper, we introduced a numerical algorithm for the solution of second order SPBVPs. The algorithm was developed by using the 2-scale relation of a well-known interpolating subdivision scheme. This algorithm gives the approximate solution of second order SPBVPs with a fourth order of approximation. We presented the comparison of maximum absolute error of the solutions obtained from subdivision (i.e. our method), spline [2,3,5,6], finite difference [9,10,12,13] and Haar wavelet [7,14] algorithms. We concluded that our algorithm gives smaller maximum absolute error.