# A Simple Modification of Homotopy Perturbation Method for the Solution of Blasius Equation in Semi-Infinite Domains

Authors: M. Aghakhani; M. Suhatril; M. Mohammadhassani; M. Daie; A. Toghroli
Journal: Mathematical Problems in Engineering
 (2015)

---

## Abstract

A simple modification of the homotopy perturbation method is proposed for the solution of the Blasius equation with two different boundary conditions. Padé approximate is used to deal with the boundary condition at infinity. The results obtained from the analytical method are compared to Howarth’s numerical solution and fifth order Runge-Kutta Fehlberg method indicating a very good agreement. The proposed method is a simple and reliable modification of homotopy perturbation method, which does not require the existence of a small parameter, linearization of the equation, or computation of Adomian’s polynomials.


---

## Body

## 1. Introduction

Nonlinear differential equations are ubiquitous in science and technology. However, finding analytical solutions for this class of equations always has been a challenging task. In the recent years, several approximate methods were proposed for the analytical solution of nonlinear differential equations that do not depend on the existence of a small or large parameter in the equation. Among them, Homotopy Analysis Method (HAM) [1], Adomian Decomposition Method (ADM) [2], Variational Iteration Method (VIM) [3], Differential Transformation Method (DTM) [4], and Homotopy Perturbation Method (HPM) [5] can be mentioned. The homotopy perturbation method, first proposed by He [5], combines ease of implementation of perturbation methods with the flexibility of the homotopy analysis method. In addition, in the recent years, some novel methods for approximate solution of nonlinear differential equations emerged, such as Optimal Homotopy Asymptotic Method (OHAM) [6], Generalized Homotopy Method (GHM) [7], and (G′′/G)-expansion method [8].In the past, many scientists attempted to suggest an improvement to the homotopy perturbation method. Their studies were mainly focused on enlarging convergence radius and accelerating convergence of the solution and new suggestions for homotopy construction as well as alternation of linear and nonlinear part of homotopy based on the applied problem. Jafari and Aminataei [9] introduced a new treatment for homotopy perturbation method, which improves results from HPM. They discussed convergency of the proposed method. In order to demonstrate efficiency, accuracy, and superiority of the suggested method, the new modification was applied to some experiments in their paper. Yusufoğlu [10] purposed an alternation of HPM for exact solution of system of linear equations. He introduced an accelerating parameter to solve the system linear equation. Through embedding the accelerating parameters, convergence improved and required iterations reduced to one iteration. An effective method for convergence improvement of HPM for the solution of fractional differential equations was suggested by Hosseinnia et al. [11]. They proposed a method to select the linear part in the HPM to keep the inherent stability of fractional equations. To illustrate the improvement, Riccati fractional differential equations were solved. Results indicated the accuracy and effectiveness of the method compared to traditional HPM. A modification of HPM to solve fractional multidimensional diffusion equations is introduced by Kumar et al. [12]. In their study, Sumudu transform was utilized for transformation of partial differential equations into a new form. Dong et al. [13] obtained the solution of strongly nonlinear mixed Volterra-Fredholm integral equation by using an improvement to the homotopy perturbation method. Although the traditional HPM is divergent in these kinds of problems, the proposed method is convergent and leads to the exact solution.In the present paper, a simple and reliable modification of HPM is employed for the solution of two different forms of nonlinear Blasius equation in a semi-infinite domain. For that, we consider two forms of the Blasius equation arising in fluid flow inside the velocity boundary layer as follows.The first form of the Blasius equation is as follows:(1)d3dη3fη+12fηd2dη2fη=0,f0=0,f′0=0,f′η=1asη⟶∞.And the second form is as follows: (2)d3dη3fη+12fηd2dη2fη=0,f0=0,f′0=1,f′η=0asη⟶∞.It can be seen that two equations are the same except for boundary conditions. The first form of the equation is the well-known classical Blasius first derived by Blasius [14] and dates back about a century, which describes the velocity profile of two-dimensional viscous laminar flow over a finite flat plate. This form of the Blasius equation is the simplest form and the origin of all boundary layer equations in fluid mechanics. The second form of the equation, proposed more recently, arises in the steady free convection about a vertical flat plate embedded in a saturated porous medium [15], Laminar boundary layers at the interface of cocurrent parallel streams [16], or the flow near the leading edge of a very long, steadily operating conveyor belt [17]. Despite the fact that the Blasius equation is a simple third order nonlinear differential equation with smooth and monotonic solution, the hope that solution techniques employed to handle this equation can be extended to the difficult equations motivated several researchers to conduct abundant studies around it by means of different analytical and numerical approaches [18–21]. It is worth mentioning that despite that this problem is a century old, it is still being investigated by several authors, and some recent paper was published about it [22–26].In the recent years, the homotopy perturbation method is successfully applied to a wide variety of problems in science and technology. Analytical solution of the Blasius equation plays a very important role in the design and optimization of fluid devices. In the present paper, a simple modification of this analytical method was employed for the solution of the Blasius equation with two different forms of boundary conditions. For that, first, the equations are transferred to corresponding initial value problems with appropriate boundary conditions. Then, the initial value problems are solved by the homotopy perturbation method. Finally, the diagonal Padé transformation is used to handle the boundary condition at infinity and enlarge the convergence radius of the resulting series. For the convenience of the reader, we first introduce the homotopy perturbation method in the following section.


## 2. Homotopy Perturbation Method

The combination of the perturbation method and the homotopy method is called the homotopy perturbation method (HPM). This method lacks the limitations of the traditional perturbation methods and can take the full advantage of the traditional perturbation techniques. To illustrate the basic idea of this method [5], consider the following nonlinear differential equation:(3)Au-fr=0r∈Ω,subjected to the boundary conditions (4)Bu,∂u∂n=0r∈Γ,where A is a general differential operator, B is a boundary operator, fr is a known analytical function, and Γ is the boundary of the domain. ∂u/∂n denotes differentiation along the normal drawn outwards from Ω. The operator A, generally speaking, can be divided into two parts of L and N, where L is the linear part, while N is a nonlinear one. Hence, (3) can be rewritten as follows:(5)Lu+Nu-fr=0r∈Ω.We can construct the following homotopy which does not require that (3) has a small parameter:(6)Hν,p=1-pLν-Lu0+pAν-fr=0,in which p∈[0,1] is called homotopy parameter and u0 is an initial approximation of  (3) satisfying the given conditions. From (3), we have(7)Hν,0=Lν-Lu0=0,Hν,1=Aν-fr=0.Solution of (3) can be written as a power series in p as follows:(8)ν=ν0+pν1+p2ν2+⋯.Setting p=1 results in the approximate solution in the form of(9)V=limp→1⁡ν=ν0+ν1+ν2+⋯.The convergence of (9) has been proven by He [27].


## 3. Application of HPM-Padé to the First Form of the Blasius Equation

In this section, we consider the following form of the Blasius equation:(10)d3dη3fη+12fηd2dη2fη=0,f0=0,f′0=0,f′η=1asη⟶∞.Due to the boundary condition at infinity, direct application homotopy perturbation method is not appropriate for this equation and will not yield the desired results. Consequently, a modification in the equation or boundary condition is necessary. Equation (10) is reformulated as follows:(11)d3dη3fη+12fηd2dη2fη=0,f0=0,f′0=0,f′′0=σ, where σ is an unknown to be determined later. This modification of boundary conditions converts (10) from a boundary value problem to the corresponding initial value problem, which can be easily handled by HPM. This procedure is similar to the shooting method for numerical solution of boundary value problems in which the given boundary value problem is converted into an initial value problem.Following the outline given in Section2, a homotopy is constructed for (11):(12)Hf,p=1-pf′′′+pf′′′+12ff′′=0.It is assumed that (12) has a solution in the following form:(13)ν=f0+pf1+p2f2+p3f3⋯.Substituting f from (13) into (12) and equating terms with identical powers of p, we have(14)p0:f0′′=0,f00=0,f0′0=0,f0′′0=σ,p1:f1′′+12f0f0′′=0,f10=0,f1′0=0,f1′′0=0,p2:f2′′+12f0f1′′+f1f0′′=0,f20=0,f2′0=0,f2′′0=0,p3:f3′′+12f0f2′′+f1f1′′+f2f0′′=0,f30=0,f3′0=0,f3′′0=0,and ….Equations (14) can be rapidly solved, resulting in the following:(15)f0=12ση2,f1=-1240σ2η5,f2=114583333316800000000000σ3η8,f3=-232514881198000000000000σ4η11,and ….Therefore, according to (9), the solution of (12) reads(16)fη≔12ση2-1240σ2η5+114583333316800000000000σ3η8-232514881198000000000000σ4η11+1819996539100000000000000σ5η14+⋯.The unknown initial curvature in (16), σ, can be determined by imposing the boundary condition of (10) at infinity, that is, f′(∞)=1. For that, Padé approximants of (16) were formed. Padé approximants have the advantage of converting a function f(η) into a rational function in order to obtain more information about f(η). A [M/N] Padé approximate to function f(η) is the quotient of two polynomials, P(η)M and Q(η)N of degrees M and N, respectively. It is well known that if a function is free of singularities on the real axis, the Padé approximant will usually converge on the entire real axis [30]. Following the procedure suggested by Boyd [31], in order to determine the unknown parameter σ in (16), the boundary condition at infinity was imposed to the diagonal Padé approximant. Then, the roots of Padé approximant are used for the calculation of the unknown value. The result of this calculation is presented in Table 1. It can be seen from this table that the value obtained for σ=f′′(0) agrees very well with the numerical results. In Tables 2–4, f, f′, and f′′ obtained from 12th order HPM-Padé approximation are compared with Howarth’s numerical solution [28]. Furthermore, as it can be seen from Tables 2–4, purposed modification of HPM-Padé method is more accurate than the variational iteration method by He [29] and is valid for a wider range of the solution domain.Table 1
Root of the Padé approximates forσ=f′′(0), the first form of the Blasius equation.
Order
f
′
′
(
0
)
Error
[
2,2
]
0.33548
3.42
E
-
3
[
4,4
]
0.33267
6.13
E
-
4
[
6,6
]
0.33213
7.27
E
-
5
[
8,8
]
0.332052
5.34
E
-
6
[
10,10
]
0.3320549
2.44
E
-
6
Numerical [28]
0.33205734Table 2
Comparison betweenf(η) obtained from HPM-Padé with VIM and numerical method, first form of the Blasius equation.
η
Exact,Howarth [28]
VIM,He [29]
HPM-Padé(present)
Relative error (%)(VIM)
Relative error (%)(HPM-Padé)
0
0
0
0
0
0
1
0.16577
0.19319
0.16557
16.68
0.12
2
0.65003
0.6794
0.65001
4.51
0.003
3
1.39682
1.39106
1.39679
0.41
0.002
4
2.30576
2.24573
2.30572
2.6
0.001
5
3.28329
3.17748
3.28309
3.22
0.006
6
4.27964
4.14688
4.27767
3.10
0.04
7
5.27926
5.13359
5.26736
2.76
0.22Table 3
Comparison betweenf′(η) obtained from HPM-Padé with VIM and numerical method, first form of the Blasius equation.
η
Exact,Howarth [28]
VIM,He [29]
HPM-Padé(present)
Relative error (%)(VIM)
Relative error (%)(HPM-Padé)
0
0
0
0
0
0
1
0.32979
0.35064
0.32977
6.32
0.006
2
0.62977
0.61218
0.62976
2.79
0.001
3
0.84605
0.7964
0.84603
5.87
0.002
4
0.95552
0.90185
0.95551
5.62
0.001
5
0.99150
0.95523
0.99152
3.62
0.002
6
0.99868
0.98032
0.99883
1.84
0.015
7
0.99992
0.99158
0.99943
0.83
0.48Table 4
Comparison betweenf′′(η) obtained from HPM-Padé with VIM and numerical method, first form of the Blasius equation.
η
Exact,Howarth [28]
VIM,He [29]
HPM-Padé(present)
Relative error (%)(VIM)
Relative error (%)(HPM-Padé)
0
0.33206
0.5436
0.33205
63.7
0.003
1
0.32301
0.27141
0.32300
15.97
0.003
2
0.26675
0.22748
0.26675
14.72
0
3
0.16136
0.14117
0.16135
12.52
0.006
4
0.06424
0.07469
0.06422
16.27
0.031
5
0.01591
0.036
0.01586
126.27
0.314
6
0.0024
0.01645
0.00110
—
—
7
0.00022
0.00723
0.00060
—
—


## 4. Application of HPM-Padé to the Second Form of the Blasius Equation

The second form of the Blasius equation is as follows:(17)d3dη3fη+12fηd2dη2fη=0,f0=0,f′0=1,f′η=0asη⟶∞.Similar to the previous section, in order to solve (17), this equation was modified as follows:(18)d3dη3fη+12fηd2dη2fη=0,f0=0,f′0=1,f′′0=α.In this equation,α is unknown initial curvature of (17) and will be determined later by the application of the boundary condition of (17) at infinity. According to the instructions in Section 2, by embedding p∈[0,1], a homotopy was constructed for (17) as follows:(19)Hf,p=1-pf′′′+pf′′′+12ff′′=0,assuming that this equation has a solution in the following form:(20)ν=f0+pf1+p2f2+p3f3⋯.Substituting the above equation in (19) and rearranging will result in the following:(21)p0:f0′′=0,f00=0,f0′0=0,f0′′0=α,p1:f1′′+12f0f0′′=0,f10=0,f1′0=0,f1′′0=0,p2:f2′′+12f0f1′′+f1f0′′=0,f20=0,f2′0=0,f2′′0=0,p3:f3′′+12f0f2′′+f1f1′′+f2f0′′=0,f30=0,f3′0=0,f3′′0=0,and ….Ordinary differential equations (21) can be solved easily. The result is as follows:(22)f0=12αη2+η,f1=1240α2η5-148αη4,f2=11161280α3η8+1120160α2η7+1960αη6,f3=-1107520α25198α3η11+2518α2η10+439αη9+5η8,and ….Consequently, the solution of (18) can be written as follows:(23)fη≔η+12αη2-148αη4+1240α2η5+1960αη6+1120160α2η7+⋯.In order to determine the unknown initial curvatureα in (23), boundary conditions of (17) at infinity, f′(η)=0, must be applied. For that, Padé approximants of (23) which enlarge convergence radius of the solution were used. Then, α was determined from f′(η)=0 to the Padé approximants. Initial curvature of (23) obtained from this method is compared to fifth order Runge-Kutta Fehlberg numerical method in Table 5. It is worth mentioning that examining the behavior of (17) reveals that its initial curvature must be negative. Therefore, the negative root of Padé approximants is selected. In Tables 6–8, the result of 8th order HPM-Padé approximation is presented against that of exact (numerical) method. It can be seen that there is a good agreement between the results of the proposed method and numerical solution.Table 5
Root of the Padé approximates forα=f′′(0), the second form of the Blasius equation.
Order
f
′
′
(
0
)
Error
[
4,4
]
-
0.52270
7.8
E
-
2
[
6,6
]
-
0.48944
4.5
E
-
2
[
8,8
]
-
0.46407
2.0
E
-
2
[
12,12
]
-
0.44654
2.7
E
-
3
[
16,16
]
-
0.44565
1.9
E
-
4
[
20,20
]
-
0.44372
1.9
E
-
5
5th order Runge-Kutta Fehlberg (present)
-
0.44374Table 6
Comparison betweenf(η) obtained from HPM-Padé with the numerical method, second form of the Blasius equation.
η
Numerical(5th order Runge-Kutta Fehlberg)
HPM-Padé(present)
Relative error (%)(HPM-Padé)
0
0
0
0
1
0.786198
0.78620
0.00047
2
1.218546
1.21855
0.00053
3
1.432728
1.43273
0.0002
4
1.533086
1.53308
0.00022
5
1.578851
1.57884
0.00029
6
1.599437
1.59945
0.0011
7
1.612470
1.61280
0.002Table 7
Comparison betweenf′(η) obtained from HPM-Padé with the numerical method, second form of the Blasius equation.
η
Numerical(5th order Runge-Kutta Fehlberg)
HPM-Padé(present)
Relative error (%)(HPM-Padé)
0
1
1
0
1
0.587153
0.587153
0.000021
2
0.301784
0.301783
0.000056
3
0.144016
0.144016
0.00012
4
0.066244
0.066243
0.0013
5
0.029956
0.029949
0.00025
6
0.013469
0.013434
0.025
7
0.006119
0.006005
1.8Table 8
Comparison betweenf′′(η) obtained from HPM-Padé with the numerical method, second form of the Blasius equation.
η
Numerical(5th order Runge-Kutta Fehlberg)
HPM-Padé(present)
Relative error (%)(HPM-Padé)
0
-
0.443749
-
0.443748
0.00016
1
-
0.358313
-
0.358312
0.00064
2
-
0.214505
-
0.214505
0.00015
3
-
0.109834
-
0.109834
0.00017
4
-
0.052157
-
0.052159
0.0043
5
-
0.023906
-
0.023922
0.067
6
-
0.010736
-
0.010800
0.59
7
-
0.0466585
-
0.048415
3.7


## 5. Conclusion

In the present paper, a simple modification of the homotopy perturbation method is proposed for the solution of the Blasius equation in semi-infinite domains. The equation in semi-infinite domain is transferred into equivalent initial value problems, which results in appearance of an unknown coefficient. In order to determine the unknown coefficient, the boundary condition of the problem at infinity is imposed to Padé approximant of the solution. The results are in very good agreement with numerical and previous data available in the literature.



---

*Source: 671527-2015-11-08.xml*