An Iteration Method Generating Analytical Solutions for Blasius Problem

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.

1. Introduction

We consider the generalized Blasius’ equation 𝑦^′′′(𝑥)+𝛼𝑦(𝑥)𝑦^′′(𝑥)=0,0≤𝑥<∞, (1.1) where 𝛼=1/2 or 𝛼=1, with boundary conditions 𝑦(0)=0,𝑦^′(0)=𝜇,𝑦^′(∞)=1. (1.2) 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 𝑆(𝑥)=𝜅2𝑥²−𝜅²𝑥240⁵+11𝜅161280³𝑥⁸−5𝜅4257792⁴𝑥¹¹+⋯, (1.3) where 𝜅=𝑦^′′(0)≈0.3320573362. The Blasius series, however, converges for |𝑥|<𝜌=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 𝑦(0)=0,𝑦^′(0)=𝜇,𝑦^′′(0)=𝜅, (2.1) 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, for 𝑦=𝑦(𝑥) the Blasius equation (1.1) becomes 𝑦^′′+𝛼𝑦𝑦^′′=𝛼𝑦′². (2.2) From the boundary conditions in (2.1), it follows that 𝑦^′′+𝛼𝑦𝑦′=𝛼^𝑥₀𝑦^′(𝑡)²𝑑𝑡+𝜅∶=𝐴𝑥,𝑦. (2.3) This can be represented by 𝛼𝑦′+2𝑦^2′=𝐴𝑥,𝑦, (2.4) which implies 𝑦′=^𝑥₀𝐴𝛼𝑡,𝑦𝑑𝑡−2𝑦²+𝜇∶=𝐵𝑥,𝑦. (2.5) In the result, we have 𝑦(𝑥)=^𝑥₀𝐵𝑡,𝑦𝑑𝑡. (2.6)

If we denote by 𝑦_𝑛 the 𝑛th iterate solution and substitute it into the right hand side of (2.6), we have an iteration formula 𝑦_𝑛+1(𝑥)=^𝑥₀𝐵𝑡,𝑦_𝑛(𝑡)𝑑𝑡,𝑛≥0, (2.7) where 𝐵𝑥,𝑦_𝑛(𝑥)=^𝑥₀𝐴𝑡,𝑦_𝑛(𝑡)𝛼𝑑𝑡−2𝑦_𝑛(𝑥)²+𝜇. (2.8) From (2.3) and (2.7), the function 𝐴(𝑥,𝑦_𝑛) can be represented by 𝐴𝑥,𝑦_𝑛(𝑥)=𝛼^𝑥₀𝐵𝑡,𝑦_𝑛−1(𝑡)²𝑑𝑡+𝜅, (2.9) for 𝑛≥1 with 𝐴(𝑡,𝑦₀)≡𝜅. Referring to the boundary conditions in (2.1), we may take the initial solution as 𝑦₀(𝑥)=𝜅2𝑥²+𝜇𝑥. (2.10) The proposed method can be summarized by the following algorithm.

Algorithm A. We have the following steps.Step 1. Set initial guesses 𝑦(𝑥)∶=𝑦₀(𝑥),𝐵𝑥,𝑦∶=𝑦^′₀(𝑥). (2.11)Step 2. For a large integer 𝐿>0, perform the iteration (2.7)–(2.9) using symbolic computations ()𝐴𝑙∶=0,while𝑙∶=𝑙+1≤𝐿,𝑥,𝑦∶=𝛼^𝑥₀𝐵𝑡,𝑦(𝑡)²𝐵𝑑𝑡+𝜅,𝑥,𝑦∶=^𝑥₀𝐴(𝑡)𝛼𝑡,𝑦𝑑𝑡−2𝑦(𝑥)²𝑦(𝑥)+𝜇,∶=^𝑥₀𝐵𝑡,𝑦(𝑡)𝑑𝑡. (2.12)
It should be noted that by performing this algorithm, we can also obtain the approximates 𝐵(𝑥,𝑦_𝑛) to the velocity 𝑦′(𝑥)=𝐵(𝑥,𝑦).

3. Analytical Solutions

Performing the above algorithm by using the symbolic calculation software Mathematica, we have the successive approximate solutions below 𝑦₁(𝑥)=𝜅𝑥²2−𝛼𝑘²𝑥⁵+120𝑥−𝛼𝑘𝑥⁴𝑦24𝜇,₂(𝑥)=𝜅𝑥²2−𝛼𝑘²𝑥⁵+12011𝛼²𝑘³𝑥⁸−𝛼40320³𝑘⁴𝑥¹¹+712800𝑥−𝛼𝑘𝑥⁴+2411𝛼²𝑘²𝑥⁷−𝛼5040³𝑘³𝑥¹⁰𝜇+𝛼64800₂𝑘𝑥⁶−𝛼240³𝑘²𝑥⁹𝜇24192²,𝑦₃(𝑥)=𝜅𝑥²2−𝛼𝑘²x⁵+12011𝛼²𝑘³𝑥⁸−403205𝛼³𝑘⁴𝑥¹¹[phone]𝛼⁴𝑘⁵𝑥¹⁴−871782912005449𝛼⁵𝑘⁶𝑥¹⁷[phone]𝛼⁶𝑘⁷𝑥²⁰−𝛼8935557120000⁷𝑘⁸𝑥²³[phone]𝑥−𝛼𝑘𝑥⁴+2411𝛼²𝑘²𝑥⁷−50405𝛼³𝑘³𝑥¹⁰[phone]𝛼⁴𝑘⁴𝑥¹³−62270208005449𝛼⁵𝑘⁵𝑥¹⁶[phone]𝛼⁶𝑘⁶𝑥¹⁹−𝛼446777856000⁷𝑘⁷𝑥²²𝜇+𝛼2133952128000₂𝑘𝑥⁶−24043𝛼³𝑘²𝑥⁹[phone]𝛼⁴𝑘³𝑥¹²−1596672001147𝛼⁵𝑘⁴𝑥¹⁵[phone]𝛼⁶𝑘⁵𝑥¹⁸−10260214272000197𝛼⁷𝑘⁶𝑥²¹𝜇49145564160000²−𝛼³𝑘𝑥⁸−268823𝛼⁴𝑘²𝑥¹¹[phone]𝛼⁵𝑘³𝑥¹⁴−26417664001829𝛼⁶𝑘⁴𝑥¹⁷+𝛼414554112000⁷𝑘⁵𝑥²⁰𝜇66189312000³−𝛼⁵𝑘²𝑥¹³−329472017𝛼⁶𝑘³𝑥¹⁶+𝛼3251404800⁷𝑘⁴𝑥¹⁹𝜇47259168768⁴. (3.1)

For the case of 𝛼=1/2, we have𝑦₁(𝑥)=𝜅𝑥²2−𝑘²𝑥⁵+240𝑥−𝑘𝑥⁴𝑦48𝜇,₂(𝑥)=𝜅𝑥²2−𝑘²𝑥⁵+24011𝑘³𝑥⁸−𝑘161280⁴𝑥¹¹+5702400𝑥−𝑘𝑥⁴+4811𝑘²𝑥⁷−𝑘20160³𝑥¹⁰𝜇+518400𝑘𝑥⁶−𝑘960²𝑥⁹𝜇193536²,𝑦₃(𝑥)=𝜅𝑥²2−𝑘²𝑥⁵+24011𝑘³𝑥⁸−1612805𝑘⁴𝑥¹¹[phone]𝑘⁵𝑥¹⁴−13948526592005449𝑘⁶𝑥¹⁷+12507689779200083𝑘⁷𝑥²⁰−𝑘571875655680000⁸𝑥²³+6282355064832000𝑥−𝑘𝑥⁴+4811𝑘²𝑥⁷−201605𝑘³𝑥¹⁰[phone]𝑘⁴𝑥¹³−996323328005449𝑘⁵𝑥¹⁶[phone]𝑘⁶𝑥¹⁹−𝑘28593782784000⁷𝑥²²𝜇+273145872384000𝑘𝑥⁶−96043𝑘²𝑥⁹[phone]𝑘³𝑥¹²−25546752001147𝑘⁴𝑥¹⁵+25360957440014057𝑘⁵𝑥¹⁸−656653713408000197𝑘⁶𝑥²¹𝜇6290632212480000²−𝑘𝑥⁸−2150423𝑘²𝑥¹¹[phone]𝑘³𝑥¹⁴−845365248001829𝑘⁴𝑥¹⁷+𝑘26531463168000⁵𝑥²⁰𝜇8472231936000³−𝑘²𝑥¹³−10543104017𝑘³𝑥¹⁶+𝑘208089907200⁴𝑥¹⁹𝜇6049173602304⁴. (3.2) One can see that the result is consistent with the known series solution [11, 12].

In particular, when 𝛼=1 and 𝜇=0, it follows that 𝑦₁(𝑥)=𝜅𝑥²2−𝑘²𝑥⁵,𝑦120₂(𝑥)=𝜅𝑥²2−𝑘²𝑥⁵+12011𝑘³𝑥⁸−𝑘40320⁴𝑥¹¹,𝑦712800₃(𝑥)=𝜅𝑥²2−𝑘²𝑥⁵+12011𝑘³𝑥⁸−403205𝑘⁴𝑥¹¹[phone]𝑘⁵𝑥¹⁴−871782912005449𝑘⁶𝑥¹⁷[phone]𝑘⁷𝑥²⁰−𝑘8935557120000⁸𝑥²³,𝑦49080898944000₄(𝑥)=𝜅𝑥²2−𝑘²𝑥⁵+12011𝑘³𝑥⁸−403205𝑘⁴𝑥¹¹[phone]𝑘⁵𝑥¹⁴−290594304002173649𝑘⁶𝑥¹⁷+35568742809600013722337𝑘⁷𝑥²⁰−11585247657984000027184438601𝑘⁸𝑥²³+1292600836944248832000012320831753849𝑘⁹𝑥²⁶𝑘403291461126605635584000000+⋯−¹⁶𝑥⁴⁷.463172433275878342410240000000 (3.3) In this case, 𝜅≈0.4695999883. For comparison, we refer to another analytical solution obtained by the Adomian decomposition method as follows: 𝑦^A₁(𝑥)=𝜅𝑥²2−𝑘²𝑥⁵,𝑦120_A2(𝑥)=𝜅𝑥²2−𝑘²𝑥⁵+12011𝑘³𝑥⁸,𝑦40320_A3(𝑥)=𝜅𝑥²2−𝑘²𝑥⁵+12011𝑘³𝑥⁸−403205𝑘⁴𝑥¹¹,𝑦532224_A4(𝑥)=𝜅𝑥²2−𝑘²𝑥⁵+12011𝑘³𝑥⁸−403205𝑘⁴𝑥¹¹[phone]𝑘⁵𝑥¹⁴,𝑦29059430400_A5(𝑥)=𝜅𝑥²2−𝑘²𝑥⁵+12011𝑘³𝑥⁸−403205𝑘⁴𝑥¹¹[phone]𝑘⁵𝑥¹⁴−290594304001272379𝑘⁶𝑥¹⁷.118562476032000 (3.4) This solution is based on 𝑦^A_𝑛(𝑥)=𝑦₀(𝑥)+^𝑛_𝑘=1𝑢_𝑘(𝑥),𝑛=1,2,…, (3.5) with 𝑢_𝑘+1(𝑥)=−𝐿⁻¹𝐴_𝑘(𝑥),𝑘=0,1,2,…, (3.6)𝑢₀(𝑥)=𝑦₀(𝑥), and Adomian polynomial 𝐴_𝑘(𝑥) generated by the formula [6] 𝐴_𝑘(𝑥)=1𝑑𝑘!𝑑𝜆^𝑘𝑘𝑘_𝑗=0𝜆^𝑗𝑢_𝑗(𝑥)^𝑘_𝑗=0𝜆^𝑗𝑢_𝑗^′′(𝑥)_𝜆=0, (3.7) where 𝐿⁻¹ is an inverse operator of 𝐿=𝑑³/𝑑𝑥³. Comparing the formulas in (3.3) and (3.4), one can see that the presented analytical solution 𝑦_𝑛(𝑥) has more terms than 𝑦^A_𝑛(𝑥) in each 𝑛th iteration. In other words, 𝑦_𝑛(𝑥)=𝑦^A_𝑛(𝑥)+^{2𝑛−𝑛−1}_𝑗=1𝑑_𝑗𝑥^{3(𝑛+𝑗)+2}, (3.8) for any integer 𝑛≥2. In practice, Figure 1 depicts that the presented solutions 𝑦_𝑛(𝑥) and their derivatives 𝑦^′_𝑛(𝑥), 𝑛=1,3,5 approximate exact ones better than 𝑦^A_𝑛(𝑥) and (𝑦^A_𝑛)^′(𝑥). Therein, we chose the initial solution 𝑦₀(𝑥)=(𝜅/2)𝑥² as given in (2.10) and took a numerical solution for the exact solution which is denoted by 𝑦^∗(𝑥). Moreover, Table 1 includes numerical results of the errors and the CUP times spent in computations for the presented solution 𝑦_𝑛(𝑥) compared with those of 𝑦^A_𝑛(𝑥). The 𝐿_∞ 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 𝐿₂ error means ‖𝑦_𝑛−𝑦^∗‖₂ 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.

(a) Presented solutions

(b) Adomian’s solutions

Figure 1 

Graphs of the presented solutions 𝑦𝑛(𝑥) and their derivatives 𝑦′𝑛(𝑥) in (a) and those of Adomian’s solutions 𝑦A𝑛(𝑥) and their derivatives (𝑦A𝑛)′(𝑥) in (b).

By 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: 𝑦(𝑥)=^∞_𝑘=0(−1)^𝑘𝑝_𝑘𝜅^𝑘+1()!𝑥3𝑘+2^3𝑘+2, (3.9) where 𝑝₀=1 and 𝑝_𝑘=^𝑘−1_𝑗=0𝑝3𝑘−13𝑗^𝑗𝑝_{𝑘−𝑗−1},𝑘≥1. (3.10) Theoretical convergence analysis with extended application of the presented method to more general problems is left for further works.

Acknowledgments

The author would like to show his gratitude to the Department of Mathematics at the University of British Columbia, where the author worked as a visiting scholar for a year. Particularly, the author is heartily thankful to Professor Anthony Peirce for his favor and kind assistance. In addition, the author shows his sincere gratitude to reviewers for their helpful comments and valuable suggestions on the first draft of this paper. This research was supported by Basic Science Research Program through the National Research Foundation of Korea (NRF) funded by the Ministry of Education, Science, and Technology (20110006106).