# An Iteration Method Generating Analytical Solutions for Blasius Problem

Authors: Beong In Yun
Journal: Journal of Applied Mathematics
 (2011)

---

## Abstract

We derive a new iteration method for finding solution of the generalized Blasius problem. This method results in the analytical series solutions which are consistent with the existing series solutions for some special cases.


---

## Body

## 1. Introduction

We consider the generalized Blasius’ equation(1.1)y′′′(x)+αy(x)y′′(x)=0,0≤x<∞,
where α=1/2 or α=1, with boundary conditions
(1.2)y(0)=0,y′(0)=μ,y′(∞)=1.
This problem describes the boundary layer flow over a moving plate with constant velocity μ. For a special case of α=1/2 and μ=0, the series solution of the Blasius problem becomes
(1.3)S(x)=κ2x2-κ2240x5+11161280κ3x8-54257792κ4x11+⋯,
where κ=y′′(0)≈0.3320573362. The Blasius series, however, converges for |x|<ρ=5.6900380545. In the literature [1–3], it was shown that the limitation can be overcome by Padé approximants or an Euler-accelerated series.Lots of analytical methods such as Adomian decomposition methods [4–6], variational iteration methods [7–11], and homotopy analysis methods [12–14] have been proposed.


## 2. Derivation of an Iteration Formula

We develop a new iteration method to find the analytical series solution of the Blasius problem (1.1) subject to the boundary condition
(2.1)y(0)=0,y′(0)=μ,y′′(0)=κ,
where the curvature κ of the solution is assumed to be known. It should be noted that in order to make the problem easy to be solved, we consider the one point boundary conditions in (2.1) instead of the two-point boundary conditions in (1.2).First, fory=y(x) the Blasius equation (1.1) becomes
(2.2)(y′′+αyy′)′=α(y′)2.
From the boundary conditions in (2.1), it follows that
(2.3)y′′+αyy′=α∫0x{y′(t)}2dt+κ:=A(x,y).
This can be represented by
(2.4)(y′+α2y2)′=A(x,y),
which implies
(2.5)y′=∫0xA(t,y)dt-α2y2+μ:=B(x,y).
In the result, we have
(2.6)y(x)=∫0xB(t,y)dt.If we denote byyn the nth iterate solution and substitute it into the right hand side of (2.6), we have an iteration formula
(2.7)yn+1(x)=∫0xB(t,yn(t))dt,n≥0,
where
(2.8)B(x,yn(x))=∫0xA(t,yn(t))dt-α2yn(x)2+μ.
From (2.3) and (2.7), the function A(x,yn) can be represented by
(2.9)A(x,yn(x))=α∫0x{B(t,yn-1(t))}2dt+κ,
for n≥1 with A(t,y0)≡κ. Referring to the boundary conditions in (2.1), we may take the initial solution as
(2.10)y0(x)=κ2x2+μx.
The proposed method can be summarized by the following algorithm.Algorithm A.
We have the following steps.Step 1.
Set initial guesses(2.11)y(x):=y0(x),B(x,y):=y0′(x).Step 2.
For a large integerL>0, perform the iteration (2.7)–(2.9) using symbolic computations
(2.12)l:=0,while(l:=l+1)≤L,A(x,y):=α∫0x{B(t,y(t))}2dt+κ,B(x,y):=∫0xA(t,y(t))dt-α2y(x)2+μ,y(x):=∫0xB(t,y(t))dt.
It should be noted that by performing this algorithm, we can also obtain the approximatesB(x,yn) to the velocity y′(x)=B(x,y).


## 3. Analytical Solutions

Performing the above algorithm by using the symbolic calculation softwareMathematica, we have the successive approximate solutions below
(3.1)y1(x)=κx22-αk2x5120+(x-αkx424)μ,y2(x)=κx22-αk2x5120+11α2k3x840320-α3k4x11712800+(x-αkx424+11α2k2x75040-α3k3x1064800)μ+(α2kx6240-α3k2x924192)μ2,y3(x)=κx22-αk2x5120+11α2k3x840320-5α3k4x11532224+10033α4k5x1487178291200-5449α5k6x173908653056000+83α6k7x208935557120000-α7k8x2349080898944000+(x-αkx424+11α2k2x75040-5α3k3x1048384+10033α4k4x136227020800-5449α5k5x16229920768000+83α6k6x19446777856000-α7k7x222133952128000)μ+(α2kx6240-43α3k2x9120960+1157α4k3x12159667200-1147α5k4x157925299200+14057α6k5x1810260214272000-197α7k6x2149145564160000)μ2-(α3kx82688-23α4k2x112217600+967α5k3x142641766400-1829α6k4x17414554112000+α7k5x2066189312000)μ3-(α5k2x133294720-17α6k3x163251404800+α7k4x1947259168768)μ4.For the case ofα=1/2, we have(3.2)y1(x)=κx22-k2x5240+(x-kx448)μ,y2(x)=κx22-k2x5240+11k3x8161280-k4x115702400+(x-kx448+11k2x720160-k3x10518400)μ+(kx6960-k2x9193536)μ2,y3(x)=κx22-k2x5240+11k3x8161280-5k4x114257792+10033k5x141394852659200-5449k6x17125076897792000+83k7x20571875655680000-k8x236282355064832000+(x-kx448+11k2x720160-5k3x10387072+10033k4x1399632332800-5449k5x167357464576000+83k6x1928593782784000-k7x22273145872384000)μ+(kx6960-43k2x9967680+1157k3x122554675200-1147k4x15253609574400+14057k5x18656653713408000-197k6x216290632212480000)μ2-(kx821504-23k2x1135481600+967k3x1484536524800-1829k4x1726531463168000+k5x208472231936000)μ3-(k2x13105431040-17k3x16208089907200+k4x196049173602304)μ4.
One can see that the result is consistent with the known series solution [11, 12].In particular, whenα=1 and μ=0, it follows that
(3.3)y1(x)=κx22-k2x5120,y2(x)=κx22-k2x5120+11k3x840320-k4x11712800,y3(x)=κx22-k2x5120+11k3x840320-5k4x11532224+10033k5x1487178291200-5449k6x173908653056000+83k7x208935557120000-k8x2349080898944000,y4(x)=κx22-k2x5120+11k3x840320-5k4x11532224+9299k5x1429059430400-2173649k6x17355687428096000+13722337k7x20115852476579840000-27184438601k8x2312926008369442488320000+12320831753849k9x26403291461126605635584000000+⋯-k16x47463172433275878342410240000000.
In this case, κ≈0.4695999883. For comparison, we refer to another analytical solution obtained by the Adomian decomposition method as follows:
(3.4)y1A(x)=κx22-k2x5120,y2A(x)=κx22-k2x5120+11k3x840320,y3A(x)=κx22-k2x5120+11k3x840320-5k4x11532224,y4A(x)=κx22-k2x5120+11k3x840320-5k4x11532224+9299k5x1429059430400,y5A(x)=κx22-k2x5120+11k3x840320-5k4x11532224+9299k5x1429059430400-1272379k6x17118562476032000.
This solution is based on
(3.5)ynA(x)=y0(x)+∑k=1nuk(x),n=1,2,…,
with
(3.6)uk+1(x)=-L-1(Ak(x)),k=0,1,2,…,u0(x)=y0(x), and Adomian polynomial Ak(x) generated by the formula [6]
(3.7)Ak(x)=1k![dkdλk(∑j=0kλjuj(x))(∑j=0kλjuj′′(x))]λ=0,
where L-1 is an inverse operator of L=d3/dx3. Comparing the formulas in (3.3) and (3.4), one can see that the presented analytical solution yn(x) has more terms than ynA(x) in each nth iteration. In other words,
(3.8)yn(x)=ynA(x)+∑j=12n-n-1djx3(n+j)+2,
for any integer n≥2. In practice, Figure 1 depicts that the presented solutions yn(x) and their derivatives yn′(x), n=1,3,5 approximate exact ones better than ynA(x) and (ynA)′(x). Therein, we chose the initial solution y0(x)=(κ/2)x2 as given in (2.10) and took a numerical solution for the exact solution which is denoted by y*(x). Moreover, Table 1 includes numerical results of the errors and the CUP times spent in computations for the presented solution yn(x) compared with those of ynA(x). The L∞ error indicates the maximum error for the 50 nodes selected in the interval (0,ρ), where ρ is a radius of convergence of the series solution given in the literature [2, 15]. In fact, ρ≈4.02 for α=1 and ρ≈5.69 for α=1/2. The L2 error means ∥yn-y*∥2 over the same interval.Table 1
Comparison of the numerical results obtained by the presented method and the Adomian’s decomposition method forα=1,1/2 with μ=0.
nPresented solution,yn(x)Adomian’s solution,ynA(x)L∞ errorL2 errorCPU time (sec)L∞ errorL2 errorCPU time (sec)20.6900.3180.140.9950.4450.09α=140.3000.1101.031.0800.3820.3760.0780.0256.371.1100.3351.4880.0130.00426.51.1180.2995.7420.9760.5340.131.4070.7480.11α=1/240.4240.1851.011.5270.6430.4160.1110.0416.711.5700.5641.5580.0190.00627.41.5810.5045.93Graphs of the presented solutionsyn(x) and their derivatives yn′(x) in (a) and those of Adomian’s solutions ynA(x) and their derivatives (ynA)′(x) in (b).
(a)
Presented solutions(b)
Adomian’s solutionsBy the numerical performance, we can surmise the convergence of the algorithm proposed in this work, and the rate of convergence is better than that of the Adomian decomposition method though it spends more CPU time as shown above. In addition, for example, for the case ofα=1 and μ=0, we may guess that the presented method has the same radius of convergence, ρ=4.0234644935 as the well known Blasius’ series [2] as follows:
(3.9)y(x)=∑k=0∞(-1)kpkκk+1(3k+2)!x3k+2,
where p0=1 and
(3.10)pk=∑j=0k-1(3k-13j)pjpk-j-1,k≥1.
Theoretical convergence analysis with extended application of the presented method to more general problems is left for further works.



---

*Source: 925649-2011-08-13.xml*