Application of the perturbation iteration method to boundary layer type problems

The recently developed perturbation iteration method is applied to boundary layer type singular problems for the first time. As a preliminary work on the topic, the simplest algorithm of PIA(1,1) is employed in the calculations. Linear and nonlinear problems are solved to outline the basic ideas of the new solution technique. The inner and outer solutions are determined with the iteration algorithm and matched to construct a composite expansion valid within all parts of the domain. The solutions are contrasted with the available exact or numerical solutions. It is shown that the perturbation–iteration algorithm can be effectively used for solving boundary layer type problems.

The classification of such algorithms for algebraic equations turns out to be quite original. Motivated by the successful results, the systematic and algorithmic way of incorporating the methods were applied to first order ordinary differential equations in the pioneering work of Pakdemirli et al. (2011). The extension of the method to differential equations is non-trivial as the dependent variable and its derivatives are considered as separate variables inspired from the Lie Group theory applied to differential equations. The method is called as "the perturbation iteration method" to distinguish it from the past literature on the so-called "iteration-perturbation methods" (He 2001;Mickens 1987Mickens , 2005Mickens , 2006 which are not systematic approaches and does not produce general algorithms valid for various types of differential equations. The new perturbation-iteration method is directly applicable in a systematic algorithmic way, does not require special transformations or ad hoc assumptions. Second order differential equations were treated by the same method by Aksoy and Pakdemirli (2010) for Bratu type equations. Aksoy et al. (2012) further solved some nonlinear heat transfer equations. Dolapci et al. (2013) applied the method to Fredholm and Volterra integral equations. Şenol et al. (2013) treated the first order differential equation systems using the method. Pakdemirli (2013) reviewed the mentioned work. The basic idea behind the new method is to construct an approximate solution as in the classical perturbation methods, but then to iterate over this approximate solution to converge to the real solution of the problem. The iterations indeed enable one to converge to the real solutions for large perturbation parameters also. The mentioned work contains only applications to the regular perturbation problems which do not have secular terms in the analysis. In a very recent work, Pakdemirli (2015) successfully applied the method to problems with secularities also. The work considers the modifications of the method for non-regular problems containing secular terms which possess slow convergence problems and truncation leads to unphysical unbounded solutions.
So far, the method has not been applied to boundary layer type singular problems which remain an open area of research. The work presented here is an initial attempt to address the implementation of the method to such singular problems. It is well known that for such problems, the method of matched asymptotic expansions is the most reliable method. In the method of matched asymptotic expansions, a solution valid outside the boundary layer is constructed first, which is the regular solution of the original problem, named the outer solution. Then another solution inside the boundary layer region where abrupt changes in the solutions is observed is constructed, called the inner solution. The outer and inner solutions are matched over the overlapping region and a composite expansion valid throughout the whole domain is constructed as a final step.
In the previous work, the systematic way of constructing perturbation-iteration algorithms are named as PIA(n,m), n representing the number of correction terms in the perturbation expansion and m representing the order of highest derivatives in the Taylor series expansion. The simplest algorithm is the PIA(1,1) with one correction term in the perturbation expansion and first order derivatives in the Taylor series expansion. For consistency and solvable iteration equations, n should be equal or less than m (n ≤ m). As n or m increases, the convergence became faster but the algebraic complexity of the equations may end up with analytically unsolvable equations. In the previous work on differential equations, solutions are presented for PIA(1,1), PIA(1,2), PIA(1,3) and PIA(2,2) algorithms. In this pioneering and preliminary work on boundary layers, only the PIA(1,1) algorithms are applied to boundary layer type problems. Linear and nonlinear problems are treated with the PIA(1,1) method and solutions are contrasted with the available exact or numerical solutions of the original problem. For PIA(1,1), since the order of equations to be solved is less than the number of boundary conditions, the outer and inner solutions must be found by separate iterations and matched together similar to the method of matched asymptotic expansions.
Finally, as a last comment, the Lie Group methods of constructing exact solutions to differential equations cannot be applied to the problems considered here. The main reason is that the problems are defined over finite domains. For such domains, the symmetries of the equations do not remain stable after the boundary conditions are applied. To obtain group invariant solutions (similarity solutions) the equations as well as the boundary conditions and boundaries should remain invariant under the transformations especially for nonlinear problems (Bluman and Kumei 1989). For infinite or semi-infinite domains, this is more likely to occur but for finite domains, the symmetries are lost. For the singular problems with finite domains, the perturbation solutions, especially the method of matched asymptotic expansions is one of the most reliable analytical method.
Here, an alternative approach which is a modification of the perturbation-iteration method to handle such problems is developed.

Theory of the perturbation iteration algorithm for second order differential equations
In this section, PIA(1,1) algorithm is derived for second order differential equations. In accordance, only one correction term in the perturbation expansion and Taylor expansions up to first order derivatives are taken in constructing the iteration procedure. Consider a general second order differential equation, with y = y(x) and ε is the perturbation parameter. Note that a parameter such as ε may not exist at all in the original equation. In that case, the parameter should carefully be inserted into the appropriate part of the equation. The role of the perturbation parameter is to simplify the resulting iteration equations which cannot be solved otherwise. If the parameter is artificially introduced, it can be taken as 1 at the end. Only one correction term is taken in the perturbation expansion Upon substitution of (2) into (1) and expanding in a Taylor series with first order derivatives only yields where F y = ∂F ∂y , F y ′ = ∂F ∂y ′ , F y ′′ = ∂F ∂y ′′ , F ε = ∂F ∂ε and all derivatives are evaluated at ε = 0. It is readily observed that the above equation is a variable coefficient non-homogenous linear second order differential equation with respect to the unknown (y c ) n in its most general form. Starting with an initial guess y 0 , first (y c ) 0 is calculated from (3) and then substituted into (2) to calculate y 1 . The iteration procedure is repeated using (3) and (2) until a satisfactory result is obtained. In solving boundary layer problems, for outer solution, the original equation is iterated. However for constructing the inner solution, first the equation is expressed in terms of the boundary layer variable and then the iteration is performed over the transformed equation. The last step is to match the solutions and construct a composite expansion valid within all the domain. This last step is similar to the calculations of the method of matched asymptotic expansions. If the inner and outer solutions cannot be matched, possible reasons may be the wrong location of the boundary layer or a bad initial guess for the iterations.

Applications to linear differential equations
Two linear second order constant coefficient problems, one homogenous and the other non-homogenous are treated using the PIA(1,1) algorithm.

Example Problem 1
Consider the singular perturbation problem For the problem, an outer solution and an inner solution will be sought and both solutions will be matched to construct a composite solution.

The outer solution
To find the outer solution, the original equation is taken in the analysis Before starting, one of the important issues is the location of the boundary layer. This requires knowledge of the physics of the problem. For variable coefficient linear second order singular equations, a systematic way of locating the boundary condition indeed exists (Nayfeh 1981). For the problem considered, the boundary layer is located near x ≈ 0. Therefore, the outer solution is not expected to satisfy the condition at x = 0. Substituting (5) into (3) and re-arranging, the iteration equation is Note that due to the small term multiplying the highest derivative, the original equation yields a first order iteration equation which leads to inconsistencies if both conditions are forced to satisfy the solution. This justifies an outer and inner solutions to be constructed separately. Inspired by the condition at the right, a simple initial guess is suggested to start the iteration. Substituting this initial guess into (6) and solving for (y c ) 0 yields Hence, the first iteration result is and imposing the right hand side boundary condition yields Progressing in a similar way, the second iteration result is This result can also be found with the method of matched asymptotic expansions. However, with PIA(1,1), the same solution is retrieved from a very simple initial guess. This solution does not satisfy the boundary condition at the left hand side i.e. y 2 (0) � = 0.

The inner solution
The inner solution is valid within a small region called the boundary layer where the solution has a sharp turn. The first step is to determine the boundary layer variable by stretching the coordinate inside this layer Substituting the transformed variable into the original equation yields where Y = Y (ξ ). Balancing should be performed for the terms. There are two options: Balancing the first and second term or the first and third term. Balancing the second and third term is not an option because ν = 0 and the original equation is retrieved with no stretching.

First and second term
For this case 1 − 2ν = −ν ⇒ ν = 1. Equation (13) after multiplication with ε is The goal is to retain as much terms as one can at the leading order. In this case, two terms remain in the first order and the choice is admissible.

First and third term
Balancing the first and third term requires 1 − 2ν = 0 ⇒ ν = 1/2. Equation (13) is Since the middle term is the only leading dominant term, this choice is not an admissible choice and hence discarded. To summarize, the distinguished limit is ν = 1 and the boundary layer variable is and the transformed equation to be solved is Equation (3) takes the special form then Using (18) and (2) and applying the boundary condition at the left, i.e. Y(0) = 0, an iteration process can be constructed. Starting from a very simple guess The successive two iteration are Note that this solution is valid in the neighborhood of x = 0 and therefore is not expected to satisfy the boundary condition at the right hand side. The equations solved are second order while the conditions are less than the required value of 2. The undetermined coefficients will be determined by the matching conditions.

Matching
The inner and outer solutions should smoothly combine in the overlapping region which requires solutions to be matched. Employing the Van Dyke's (1975) matching principle in the overlapping region. In accordance, the outer expansion is written in terms of the inner variable, approximated and equated to the inner expansion written in terms of the outer variable and approximated. Hence and then approximation is taken up to two terms for fixed ξ Returning back to the original variable x now For the right hand side of the equation in the matching, the inner expansion is represented in terms of the outer variable x and approximated up to two terms for fixed x Note that for x kept fixed and not negligibly small, e −2x/ε is an exponentially small term which can be neglected. Equating (27) and (25) The second term is an exponentially small term which can be neglected. Under the approximations, the exact solution is Approximating further the exponential terms with the perturbation parameter the result is This approximation is comparable with the perturbation iteration solution if further is taken in (40).

Example Problem 2
Consider the non-homogenous linear boundary value problem

The outer solution
To find the outer solution, the original equation is used The boundary layer is located near x ≈ 0. Therefore, the outer solution is not expected to satisfy the condition at x = 0. Substituting (43)

into (3) and rearranging, the iteration equation is
Choose a simple initial guess compatible with the boundary condition at the right to start the iteration. The first iteration solution satisfying the boundary condition at the right is This solution is not expected to satisfy the boundary condition at the left hand side i.e. y 1 (0) � = 0.

The inner solution
The inner solution is valid within a small region called the boundary layer where the solution has a sharp turn. It turns out that the boundary layer variable by stretching the coordinate inside this layer is and the transformed equation is Note that this solution is valid in the neighborhood of x = 0 and therefore is not expected to satisfy the boundary condition at the right hand side. The undetermined coefficient will be determined by the matching conditions.

Matching
The matching of solutions requires in the overlapping region. In accordance, the outer expansion is written in terms of the inner variable, approximated and equated to the inner expansion written in terms of the outer variable and approximated. Hence The inner expansion is represented in terms of the outer variable x and approximated Note that for x kept fixed and not negligibly small, e −2x/ε is an exponentially small term which can be neglected. Equating (54) and (55) gives the undetermined constant

Hence the inner solution in terms of the original variable is
The composite solution is or which is a valid approximation throughout the whole domain. This solution satisfies both of the boundary conditions and can also be retrieved by the method of matched asymptotic expansions.The exact solution of the problem is Neglecting the exponentially small term e −2/ε and O(ε) terms, solution (59) is obtained, hence the solution is a valid approximation of the exact solution.

Applications to non-linear differential equations
A first order and a second order nonlinear boundary value problems are solved in this section using the method.

Example Problem 3
Consider the first order non-linear boundary value problem

For the outer solution, the original equation is used
The boundary layer is assumed to be located near x ≈ 0. If there arises any inconsistency in the solutions, the assumption may be revised. The outer solution is not expected to satisfy the condition at x = 0. Substituting (62) into (3) and re-arranging, the iteration equation is Since there is no boundary condition at the right hand side, the initial guess might be selected so that the iteration equation takes the simplest form. Hence n − 4 + 2y n ε(y c ) n + εy ′ n = 0, n = 0, 1, 2, . . .
(64) y 0 = 2 The first iteration solution is which is the same solution with the starting value. Therefore the outer expansion seems to be converged.

The inner solution
The boundary layer variable by stretching the coordinate inside this layer is and the transformed equation is  SpringerPlus (2016) 5:208 The exact solution is Solution (75) is a very good approximation of the exact solution (76). In Fig. 1, ε = 0.1, and in Fig. 2, ε = 100 are taken. For small or large perturbation parameters, the agreement is excellent.

Example Problem 4
Consider the second order non-linear boundary value problem (Cole 1968) The outer solution

For the outer solution, the original equation is used
The boundary layer is located near x ≈ 0. The outer solution is not expected to satisfy the condition at x = 0. Substituting (78) into (3) and rearranging, the iteration equation is y 0 = 3 which is inspired from the right hand side boundary condition would not be a good initial guess. Rather, one has to look for a trivial function which most simplifies the iteration equation for n = 0. Such an initial guess is The first and second iteration solutions satisfying the right hand side condition is which shows that the outer solution converged.

The inner solution
The boundary layer variable by stretching the coordinate inside this layer is and the transformed equation is The iteration Eq. (3) takes the form Starting from a simple guess which makes the most simplification in the iteration equation, the first iteration is Matching this solution with the first iteration of the outer solution makes a simplification and the constant c 1 turns out to be zero. Hence, the final first iteration is (78) F (y, y ′ , y ′′ , ε) = εy ′′ + yy ′ − y = 0 (79) y n y ′ n − y n + (y ′ n − 1)ε(y c ) n + y n ε(y ′ c ) n + εy ′′ n = 0, n = 0, 1, 2, . . .
(86) Y 0 = 2 (87) Y 1 = (2 + c 1 )(1 − e −2ξ ) + εξ and nonlinear problems are treated with the method and results are compared with the available analytical or numerical solutions. It is shown that the method can be effectively applied to singular problems producing admissible analytical approximate solutions. The author declares that he has no competing interests