- Email: [email protected]
An optimal homotopy asymptotic method applied to the steady flow of a fourth-grade fluid past a porous plate
Applied Mathematics Letters 22 (2009) 245–251
Contents lists available at ScienceDirect
Applied Mathematics Letters journal homepage: www.elsevier.com/locate/aml
An optimal homotopy asymptotic method applied to the steady flow of a fourth-grade fluid past a porous plate Vasile Marinca ∗ , Nicolae Herişanu, Constantin Bota, Bogdan Marinca Politehnica University of Timişoara, Bd. Mihai Viteazu 1, 300222 Timişoara, Romania
article
info
Article history: Received 15 August 2007 Received in revised form 20 January 2008 Accepted 24 March 2008 Keywords: Optimal Homotopy Asymptotic Method (OHAM) Fourth-grade fluid Porous plate
a b s t r a c t A new analytic approximate technique for addressing nonlinear problems, namely the Optimal Homotopy Asymptotic Method (OHAM), is proposed and used in an application to the steady flow of a fourth-grade fluid. This approach does not depend upon any small/large parameters. This method provides us with a convenient way to control the convergence of approximation series and adjust convergence regions when necessary. The series solution is developed and the recurrence relations are given explicitly. The results reveal that the proposed method is effective and easy to use. © 2008 Elsevier Ltd. All rights reserved.
1. Introduction Considerable attention has been paid in recent years to problems of flow of non-Newtonian fluids. The governing equations are very complicated and highly nonlinear as compared to those for Newtonian fluids. There are few analytic solutions of the equations involving Newtonian fluids and such solutions become even rarer when equations for nonNewtonian fluids are considered. There is no single model available which clearly exhibits the properties of all nonNewtonian fluids. The aim of the present work is to solve the problem of the steady flow of a fourth-grade fluid past a porous plate using an analytical technique, OHAM. We consider the fourth-grade fluid in the present study for generality. The governing differential equation is highly nonlinear and it is not easy to obtain an analytic solution using traditional methods. Also, the order of the differential equation in the case of fourth-grade fluid is higher than that of the Navier–Stokes equation. The no-slip boundary condition is sufficient for a Newtonian fluid but for a fourth-grade one it may not be sufficient and, therefore, one needs additional conditions at the boundary. The aim here is just to venture further into the regime of nonlinear fluids. In view of this we have established an analytic solution based on the OHAM. 2. Mathematical formulation Consider the steady flow of a fourth-grade fluid past a porous plate. If T is the Cauchy stress tensor, V is the velocity field, ρ is the fluid density, we choose the x-axis parallel to the plate and the y-axis normal to it, the velocity field depending only on y. The flow past a porous plate with suction or injection and with a uniform stream at infinity has a twodimensional structure. The complete set of governing equations of the fourth-grade fluid consists of the incompressibility conditions [1–3] div V = 0
∗ Corresponding author. E-mail address: [email protected] (V. Marinca). 0893-9659/$ – see front matter © 2008 Elsevier Ltd. All rights reserved. doi:10.1016/j.aml.2008.03.019
(1)
V. Marinca et al. / Applied Mathematics Letters 22 (2009) 245–251
246
and the momentum equation dV dt
=
1
ρ
div T.
(2)
After other simple manipulations, from Eqs. (1) and (2) we obtain 2 2 5 d4 u du d3 u du d u d2 u 3d u − ρV0 − γ V = µ 2 − α1 V0 3 + β1 V02 4 + 6(β2 + β3 ) 1 0 dy dy dy dy dy dy2 dy5 " 2 # 2 du d u d − 2(3γ2 + γ3 + γ4 + γ5 + 3γ7 + γ8 )V0 dy dy dy2
(3)
where αi (i = 1, 2), βi (i = 1, 2, 3), γi (i = 1, 2, . . . , 8) are material constants and µ is the coefficient of viscosity, V0 > 0 is the suction velocity and V0 < 0 corresponds to the injection velocity. The boundary conditions on u are u(0) = 0;
u(y) → u0
as y → ∞
(4)
where u0 is the mainstream velocity. Introducing the quantities u¯ =
β¯ =
u u0
y¯ =
;
u0 y v
6(β2 + β3 )u40 V02
ρ
v3
v0
α¯1 =
α1 u20 V02 ; ρv2
β¯ 1 =
;
v¯0 =
;
γ¯ = 2(2γ2 + γ3 + γ4 + γ5 + 3γ7 + γ8 )
u0
;
β1 u40 V04 ; ρv3
γ¯ 1 =
γ1 u60 V06 ; ρv4
(5)
u60 V04
ρv 4
and a new variable
τ = v0 y
(6)
the boundary value problem takes the form d2 u dτ 2
+
du dτ
u(0) = 0;
− α1
d3 u dτ 3
+ β1
d4 u dτ 4
− γ1
d5 u dτ 5
+β
du dτ
2
d2 u dτ 2
− γ 2
du
d2 u
dτ
dτ 2
!2
+
du dτ
2
d3 u dτ 3
=0
u(∞) = 1.
(7) (8)
For simplicity we omitted the bars of the non-dimensional quantities. The fifth-order differential equation (7) subject to the boundary conditions (8) can be treated using OHAM. 3. Basic idea of the Optimal Homotopy Asymptotic Method (OHAM) We apply OHAM to the following differential equation: L(u(τ)) + g(τ) + N(u(τ)) = 0,
B (u) = 0
(9)
where L is a linear operator, τ denotes an independent variable, u(τ) is an unknown function, g(τ) is a known function, N(u(τ)) is a nonlinear operator and B is a boundary operator. By means of OHAM one first constructs a family of equations:
(1 − p)[L(φ(τ, p)) + g(τ)] = H(p)[L(φ(τ, p)) + g(τ) + N(φ(τ, p))],
B(φ(τ, p)) = 0
(10)
where p ∈ [0, 1] is an embedding parameter, H(p) is a nonzero auxiliary function for p 6= 0 and H(0) = 0, φ(τ, p) is an unknown function. Obviously, when p = 0 and p = 1 it holds that
φ(τ, 0) = u0 (τ),
φ(τ, 1) = u(τ).
(11)
Thus, as p increases from 0 to 1, the solution φ(τ, p) varies from u0 (τ) to the solution u(τ), where u0 (τ) is obtained from Eq. (10) for p = 0: L(u0 (τ)) + g(τ) = 0,
B (u0 ) = 0.
(12)
We choose the auxiliary function H(p) in the form H(p) = pC1 + p2 C2 + · · ·
where C1 , C2 , . . . are constants which can be determined later. Expanding φ(τ, p) in a series with respect to p, one has X φ(τ, p, Ci ) = u0 (τ) + uk (τ, Ci )pk , i = 1, 2, . . . . k≥1
(13)
(14)
V. Marinca et al. / Applied Mathematics Letters 22 (2009) 245–251
247
Substituting Eq. (14) into Eq. (10), collecting the same powers of p, and equating each coefficient of p to zero, we obtain L(u1 (τ)) = C1 N0 (u0 (τ)), B (u1 ) = 0
(15)
L(uk (τ) − uk−1 (τ)) = Ck N0 (u0 (τ))
+
k−1 X
Ci [L(uk−i (τ)) + Nk−i (u0 (τ), u1 (τ), . . . , uk−i (τ))],
k = 2, 3, . . . , B (uk ) = 0, Ck = K (t)
(16)
i =1
where Ni , i ≥ 0, are the coefficients of pi in the nonlinear operator N: N(u(τ)) = N0 (u0 (τ)) + pN1 u0 (τ), u1 (τ) + p2 N2 (u0 (τ), u1 (τ), u2 (τ)) + · · ·
(17)
and K (t) is a function which will be defined later. It should be emphasized that the uk for k ≥ 0 are governed by the linear Eq. (12), (15) and (16) with the linear boundary conditions that come from original problem, which can be easily solved. The convergence of the series (14) depends upon the auxiliary constants C1 , C2 , . . .. If it is convergent at p = 1, one has X u(τ, Ci ) = u0 (τ) + uk (τ, Ci ). (18) k≥1
Generally speaking, the solution of Eq. (9) can be determined approximately in the form u˜ (m) = u0 (τ) +
m X
uk (τ, Ci ).
(19)
k=1
We note that the last coefficient Cm can be function of τ . Substituting Eq. (19) into Eq. (9), there results the following residual: R(τ, Ci ) = L(u˜ (m) (τ, Ci )) + g(τ) + N(u˜ (m) (τ, Ci )).
(20)
If R(τ, Ci ) = 0 then u˜ (m) (τ, Ci ) happens to be the exact solution. Generally such a case will not arise for nonlinear problems, but we can minimize the functional Z b J(C1 , C2 , . . . , Cn ) = R2 (τ, C1 , C2 , . . . Cm )dτ (21) a
where a and b are two values, depending on the given problem. The unknown constants Ci (i = 1, 2, . . . m) can be identified from the conditions
∂J ∂J = = · · · = 0. ∂C1 ∂C2
(22)
With these constants known, the approximate solution (of order m) (19) is well determined. The constants C1 , C2 , . . . can be determined in other forms (collocation, Ritz method etc.). It is easy to observe that the so-called Homotopy Perturbation Method (HPM) proposed by He [4,5] is only a special case of Eq. (10) when H(p) = −p, and the Homotopy Analysis Method (HAM) proposed by Liao [6,7] is also another special case of Eq. (10) when H(p) = ph¯ and F (r¯, τ) = 1, where the parameter h¯ is determined from so-called “h¯ -curves”. It can be observed that the method proposed in this work generalizes these two methods using the special (more general) auxiliary function H(p). Another important feature of the proposed method (OHAM) is that using Eq. (30), a minimization of errors is obtained. 4. Application of OHAM to the steady flow of fourth-grade fluid According to Eqs. (7) and (8), we choose the linear operator L(φ(τ, p)) =
∂2 φ (τ, p) ∂φ(τ, p) + ∂τ 2 ∂τ
(23)
and we define a nonlinear operator N(φ(τ, p)) = −α1
∂3 φ(τ, p) ∂4 φ(τ, p) ∂5 φ(τ, p) ∂φ(τ, p) + β1 − γ1 +β 3 4 ∂τ ∂τ ∂τ 5 ∂τ
∂φ(τ, p) ∂2 φ(τ, p) − γ 2 ∂τ ∂τ2
!2
∂φ(τ, p) + ∂τ
2
2
∂2 φ(τ, p) ∂τ 2
!
∂3 φ(τ, p) = 0. ∂τ3
(24)
The initial conditions are
φ(0, p) = 0,
φ(∞, p) = 1.
(25)
V. Marinca et al. / Applied Mathematics Letters 22 (2009) 245–251
248
Zeroth-order problem given by Eq. (12) (g(τ) = 0): d2 u0 (τ) dτ 2
du0 (τ)
+
dτ
= 0,
u0 (0) = 0,
u0 (∞) = 1.
(26)
It is obtained that u0 (τ) = 1 − e−τ .
(27)
First-order problem given by Eq. (15):
∂2 u1 (τ) ∂u1 (τ) ∂3 u0 (τ) ∂4 u0 (τ) ∂5 u0 (τ) ∂u0 (τ) = C1 − α1 + + β1 − γ1 +β 2 3 4 5 ∂τ ∂τ ∂τ ∂τ ∂τ ∂τ
∂u0 (τ) ∂2 u0 (τ) − γ 2 ∂τ ∂τ2 u1 (0) = 0,
!2
∂u0 (τ) + ∂τ
2
2
∂2 u0 (τ) ∂τ 2
!
∂3 u0 (τ) ∂τ3
(28)
u1 (∞) = 0.
The solution of Eq. (28) is given by u1 (τ) =
C1
6
(β + 3γ) e−τ − e−3τ + C1 (α1 + β1 + γ1 )τ e−τ .
(29)
Second-order problem given by Eq. (16) for k = 2:
∂2 u1 (τ) ∂u1 (τ) ∂3 u1 (τ) ∂2 u2 (τ) ∂u2 (τ) ∂2 u1 (τ) ∂u1 (τ) ∂4 u1 (τ) = + C − α + + + + β 1 1 1 ∂τ2 ∂τ ∂τ 2 ∂τ ∂τ 2 ∂τ ∂τ3 ∂τ 4 (
∂5 u1 (τ) ∂u0 (τ) ∂2 u0 (τ) ∂u1 (τ) ∂u0 (τ) + − γ1 +β 2 5 ∂τ ∂τ ∂τ 2 ∂τ ∂τ "
2
∂u0 (τ) ∂2 u0 (τ) ∂2 u1 (τ) ∂u0 (τ) ∂3 u0 (τ) ∂u1 (τ) ∂u0 (τ) +4 + 2 + ∂τ ∂τ 2 ∂τ 2 ∂τ ∂τ 3 ∂τ ∂τ
+ K (τ)
(
#
2
∂3 u0 (τ) ∂τ 3
∂2 u1 (τ) ∂2 u0 (τ) − γ 2 2 ∂τ ∂τ 2
!2
∂u1 (τ) ∂τ
#)
∂2 u0 (τ) ∂3 u0 (τ) ∂u0 (τ) ∂4 u0 (τ) − α + + β 1 1 ∂τ 2 ∂τ ∂τ 3 ∂τ4
∂5 u0 (τ) ∂u0 (τ) − γ1 +β ∂τ5 ∂τ u2 (0) = 0, u2 (∞) = 0.
2
∂2 u0 (τ) ∂τ2
∂2 u0 (τ) ∂u0 (τ) −γ 2 ∂τ2 ∂τ
!2
∂u0 (τ) + ∂τ
2
∂3 u0 (τ) ∂τ 3
#)
(30)
Choosing K (τ) = C2 + C3 e−2τ (C2 and C3 are constants), Eq. (30) has the solution u2 (τ) =
(
"
(α1 + β1 + γ1 ) C12 + C1 + C2 + C12 (2α1 + 3β1 + 4γ1 ) − )
−
C12
+
α1 + β1 + γ1 h
−
β + 3γ
2
(α1 + β1 + λ1 )2 τ2 e−τ +
6 2
(
β + 3γ 6 i
C12 (4β + 15γ) + C3 −
(α1 + β1 + γ1 )τe−3τ +
β + 3γ 40
C12 + C1+ C2 + C12
5(β + 3γ)
C12
9 2
C12
6
α1 +
#
(β + 3γ) τ 27 2 )
(α1 + β1 + γ1 ) 12
β1 +
81 2
γ1 −
e−τ − e−3τ
(5βC12 + 5γ C12 + 2C3 )(e−τ − e−5τ ).
β + 3γ
2
(31)
Substitution of Eqs. (27), (29) and (31) into Eq. (19) yields the second-order approximate solution (m = 2) for Eqs. (7) and (8): u˜ (2) (τ) = 1 + (Aτ 2 + Bτ + C )e−τ + (Dτ + E)e−3τ + F e−5τ
(32)
V. Marinca et al. / Applied Mathematics Letters 22 (2009) 245–251
249
Fig. 1. Comparison between the results obtained from OHAM and Sajid’s results [1] obtained from HAM in case (a): α1 = 0.5, β1 = 1, γ1 = 0.3, β = 1, γ = 1.
Fig. 2. Comparison between the results obtained from OHAM and Sajid’s results [1] obtained from HAM in case (b): α1 = −0.5, β1 = 1, γ1 = 0.3,
β = 0.5, γ = 0.5.
where A=− B=
C12
2
(α1 + β1 + γ1 )
C =
β + 3γ
6
− D=
(α1 + β1 + γ1 )2
C12
2
E = −
− F=−
5C12 12
" C12
+ 2C1 + C2 +
C12 + 2C1 + C2 + C12
9 2
C12
(2α1 + 3β1 + 4γ1 ) −
α1 +
(β + 3γ)(α1 + β1 + γ1 ) +
27 2
β + 3γ 40
β1 +
81 2
γ1 −
C12
6
(β + 3γ)
β + 3γ
2
+
#
α1 + β1 + γ1 h
(5βC12 + 5γ C12 + 2C3 ) − 1
(β + 3γ)(α1 + β1 + γ1 ) β + 3γ 6
C12 + 2C1 + C2 + C12
α1 + β1 + γ1 h 6
(β + 3γ) 40
9 2
α1 + i
27
C12 (4β + 15γ) + C3 +
(5βC12 + 5γ C12 + 2C3 ).
β1 +
81
2 2 5(β + 3γ) 12
γ1 −
β + 3γ
2
C12 (α1 + β1 + γ1 )
6
C12 (4β + 15γ) + C3
i
(33)
V. Marinca et al. / Applied Mathematics Letters 22 (2009) 245–251
250
Fig. 3. Residual R(τ) given by (34) for case (a).
Fig. 4. Residual R(τ) given by (34) for case (b).
Substituting the approximate solution of the second-order (32) into Eq. (20) yields the residual and the functional J, respectively: !2 2 2 2 d3 u˜ d u˜ d3 u˜ d2 u˜ du˜ d4 u˜ d5 u˜ du˜ du˜ d2 u˜ du˜ − α1 3 + β1 4 − γ1 5 + β R(τ, C1 , C2 , C3 ) = + − −γ 2 + (34) dτ 2 dτ dτ dτ dτ dτ dτ 2 dτ dτ 2 dτ dτ 3 Z ∞ J(C1 , C2 , C3 ) = R2 (τ, C1 , C2 , C3 )dτ. (35) 0
The constants C1 , C2 , C3 result from the conditions (22). In the particular cases: (a) α1 = −0.5; β1 = 1; γ1 = 0.3; β = 1; γ = 1, it is obtained that C1 = −0.12454086; C2 = −0.14906175; C3 = −0.08325414, and the approximate solution of second order in the form u(τ) = 1 + (−0.00496333673 · τ 2 − 0.274671983803 · τ − 1.027435052878)e−τ
+ (0.02481668364 · τ + 0.02629465182)e−3τ + 0.001140401055e−5τ .
(36)
(b) α1 = −0.5; β1 = 1; γ1 = 0.3; β = 0.5; γ = 0.5, it is obtained that C1 = −0.08771268; C2 = −0.06229223; C3 = −0.01754022 and the approximate solution of second order u(τ) = 1 + (−0.00246192499 · τ 2 − 0.166375477221 · τ − 1.016783012456)e−τ
+ (0.006154812473 · τ + 0.01695236848)e−3τ − 0.000169356021e−5τ .
(37)
Figs. 1 and 2 show a comparison between the results obtained using OHAM and the results obtained by Sajid et al. [1] who employed HAM for the same problem. One can see from Figs. 3 and 4 the accuracy of the solution obtained by the present method, which in both cases is very good. The residual R(τ) has maximum magnitude 0.06 in case (a) and 0.03 in case (b), which proves the accuracy of the approximate solution. 5. Conclusions In this paper, a new technique is proposed for solving a nonlinear equation. This procedure is explicit, effective and has a distinct advantage over usual approximation methods in that the approximate solution obtained here is valid not only for weakly nonlinear equations, but also for strongly nonlinear ones. The convergence and low error are remarkable. In this method we control the convergence using a number k of auxiliary constants C1 , C2 , . . . , Ck which are optimally determined. In this respect the functional J given by (21) must be minimized, which means that conditions (22) must be satisfied. This
V. Marinca et al. / Applied Mathematics Letters 22 (2009) 245–251
251
condition allows the determination of the parameters C1 , C2 . . . , Ck and with these known, the approximate solution is well determined. The error obtained for the solution is remarkably low and this is proved by Figs. 3 and 4, where the residual R(τ) given by (34) is plotted for cases (a) and (b). References [1] M. Sajid, T. Hayat, S. Asghar, On the analytic solution of the steady flow of a fourth grade fluid, Phys. Lett. A 355 (2006) 18–26. [2] T. Hayat, M. Sajid, On analytic solution for thin film flow of fourth grade fluid down a vertical cylinder, Phys. Lett. A 361 (2007) 316–322. [3] K.R. Rajagopal, On the boundary conditions for fluids of the differential type, in: A. Sequira (Ed.), J. Navier–Stokes Equation and Related Non-linear Problems, Plenum New York, 1995, pp. 273–278. [4] J.H. He, Homotopy perturbation technique, Comput. Methods Appl. Mech. Eng. 178 (1999) 257–262. [5] J.H. He, Homotopy perturbation method for solving boundary value problems, Phys. Lett. A 350 (1–2) (2006) 1187–1193. [6] S.J. Liao, A second-order approximate analytical solution of a simple pendulum by the process analysis method, ASME, J. Appl. Mech. 59 (1992) 970–975. [7] S.J. Liao, A.T. Chwang, Application of homotopy analysis method in nonlinear oscillations, ASME J. Appl. Mech. 65 (1998) 914–922.
Contents lists available at ScienceDirect
Applied Mathematics Letters journal homepage: www.elsevier.com/locate/aml
An optimal homotopy asymptotic method applied to the steady flow of a fourth-grade fluid past a porous plate Vasile Marinca ∗ , Nicolae Herişanu, Constantin Bota, Bogdan Marinca Politehnica University of Timişoara, Bd. Mihai Viteazu 1, 300222 Timişoara, Romania
article
info
Article history: Received 15 August 2007 Received in revised form 20 January 2008 Accepted 24 March 2008 Keywords: Optimal Homotopy Asymptotic Method (OHAM) Fourth-grade fluid Porous plate
a b s t r a c t A new analytic approximate technique for addressing nonlinear problems, namely the Optimal Homotopy Asymptotic Method (OHAM), is proposed and used in an application to the steady flow of a fourth-grade fluid. This approach does not depend upon any small/large parameters. This method provides us with a convenient way to control the convergence of approximation series and adjust convergence regions when necessary. The series solution is developed and the recurrence relations are given explicitly. The results reveal that the proposed method is effective and easy to use. © 2008 Elsevier Ltd. All rights reserved.
1. Introduction Considerable attention has been paid in recent years to problems of flow of non-Newtonian fluids. The governing equations are very complicated and highly nonlinear as compared to those for Newtonian fluids. There are few analytic solutions of the equations involving Newtonian fluids and such solutions become even rarer when equations for nonNewtonian fluids are considered. There is no single model available which clearly exhibits the properties of all nonNewtonian fluids. The aim of the present work is to solve the problem of the steady flow of a fourth-grade fluid past a porous plate using an analytical technique, OHAM. We consider the fourth-grade fluid in the present study for generality. The governing differential equation is highly nonlinear and it is not easy to obtain an analytic solution using traditional methods. Also, the order of the differential equation in the case of fourth-grade fluid is higher than that of the Navier–Stokes equation. The no-slip boundary condition is sufficient for a Newtonian fluid but for a fourth-grade one it may not be sufficient and, therefore, one needs additional conditions at the boundary. The aim here is just to venture further into the regime of nonlinear fluids. In view of this we have established an analytic solution based on the OHAM. 2. Mathematical formulation Consider the steady flow of a fourth-grade fluid past a porous plate. If T is the Cauchy stress tensor, V is the velocity field, ρ is the fluid density, we choose the x-axis parallel to the plate and the y-axis normal to it, the velocity field depending only on y. The flow past a porous plate with suction or injection and with a uniform stream at infinity has a twodimensional structure. The complete set of governing equations of the fourth-grade fluid consists of the incompressibility conditions [1–3] div V = 0
∗ Corresponding author. E-mail address: [email protected] (V. Marinca). 0893-9659/$ – see front matter © 2008 Elsevier Ltd. All rights reserved. doi:10.1016/j.aml.2008.03.019
(1)
V. Marinca et al. / Applied Mathematics Letters 22 (2009) 245–251
246
and the momentum equation dV dt
=
1
ρ
div T.
(2)
After other simple manipulations, from Eqs. (1) and (2) we obtain 2 2 5 d4 u du d3 u du d u d2 u 3d u − ρV0 − γ V = µ 2 − α1 V0 3 + β1 V02 4 + 6(β2 + β3 ) 1 0 dy dy dy dy dy dy2 dy5 " 2 # 2 du d u d − 2(3γ2 + γ3 + γ4 + γ5 + 3γ7 + γ8 )V0 dy dy dy2
(3)
where αi (i = 1, 2), βi (i = 1, 2, 3), γi (i = 1, 2, . . . , 8) are material constants and µ is the coefficient of viscosity, V0 > 0 is the suction velocity and V0 < 0 corresponds to the injection velocity. The boundary conditions on u are u(0) = 0;
u(y) → u0
as y → ∞
(4)
where u0 is the mainstream velocity. Introducing the quantities u¯ =
β¯ =
u u0
y¯ =
;
u0 y v
6(β2 + β3 )u40 V02
ρ
v3
v0
α¯1 =
α1 u20 V02 ; ρv2
β¯ 1 =
;
v¯0 =
;
γ¯ = 2(2γ2 + γ3 + γ4 + γ5 + 3γ7 + γ8 )
u0
;
β1 u40 V04 ; ρv3
γ¯ 1 =
γ1 u60 V06 ; ρv4
(5)
u60 V04
ρv 4
and a new variable
τ = v0 y
(6)
the boundary value problem takes the form d2 u dτ 2
+
du dτ
u(0) = 0;
− α1
d3 u dτ 3
+ β1
d4 u dτ 4
− γ1
d5 u dτ 5
+β
du dτ
2
d2 u dτ 2
− γ 2
du
d2 u
dτ
dτ 2
!2
+
du dτ
2
d3 u dτ 3
=0
u(∞) = 1.
(7) (8)
For simplicity we omitted the bars of the non-dimensional quantities. The fifth-order differential equation (7) subject to the boundary conditions (8) can be treated using OHAM. 3. Basic idea of the Optimal Homotopy Asymptotic Method (OHAM) We apply OHAM to the following differential equation: L(u(τ)) + g(τ) + N(u(τ)) = 0,
B (u) = 0
(9)
where L is a linear operator, τ denotes an independent variable, u(τ) is an unknown function, g(τ) is a known function, N(u(τ)) is a nonlinear operator and B is a boundary operator. By means of OHAM one first constructs a family of equations:
(1 − p)[L(φ(τ, p)) + g(τ)] = H(p)[L(φ(τ, p)) + g(τ) + N(φ(τ, p))],
B(φ(τ, p)) = 0
(10)
where p ∈ [0, 1] is an embedding parameter, H(p) is a nonzero auxiliary function for p 6= 0 and H(0) = 0, φ(τ, p) is an unknown function. Obviously, when p = 0 and p = 1 it holds that
φ(τ, 0) = u0 (τ),
φ(τ, 1) = u(τ).
(11)
Thus, as p increases from 0 to 1, the solution φ(τ, p) varies from u0 (τ) to the solution u(τ), where u0 (τ) is obtained from Eq. (10) for p = 0: L(u0 (τ)) + g(τ) = 0,
B (u0 ) = 0.
(12)
We choose the auxiliary function H(p) in the form H(p) = pC1 + p2 C2 + · · ·
where C1 , C2 , . . . are constants which can be determined later. Expanding φ(τ, p) in a series with respect to p, one has X φ(τ, p, Ci ) = u0 (τ) + uk (τ, Ci )pk , i = 1, 2, . . . . k≥1
(13)
(14)
V. Marinca et al. / Applied Mathematics Letters 22 (2009) 245–251
247
Substituting Eq. (14) into Eq. (10), collecting the same powers of p, and equating each coefficient of p to zero, we obtain L(u1 (τ)) = C1 N0 (u0 (τ)), B (u1 ) = 0
(15)
L(uk (τ) − uk−1 (τ)) = Ck N0 (u0 (τ))
+
k−1 X
Ci [L(uk−i (τ)) + Nk−i (u0 (τ), u1 (τ), . . . , uk−i (τ))],
k = 2, 3, . . . , B (uk ) = 0, Ck = K (t)
(16)
i =1
where Ni , i ≥ 0, are the coefficients of pi in the nonlinear operator N: N(u(τ)) = N0 (u0 (τ)) + pN1 u0 (τ), u1 (τ) + p2 N2 (u0 (τ), u1 (τ), u2 (τ)) + · · ·
(17)
and K (t) is a function which will be defined later. It should be emphasized that the uk for k ≥ 0 are governed by the linear Eq. (12), (15) and (16) with the linear boundary conditions that come from original problem, which can be easily solved. The convergence of the series (14) depends upon the auxiliary constants C1 , C2 , . . .. If it is convergent at p = 1, one has X u(τ, Ci ) = u0 (τ) + uk (τ, Ci ). (18) k≥1
Generally speaking, the solution of Eq. (9) can be determined approximately in the form u˜ (m) = u0 (τ) +
m X
uk (τ, Ci ).
(19)
k=1
We note that the last coefficient Cm can be function of τ . Substituting Eq. (19) into Eq. (9), there results the following residual: R(τ, Ci ) = L(u˜ (m) (τ, Ci )) + g(τ) + N(u˜ (m) (τ, Ci )).
(20)
If R(τ, Ci ) = 0 then u˜ (m) (τ, Ci ) happens to be the exact solution. Generally such a case will not arise for nonlinear problems, but we can minimize the functional Z b J(C1 , C2 , . . . , Cn ) = R2 (τ, C1 , C2 , . . . Cm )dτ (21) a
where a and b are two values, depending on the given problem. The unknown constants Ci (i = 1, 2, . . . m) can be identified from the conditions
∂J ∂J = = · · · = 0. ∂C1 ∂C2
(22)
With these constants known, the approximate solution (of order m) (19) is well determined. The constants C1 , C2 , . . . can be determined in other forms (collocation, Ritz method etc.). It is easy to observe that the so-called Homotopy Perturbation Method (HPM) proposed by He [4,5] is only a special case of Eq. (10) when H(p) = −p, and the Homotopy Analysis Method (HAM) proposed by Liao [6,7] is also another special case of Eq. (10) when H(p) = ph¯ and F (r¯, τ) = 1, where the parameter h¯ is determined from so-called “h¯ -curves”. It can be observed that the method proposed in this work generalizes these two methods using the special (more general) auxiliary function H(p). Another important feature of the proposed method (OHAM) is that using Eq. (30), a minimization of errors is obtained. 4. Application of OHAM to the steady flow of fourth-grade fluid According to Eqs. (7) and (8), we choose the linear operator L(φ(τ, p)) =
∂2 φ (τ, p) ∂φ(τ, p) + ∂τ 2 ∂τ
(23)
and we define a nonlinear operator N(φ(τ, p)) = −α1
∂3 φ(τ, p) ∂4 φ(τ, p) ∂5 φ(τ, p) ∂φ(τ, p) + β1 − γ1 +β 3 4 ∂τ ∂τ ∂τ 5 ∂τ
∂φ(τ, p) ∂2 φ(τ, p) − γ 2 ∂τ ∂τ2
!2
∂φ(τ, p) + ∂τ
2
2
∂2 φ(τ, p) ∂τ 2
!
∂3 φ(τ, p) = 0. ∂τ3
(24)
The initial conditions are
φ(0, p) = 0,
φ(∞, p) = 1.
(25)
V. Marinca et al. / Applied Mathematics Letters 22 (2009) 245–251
248
Zeroth-order problem given by Eq. (12) (g(τ) = 0): d2 u0 (τ) dτ 2
du0 (τ)
+
dτ
= 0,
u0 (0) = 0,
u0 (∞) = 1.
(26)
It is obtained that u0 (τ) = 1 − e−τ .
(27)
First-order problem given by Eq. (15):
∂2 u1 (τ) ∂u1 (τ) ∂3 u0 (τ) ∂4 u0 (τ) ∂5 u0 (τ) ∂u0 (τ) = C1 − α1 + + β1 − γ1 +β 2 3 4 5 ∂τ ∂τ ∂τ ∂τ ∂τ ∂τ
∂u0 (τ) ∂2 u0 (τ) − γ 2 ∂τ ∂τ2 u1 (0) = 0,
!2
∂u0 (τ) + ∂τ
2
2
∂2 u0 (τ) ∂τ 2
!
∂3 u0 (τ) ∂τ3
(28)
u1 (∞) = 0.
The solution of Eq. (28) is given by u1 (τ) =
C1
6
(β + 3γ) e−τ − e−3τ + C1 (α1 + β1 + γ1 )τ e−τ .
(29)
Second-order problem given by Eq. (16) for k = 2:
∂2 u1 (τ) ∂u1 (τ) ∂3 u1 (τ) ∂2 u2 (τ) ∂u2 (τ) ∂2 u1 (τ) ∂u1 (τ) ∂4 u1 (τ) = + C − α + + + + β 1 1 1 ∂τ2 ∂τ ∂τ 2 ∂τ ∂τ 2 ∂τ ∂τ3 ∂τ 4 (
∂5 u1 (τ) ∂u0 (τ) ∂2 u0 (τ) ∂u1 (τ) ∂u0 (τ) + − γ1 +β 2 5 ∂τ ∂τ ∂τ 2 ∂τ ∂τ "
2
∂u0 (τ) ∂2 u0 (τ) ∂2 u1 (τ) ∂u0 (τ) ∂3 u0 (τ) ∂u1 (τ) ∂u0 (τ) +4 + 2 + ∂τ ∂τ 2 ∂τ 2 ∂τ ∂τ 3 ∂τ ∂τ
+ K (τ)
(
#
2
∂3 u0 (τ) ∂τ 3
∂2 u1 (τ) ∂2 u0 (τ) − γ 2 2 ∂τ ∂τ 2
!2
∂u1 (τ) ∂τ
#)
∂2 u0 (τ) ∂3 u0 (τ) ∂u0 (τ) ∂4 u0 (τ) − α + + β 1 1 ∂τ 2 ∂τ ∂τ 3 ∂τ4
∂5 u0 (τ) ∂u0 (τ) − γ1 +β ∂τ5 ∂τ u2 (0) = 0, u2 (∞) = 0.
2
∂2 u0 (τ) ∂τ2
∂2 u0 (τ) ∂u0 (τ) −γ 2 ∂τ2 ∂τ
!2
∂u0 (τ) + ∂τ
2
∂3 u0 (τ) ∂τ 3
#)
(30)
Choosing K (τ) = C2 + C3 e−2τ (C2 and C3 are constants), Eq. (30) has the solution u2 (τ) =
(
"
(α1 + β1 + γ1 ) C12 + C1 + C2 + C12 (2α1 + 3β1 + 4γ1 ) − )
−
C12
+
α1 + β1 + γ1 h
−
β + 3γ
2
(α1 + β1 + λ1 )2 τ2 e−τ +
6 2
(
β + 3γ 6 i
C12 (4β + 15γ) + C3 −
(α1 + β1 + γ1 )τe−3τ +
β + 3γ 40
C12 + C1+ C2 + C12
5(β + 3γ)
C12
9 2
C12
6
α1 +
#
(β + 3γ) τ 27 2 )
(α1 + β1 + γ1 ) 12
β1 +
81 2
γ1 −
e−τ − e−3τ
(5βC12 + 5γ C12 + 2C3 )(e−τ − e−5τ ).
β + 3γ
2
(31)
Substitution of Eqs. (27), (29) and (31) into Eq. (19) yields the second-order approximate solution (m = 2) for Eqs. (7) and (8): u˜ (2) (τ) = 1 + (Aτ 2 + Bτ + C )e−τ + (Dτ + E)e−3τ + F e−5τ
(32)
V. Marinca et al. / Applied Mathematics Letters 22 (2009) 245–251
249
Fig. 1. Comparison between the results obtained from OHAM and Sajid’s results [1] obtained from HAM in case (a): α1 = 0.5, β1 = 1, γ1 = 0.3, β = 1, γ = 1.
Fig. 2. Comparison between the results obtained from OHAM and Sajid’s results [1] obtained from HAM in case (b): α1 = −0.5, β1 = 1, γ1 = 0.3,
β = 0.5, γ = 0.5.
where A=− B=
C12
2
(α1 + β1 + γ1 )
C =
β + 3γ
6
− D=
(α1 + β1 + γ1 )2
C12
2
E = −
− F=−
5C12 12
" C12
+ 2C1 + C2 +
C12 + 2C1 + C2 + C12
9 2
C12
(2α1 + 3β1 + 4γ1 ) −
α1 +
(β + 3γ)(α1 + β1 + γ1 ) +
27 2
β + 3γ 40
β1 +
81 2
γ1 −
C12
6
(β + 3γ)
β + 3γ
2
+
#
α1 + β1 + γ1 h
(5βC12 + 5γ C12 + 2C3 ) − 1
(β + 3γ)(α1 + β1 + γ1 ) β + 3γ 6
C12 + 2C1 + C2 + C12
α1 + β1 + γ1 h 6
(β + 3γ) 40
9 2
α1 + i
27
C12 (4β + 15γ) + C3 +
(5βC12 + 5γ C12 + 2C3 ).
β1 +
81
2 2 5(β + 3γ) 12
γ1 −
β + 3γ
2
C12 (α1 + β1 + γ1 )
6
C12 (4β + 15γ) + C3
i
(33)
V. Marinca et al. / Applied Mathematics Letters 22 (2009) 245–251
250
Fig. 3. Residual R(τ) given by (34) for case (a).
Fig. 4. Residual R(τ) given by (34) for case (b).
Substituting the approximate solution of the second-order (32) into Eq. (20) yields the residual and the functional J, respectively: !2 2 2 2 d3 u˜ d u˜ d3 u˜ d2 u˜ du˜ d4 u˜ d5 u˜ du˜ du˜ d2 u˜ du˜ − α1 3 + β1 4 − γ1 5 + β R(τ, C1 , C2 , C3 ) = + − −γ 2 + (34) dτ 2 dτ dτ dτ dτ dτ dτ 2 dτ dτ 2 dτ dτ 3 Z ∞ J(C1 , C2 , C3 ) = R2 (τ, C1 , C2 , C3 )dτ. (35) 0
The constants C1 , C2 , C3 result from the conditions (22). In the particular cases: (a) α1 = −0.5; β1 = 1; γ1 = 0.3; β = 1; γ = 1, it is obtained that C1 = −0.12454086; C2 = −0.14906175; C3 = −0.08325414, and the approximate solution of second order in the form u(τ) = 1 + (−0.00496333673 · τ 2 − 0.274671983803 · τ − 1.027435052878)e−τ
+ (0.02481668364 · τ + 0.02629465182)e−3τ + 0.001140401055e−5τ .
(36)
(b) α1 = −0.5; β1 = 1; γ1 = 0.3; β = 0.5; γ = 0.5, it is obtained that C1 = −0.08771268; C2 = −0.06229223; C3 = −0.01754022 and the approximate solution of second order u(τ) = 1 + (−0.00246192499 · τ 2 − 0.166375477221 · τ − 1.016783012456)e−τ
+ (0.006154812473 · τ + 0.01695236848)e−3τ − 0.000169356021e−5τ .
(37)
Figs. 1 and 2 show a comparison between the results obtained using OHAM and the results obtained by Sajid et al. [1] who employed HAM for the same problem. One can see from Figs. 3 and 4 the accuracy of the solution obtained by the present method, which in both cases is very good. The residual R(τ) has maximum magnitude 0.06 in case (a) and 0.03 in case (b), which proves the accuracy of the approximate solution. 5. Conclusions In this paper, a new technique is proposed for solving a nonlinear equation. This procedure is explicit, effective and has a distinct advantage over usual approximation methods in that the approximate solution obtained here is valid not only for weakly nonlinear equations, but also for strongly nonlinear ones. The convergence and low error are remarkable. In this method we control the convergence using a number k of auxiliary constants C1 , C2 , . . . , Ck which are optimally determined. In this respect the functional J given by (21) must be minimized, which means that conditions (22) must be satisfied. This
V. Marinca et al. / Applied Mathematics Letters 22 (2009) 245–251
251
condition allows the determination of the parameters C1 , C2 . . . , Ck and with these known, the approximate solution is well determined. The error obtained for the solution is remarkably low and this is proved by Figs. 3 and 4, where the residual R(τ) given by (34) is plotted for cases (a) and (b). References [1] M. Sajid, T. Hayat, S. Asghar, On the analytic solution of the steady flow of a fourth grade fluid, Phys. Lett. A 355 (2006) 18–26. [2] T. Hayat, M. Sajid, On analytic solution for thin film flow of fourth grade fluid down a vertical cylinder, Phys. Lett. A 361 (2007) 316–322. [3] K.R. Rajagopal, On the boundary conditions for fluids of the differential type, in: A. Sequira (Ed.), J. Navier–Stokes Equation and Related Non-linear Problems, Plenum New York, 1995, pp. 273–278. [4] J.H. He, Homotopy perturbation technique, Comput. Methods Appl. Mech. Eng. 178 (1999) 257–262. [5] J.H. He, Homotopy perturbation method for solving boundary value problems, Phys. Lett. A 350 (1–2) (2006) 1187–1193. [6] S.J. Liao, A second-order approximate analytical solution of a simple pendulum by the process analysis method, ASME, J. Appl. Mech. 59 (1992) 970–975. [7] S.J. Liao, A.T. Chwang, Application of homotopy analysis method in nonlinear oscillations, ASME J. Appl. Mech. 65 (1998) 914–922.