S1793962311000499 A NEW ALGORITHM FOR ... - Semantic Scholar
August 8, S1793962311000499
2011 16:37 WSPC/262-IJMSSC/S1793-9623
International Journal of Modeling, Simulation, and Scientific Computing Vol. 2, No. 3 (2011) 355–373 c World Scientific Publishing Company DOI: 10.1142/S1793962311000499
A NEW ALGORITHM FOR SOLVING NONLINEAR BOUNDARY VALUE PROBLEMS ARISING IN HEAT TRANSFER
S. S. MOTSA Mathematics Department University of Swaziland, Private Bag 4 Kwaluseni, M201, Swaziland [email protected] Received 27 May 2010 Revised 18 April 2011 Accepted 30 April 2011
In this paper, a very efficient and easy-to-use successive linearization approach for solving nonlinear differential equations is proposed. The implementation of the method is demonstrated by solving three nonlinear differential equations of different complexities arising in heat transfer. New closed form explicit analytic solutions of some of the governing nonlinear equations are obtained and compared with the results of the proposed method and with numerical solutions from the MATLAB in-built routine bvp4c. Keywords: Successive collocation.
linearization
method;
convective–radiative
fin;
spectral
1. Introduction Rectangular fins are often used to enhance heat exchange between heated surfaces and the surrounding environment. In recent years, the problem of heat transfer through extended fins has attracted the interest of many researchers in view of its numerous applications on heat exchangers, electronic components, semiconductors, solar collectors, and many other engineering applications. For most practical applications, the mathematical model of the underlying heat transfer problem consists of nonlinear differential equations whose closed form analytical solutions are either very difficult or impossible to obtain. As a result many researchers have recently devoted attention toward developing new analytical or numerical approximation algorithms or improving existing methods for solving these problems. In the last decade, several methods for solving nonlinear equations, including the Adomian decomposition method,1–4 differential transform method,5 homotopy analysis method,6–11 and homotopy perturbation method,12–14 Taylor series15,16 355
August 8, S1793962311000499
356
2011 16:37 WSPC/262-IJMSSC/S1793-9623
S. S. Motsa
and variational iteration method,17–20 have been applied to solve nonlinear heat transfer equations. In this work, a new numerical approach is proposed that is based on successive linearization of the governing equations. The applicability of the method is tested on three different nonlinear heat transfer equations of different complexity. One of the test problems is that of steady conduction in a slab with variable thermal conductivity that was previously discussed in Ref. 14 using the homotopy perturbation method (HPM) and later studied in Ref. 11 using the homotopy analysis method (HAM). It was proved in Ref. 11 that the HPM approach used in Ref. 14 is only valid for small thermal conductivity parameters (0 ≤ < 1). The second problem considered in this work is the nonlinear convective straight fin problem with temperature dependent thermal conductivity that was investigated in Refs. 5, 9, 10, 13, 15, 16 and 18. The new linearization approach was also applied on the nonlinear convective– radiative conduction problem with temperature dependent thermal conductivity that was investigated in Ref. 19 using the variational iteration method. Here, the nonlinear problems are solved to obtain the distribution of temperature using the new linearization approach and compared with results obtained using the Matlab ODE solver bvp4c, which is a finite difference approach that implements the 3-stage Lobatto III formula (see for example Ref. 21). Simple closed form analytic solutions of some of the governing nonlinear equations, which are subsequently used to validate the results of the proposed linearization method are also reported. 2. Successive Linearization Method Solution In this section, the basic idea of the proposed method of successive linearization is described and it is implemented in three test problems arising in nonlinear heat transfer. It is begun by describing the basic idea behind this method of solution. Consider a nonlinear boundary value problem of the form L[y(x), y (x), y (x)] + N [y(x), y (x), y (x)] = 0,
(1)
where y(x) is an unknown function, x is an independent variable and the primes denote ordinary differentiation with respect to x. The functions L and N represent the linear and nonlinear components of the governing equation, respectively. For illustrative purposes, Eq. (1) is to be solved for x ∈ [a, b] subjected to the boundary conditions y(a) = ya ,
y(b) = yb ,
(2)
where ya and yb are given constants. For the initial guess of the solution of Eq. (1), a function that satisfies the boundary conditions (2) is chosen. It is convenient to consider polynomial functions. Thus, a suitable initial approximation, denoted by y0 (x), that satisfies Eq. (2) would be a straight line.
August 8, S1793962311000499
2011 16:37 WSPC/262-IJMSSC/S1793-9623
New Algorithm for Solving Nonlinear BVPs in Heat Transfer
357
A function Y1 (x) is defined to represent the vertical difference between y(x) and the initial guess y0 (x) (see Fig. 1), that is Y1 (x) = y(x) − y0 (x),
or y(x) = y0 (x) + Y1 (x).
(3)
Substituting Eq. (3) in (1) gives L[Y1 , Y1 , Y1 ] + N [y0 + Y1 , y0 + Y1 , y0 + Y1 ] = −L[y0 , y0 , y0 ].
(4)
Since y0 (x) is given, solving Eq. (4) would yield an exact solution for Y1 (x). However, since the equation is nonlinear, it may not be possible to find an exact solution. Therefore, an approximate solution is needed, which is obtained by solving the linear part of the equation assuming that Y1 and its derivatives are small. This assumption enables us to use the Taylor series method to linearize the equation. If Y1 (x) is the solution of the full Eq. (4), let y1 (x) denote the solution of the linearized version of (4). Expanding (4) using Taylor series and neglecting higher order terms gives ∂N ∂N ∂N + y1 + y1 L[y1 , y1 , y1 ] + y1 ∂y1 (y0 ,y ,y ) ∂y1 (y0 ,y ,y ) ∂y1 (y0 ,y ,y ) 0
=
−L[y0 , y0 , y0 ]
−
0
0
0
N [y0 , y0 , y0 ].
0
0
(5)
Equation (5) can be written in compact form as L[y1 , y1 , y1 ] + a0,0 y1 + a1,0 y1 + a2,0 y1 = r0 (x),
(6)
subject to the boundary conditions y1 (a) = 0,
y1 (b) = 0,
(7)
where ∂N a0,0 (x) = ∂y1 (y0 ,y ,y ) 0 0 ∂N a1,0 (x) = ∂y1 (y0 ,y ,y ) 0 0 ∂N a2,0 (x) = ∂y1
(8)
(9)
(10)
(y0 ,y0 ,y0 )
r0 (x) = −L[y0 , y0 , y0 ] − N [y0 , y0 , y0 ].
(11)
Since the right-hand side of Eq. (6) is known and the left-hand side is linear, the equation can be solved for y1 (x). Assuming that the solution of the linear part (6)
August 8, S1793962311000499
358
2011 16:37 WSPC/262-IJMSSC/S1793-9623
S. S. Motsa
is close to the solution of Eq. (4), that is, Y1 (x) ≈ y1 (x), the current estimate (first-order) of the solution y(x) is y(x) ≈ y0 (x) + y1 (x).
(12)
To improve on this solution, a slack function, Y2 (x) is defined that when added to y1 (x) gives Y1 (x) (see Fig. 1), that is Y1 (x) = Y2 (x) + y1 (x).
(13)
Since y1 (x) is now known (as a solution of Eq. 6), Eq. (13) is substituted in Eq. (4) to obtain L[Y2 , Y2 , Y2 ] + N [y0 + y1 + Y2 , y0 + y1 + Y2 , y0 + y1 + Y2 ] = −L[y0 + y1 , y0 + y1 , y0 + y1 ].
(14)
Solving Eq. (14) would result in an exact solution for Y2 (x). But, since the equation is nonlinear, it may not be possible to find an exact solution. Therefore, the equation is linearized using Taylor series expansion and the resulting linear equation is solved. The solution of the linear version of Eq. (14) is denoted by y2 (x), such that Y2 (x) ≈ y2 (x). Setting Y2 (x) = y2 (x) and expanding Eq. (14), for small y2 (x) and its derivatives gives L[y2 , y2 , y2 ] + a0,1 y2 + a1,1 y2 + a2,1 y2 = r1 (x),
(15)
subject to the boundary conditions y2 (a) = 0,
y2 (b) = 0,
(16)
where ∂N a0,1 (x) = ∂y2 (y0 +y1 ,y +y ,y +y ) 0 1 0 1 ∂N a1,1 (x) = ∂y2 (y0 +y1 ,y +y ,y +y ) 0 1 0 1 ∂N a2,1 (x) = ∂y2 (y0 +y1 ,y +y ,y +y ) 0
1
0
(17)
(18)
(19)
1
r1 (x) = −L[y0 + y1 , y0 + y1 , y0 + y1 ] − N [y0 + y1 , y0 + y1 , y0 + y1 ].
(20)
After solving (15), the current (second-order) estimate of the solution y(x) is y(x) ≈ y0 (x) + y1 (x) + y2 (x).
(21)
August 8, S1793962311000499
2011 16:37 WSPC/262-IJMSSC/S1793-9623
New Algorithm for Solving Nonlinear BVPs in Heat Transfer
359
Next, Y3 (x) is defined (see Fig. 1) such that Y2 (x) = Y3 (x) + y2 (x).
(22)
Equation (22) is substituted in the nonlinear equation (14) and the linearization process described above is repeated. This process is repeated for m = 3, 4, 5, . . . , i. In general, one have Yi (x) = Yi+1 (x) + yi (x).
(23)
Thus y(x) is obtained as y(x) = Y1 (x) + y0 (x),
(24)
= Y2 (x) + y1 (x) + y0 (x),
(25)
= Y3 (x) + y2 (x) + y1 (x) + y0 (x), .. .
(26)
= Yi+1 (x) + yi (x) + · · · + y3 (x) + y2 (x) + y1 (x) + y0 (x),
(27)
= Yi+1 (x) +
i
ym (x).
m=0
The procedure for obtaining each Yi (x) is illustrated in Fig. 1 for i = 1, 2, 3.
Fig. 1.
Geometric representation of the functions Yi (i = 1, 2, 3).
(28)
August 8, S1793962311000499
360
2011 16:37 WSPC/262-IJMSSC/S1793-9623
S. S. Motsa
It is noted that when i becomes large, Yi+1 becomes increasingly smaller. Thus, for large i, the ith order solution of y(x) can be approximated by i
y(x) =
ym (x) = yi (x) +
m=0
i−1
ym (x).
(29)
m=0
Starting from a known initial guess y0 (x), the solutions for yi (x) can be obtained by successively linearizing the governing equation (1) and solving the resulting linear equation for yi (x) given that the previous guess yi−1 (x) is known. The general form of the linearized equation that is successively solved for yi (x) is given by L[yi , yi , yi ] + a0,i−1 yi + a1,i−1 yi + a2,i−1 yi = ri−1 (x),
(30)
subject to the boundary conditions yi (a) = 0,
yi (b) = 0,
(31)
where
∂N a0,i−1 (x) = ∂yi (Pi−1
m=0
ym ,
∂N a1,i−1 (x) = ∂yi (Pi−1
m=0
ym ,
∂N a2,i−1 (x) = ∂yi (Pi−1
m=0
ym ,
Pi−1
m=0
Pi−1
m=0
Pi−1
m=0
, ym
, ym
, ym
Pi−1
m=0
Pi−1
m=0
Pi−1
m=0
(32) ym )
(33) ym )
(34) ym )
i−1 i−1 i−1 i−1 i−1 i−1 ym , ym , ym − N ym , ym , ym . ri−1 (x) = −L m=0
m=0
m=0
m=0
m=0
(35)
m=0
2.1. Example 1: Nonlinear heat conduction in a slab with variable thermal conductivity The problem of one-dimensional heat conduction in a slab of finite thickness given in Refs. 11 and 14 is considered as d dθ (1 + θ) = 0, (36) dy dy θ(0) = 1,
θ(1) = 0.
The exact solution of the problem is obtained by integrating Eq. (36) once separating variables and integrating again to obtain θ + θ 2 = c1 y + c2 , 2 where c1 , c2 are integration constants. Applying the boundary conditions gives 2 2 c2 = −1 − , c3 = 1 + .
(37) then (38) (37)
(39)
August 8, S1793962311000499
2011 16:37 WSPC/262-IJMSSC/S1793-9623
New Algorithm for Solving Nonlinear BVPs in Heat Transfer
361
Thus, the exact solution of (36) is given as
θ(y) =
−1 +
1 − 2 y − 2y + 2 + 2 .
(40)
To investigate the successive linearization method (SLM) application on (36), the equation is written as L(θ, θ , θ ) + N (θ, θ , θ ) = 0,
(41)
where L(θ, θ , θ ) = θ ,
N (θ, θ , θ ) = θθ + θ2 .
(42)
Using Eqs. (30)–(35), the governing equation to be solved is found as θi + a0,i−1 θi + a1,i−1 θi + a2,i−1 θi = ri−1 (x),
(43)
subject to the boundary conditions θi (0) = 0,
θi (1) = 0,
(44)
where ∂N a0,i−1 (x) = ∂θi (Pi−1
m=0
θm ,
∂N a1,i−1 (x) = ∂θi (Pi−1
m=0 θm ,
∂N a2,i−1 (x) = ∂θi (Pi−1 ri−1 (x) = −
m=0 θm ,
i−1
n=0
θn +
Pi−1
m=0
, θm
Pi−1
m=0 θm ,
Pi−1
m=0 θm ,
i−1 n=0
θn
Pi−1
m=0
θm )
=
m=0 θm )
i−1 n=0
m=0 θm )
θn +
θm ,
(45)
m=0
= 2
Pi−1
Pi−1
i−1
i−1
θm ,
(46)
m=0
=
i−1
θm ,
(47)
m=0
i−1
2
θn
.
(48)
n=0
The initial approximation θ0 (y) must be chosen in such a way that it satisfies the boundary conditions (37). An appropriate initial guess is θ0 (y) = 1 − y.
(49)
August 8, S1793962311000499
362
2011 16:37 WSPC/262-IJMSSC/S1793-9623
S. S. Motsa
Note that Eq. (43) can be written as b1,i−1 θi + b2,i−1 θi + b3,i−1 θi = ri−1 ,
(50)
where b1,i−1 = 1 +
i−1
θn ,
b2,i−1 = 2
n=0
i−1
θn ,
b3,i−1 =
n=0
i−1
θn .
(51)
n=0
Starting from the initial approximation, θ0 , the subsequent solutions for θi , i ≥ 1 are obtained by iteratively by solving Eq. (51) subject the boundary conditions (44). Once each solution for θi (i ≥ 1) has been obtained, the approximate solutions for θ(y) are obtained as θ(y) ≈
i
θn (y),
(52)
n=0
where i is the order of SLM approximation. The coefficient parameters and the righthand side of Eq. (50) for i = 1, 2, 3, . . . , are known (from previous iterations). Thus, Eq. (50) can easily be solved using numerical methods such as finite differences, finite elements, Runge–Kutta based shooting methods or collocation methods. For this particular example, the solution for θi can be found exactly by integrating Eq. (50) subject to the boundary conditions (44) to obtain, θi =
1+
1 i−1
n=0 θm
i−1
1 + (1 − y) − θm − 2 2 n=0
i−1
2 θm
.
(53)
n=0
In this work, other equations similar to (50) are solved using the Chebyshev spectral collocation method (see for example Refs. 22–24). This method is based on approximating the unknown functions by the Chebyshev interpolating polynomials in such a way that they are collocated at the Gauss–Lobatto points defined as ξj = cos
πj , N
j = 0, 1, . . . , N,
(54)
where the variable ξ is defined as ξ = 2y − 1,
(55)
and is used to transform the physical domain of the problem [0, 1] to the space [−1, 1] on which the Chebyshev spectral collocation method can be implemented.
August 8, S1793962311000499
2011 16:37 WSPC/262-IJMSSC/S1793-9623
New Algorithm for Solving Nonlinear BVPs in Heat Transfer
363
The Chebyshev spectral collocation method involves approximating the unknown function θ(y) as a truncated series of Chebyshev polynomials of the form θi ≈
N
θi (ξk )Tk (ξj ),
j = 0, 1, . . . , N,
(56)
k=0
where Tk is the kth Chebyshev polynomial and N + 1 is the number of collocation points. Derivatives of the functions θi (y) at the collocation points are represented as N
ds θi s = Dkj θk (ξj ), s dy
(57)
k=0
where s is the order of differentiation and D = 2D, with D being the Chebyshev spectral differentiation matrix whose entries are defined as cj (−1)j+k ck z j − z k zk =− 2(1 − zk2 )
Djk =
j = k; j, k = 0, 1, . . . , N,
Dkk
k = 1, 2, . . . , N − 1,
D00 =
(58)
2N 2 + 1 = −DN N . 6
Applying the Chebyshev spectral method, by substituting Eqs. (56)–(58) in (50) gives Ai−1 Gi = Ri−1 ,
(59)
θi (ξ0 ) = θi (ξN ) = 0,
(60)
subject to
where Ai−1 = a1,i−1 D2 + a2,i−1 D + a3,i−1 , T
Gi = [θi (ξ0 ), θi (ξ1 ), . . . , θi (ξN )] , Ri−1 = [ri−1 (ξ0 ), ri−1 (ξ1 ), . . . , ri−1 (ξN )]T .
(61) (62) (63)
In the above definitions T stands for transpose and ap,i−1 (p = 1, 2, 3) denotes a diagonal matrix of size (N + 1) × (N + 1). The boundary conditions (60) are implemented by deleting the first and last rows and columns of Ai−1 , and deleting the first and last rows of Gi−1 and Ri−1 . The solutions for θi (ξ1 ), θi (ξ2 ), . . . , θi (ξN −1 ) are then obtained from Gi = A−1 i−1 Ri−1 .
(64)
August 8, S1793962311000499
364
2011 16:37 WSPC/262-IJMSSC/S1793-9623
S. S. Motsa
2.2. Example 2: Nonlinear convective–conduction fin problem The problem of one-dimensional heat conduction in convective straight fin with temperature dependent conductivity that was discussed in Refs. 1, 3, 5, 7, 9, 10, 13, 15 and 18, is considered and is given by dθ d (65) (1 + θ) − β 2 θ = 0, dy dy subject to θ (0) = 1,
θ(1) = 1.
(66)
Before applying the SLM, an exact solution for Eq. (65) is found. The transformation is introduced d2 θ du du =u , = dy 2 dy dθ
dθ = u, dy
(67)
and substituted in (65) to obtain (u2 − β 2 θ)dθ + (u + θu)du = 0.
(68)
Multiplying (68) by the integrating factor ρ = 1 + θ gives (u2 − β 2 θ)(1 + θ)dθ + u(1 + θ)2 du = 0, which is an exact differential equation with the following solution: 1 θ 1 2 θ 2 2 2 2 + + u θ + −β θ = c1 , 2 2 2 3 where c1 is an integration constant. From (67) and (70) it is found that 2θ 2 2 c2 + β θ 1 + 3 dθ = , dy 1 + θ
(69)
(70)
(71)
where c2 = 2c1 . The above equation is difficult to integrate, but it can be shown, using long division of the right-hand side of (71), that for the special case when and β are carefully selected in such a way that c2 = −
β2 , 32
(72)
the problem then reduces to (β 2 /32 )(2θ − 1)(1 + θ)2 (2θ − 1) dθ = = . dy 1 + θ 3β 2
(73)
It can be shown that the solution of the above equation is θ(y) =
β3y2 + 3 , 6
whenever =
β2 + 3 . 6
(74)
August 8, S1793962311000499
2011 16:37 WSPC/262-IJMSSC/S1793-9623
New Algorithm for Solving Nonlinear BVPs in Heat Transfer
365
For example, for the special case when = 2 and β = 3, the exact solution of (65) is θ(y) =
3y 2 + 1 . 4
(75)
When the SLM is applied on (65), using the approach described in the previous subsection, one obtains b1,i−1 θi + b2,i−1 θi + b3,i−1 θi = ri−1 ,
(76)
where b1,i−1 = 1 +
i−1
θn ,
b2,i−1 = 2
n=0
ri−1 = −
i−1
θn ,
b3,i−1 =
n=0
θn +
n=0
i−1
i−1 n=0
θn
i−1
θn +
n=0
i−1
θn − β 2 ,
(77)
n=0
i−1
2 θn
i−1
− β2
n=0
θn .
(78)
n=0
One chooses θ0 (y) =
cosh βy , cosh β
(79)
as an initial approximation. Thus, starting from the initial approximation, θ0 , the subsequent solutions for θi , i ≥ 1 are obtained by iteratively solving Eq. (76) subject to the boundary conditions θi (0) = 0,
θi (1) = 0.
(80)
Applying the Chebyshev spectral method in (76) and (80) gives Ai−1 Gi = Ri−1 ,
(81)
subject to N
DN k θi (ξk ) = 0, θi (ξ0 ) = 0,
(82)
k=0
where T
Ri−1 = [ri−1 (ξ0 ), ri−1 (ξ1 ), . . . , ri−1 (ξN )] ,
(83)
and Ai−1 and Gi−1 have the same form as Eqs. (61) and (62). After incorporating the boundary conditions (82) in (81), the solution Gi can be obtained from Gi = A−1 i−1 Ri−1 .
(84)
August 8, S1793962311000499
366
2011 16:37 WSPC/262-IJMSSC/S1793-9623
S. S. Motsa
2.3. Example 3: Nonlinear convective–radiative–conduction fin equation In this section, the method of successive linearization is applied to the nonlinear convective–radiative–conduction problem with temperature-dependent thermal conductivity that was investigated in Ref. 19. The problem extends the convection problem described in Example 2 to include a radiation term (β1 ) and the governing equations are given by d dθ (85) (1 + θ) − β 2 θ − β1 θ4 = 0, dy dy subject to θ (0) = 1,
θ(1) = 1.
(86)
When the SLM is applied on (85), using the approach described in the previous subsection, one obtains b1,i−1 θi + b2,i−1 θi + b3,i−1 θi = ri−1 ,
(87)
where b1,i−1 = 1 +
i−1
θn ,
b2,i−1 = 2
n=0
b3,i−1 =
i−1
θn − β 2 − 4β1
ri−1 = −
i−1
n=0
θn ,
(88)
n=0
i−1
n=0
i−1
3 θn
,
(89)
n=0
θn +
i−1
θn
n=0
i−1 n=0
θn +
i−1 n=0
2 θn
− β2
i−1 n=0
θ n − β1
i−1
4 θn .
n=0
(90) The rest of the analysis used to obtain θ(y) is exactly the same as the analysis given in the previous section. 3. Results and Discussion In this section, the results showing the temperature distribution and rate of heat transfer at the boundaries for the problems discussed in Examples 1 and 2 are presented. To assess the accuracy of the proposed SLM, comparison is made with some of the exact analytical solutions and the numerical solutions obtained using the MATLAB bvp4c routine that is used for solving boundary value problems. In generating the results presented in this study, N + 1 = 60 collocation points were used in the proposed SLM. Table 1 gives a comparison between the successive linearization results and the exact results for θ(y) at selected values of y obtained using (40) for different
August 8, S1793962311000499
2011 16:37 WSPC/262-IJMSSC/S1793-9623
New Algorithm for Solving Nonlinear BVPs in Heat Transfer
367
Table 1. Comparison between the successive linearization results and the exact results for θ(y) at selected values of y.
y
0th order
1st order
2nd order
3rd order
4th order
Exact
0.5
0.2 0.4 0.6 0.8
0.8000000 0.6000000 0.4000000 0.2000000
0.8285714 0.6461538 0.4500000 0.2363636
0.8284271 0.6457513 0.4494898 0.2360680
0.8284271 0.6457513 0.4494897 0.2360680
0.8284271 0.6457513 0.4494897 0.2360680
0.8284271 0.6457513 0.4494897 0.2360680
2
0.2 0.4 0.6 0.8
0.8000000 0.6000000 0.4000000 0.2000000
0.8615385 0.7090909 0.5333333 0.3142857
0.8601478 0.7041695 0.5247312 0.3062657
0.8601471 0.7041595 0.5246951 0.3062258
0.8601471 0.7041595 0.5246951 0.3062258
0.8601471 0.7041595 0.5246951 0.3062258
5
0.2 0.4 0.6 0.8
0.8000000 0.6000000 0.4000000 0.2000000
0.8800000 0.7500000 0.6000000 0.4000000
0.8770370 0.7381579 0.5750000 0.3666667
0.8770330 0.7380832 0.5745968 0.3656863
0.8770330 0.7380832 0.5745967 0.3656854
0.8770330 0.7380832 0.5745967 0.3656854
values of . It is noted from the table that the results from the two methods are in good agreement with accuracy up to 7 decimal places being achieved at third-order of approximation. This shows that the proposed linearization method converges rapidly and is very accurate. The temperature profile for nonlinear conduction problem described in Example 1 is depicted in Fig. 2 for different values of . The figure gives a comparison between the exact solution and the second-order successive linearization solution. It can be
1 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2
5 2 0.5
0.1 0
0
0.2
0.4
0.6
0.8
1
Fig. 2. The comparison of the exact solution (solid line) and second-order successive linearization solution (squares) for θ(y) using = 0.5, = 2, and = 5.
August 8, S1793962311000499
368
2011 16:37 WSPC/262-IJMSSC/S1793-9623
S. S. Motsa
Table 2. Comparison of the successive linearization results for θ (0) against the HAM results of Ref. 11 and the exact solution.
Ref. 11 (HAM)
Present results θ (0)
Exact θ (0)
Order of approximation
θ (0)
Order of approximation
0.5
1 5 10 15 20
−0.75 −0.828125 −0.833496 −0.833328 −0.833333
1 2 3 4 5
−0.83333333 −0.83333333 −0.83333333 −0.83333333 −0.83333333
−0.83333333
2
1 5 10 15 20
−0.75 −0.666992 −0.666668 −0.666667 −0.666667
1 2 3 4 5
−0.66666667 −0.66666667 −0.66666667 −0.66666667 −0.66666667
−0.66666667
5
1 5 10 15 20
−0.75 −0.5876 −0.583377 −0.583333 −0.583333
1 2 3 4 5
−0.58333333 −0.58333333 −0.58333333 −0.58333333 −0.58333333
−0.58333333
seen from the figure that there is good agreement between the two results. This shows the accuracy of the linearization method for different values of . Table 2 gives a comparison of the results for the rate of heat transfer θ (0), for the conduction problem from the current linearization approach versus the HAM of Ref. 11 and the exact solution. It can be seen from the table that the exact solution matches the current linearization results even at first-order approximation. A close approximation to the exact solution that is accurate to 6 decimal places is obtained at the 15th order using the HAM approach of Ref. 11. Consequently, it can be concluded that the current method is efficient and converges more rapidly than the standard HAM. When the solution (53) is differentiated and evaluated at y = 0, one obtains θ (0) = −
1 + (/2) . 1+
(91)
The present results given in Table 2 match the results obtained using Eq. (91). Table 3 gives a comparison between the successive linearization results and the exact results, given by Eq. (75) for θ(y) when β = 3 and = 2 in Example 2. It can be seen from the table that there is very good agreement between the exact results and the results from the proposed method. Accuracy to eight decimal places is achieved at only the third-order of approximation. This shows that the linearization method converges rapidly and is very accurate. The performance of a fin is characterized by fin efficiency. This is defined as the ratio of the fin heat transfer rate to the heat transfer rate of the entire surface if
August 8, S1793962311000499
2011 16:37 WSPC/262-IJMSSC/S1793-9623
New Algorithm for Solving Nonlinear BVPs in Heat Transfer
369
Table 3. Comparison between the successive linearization results and the exact results for θ(y), when β = 3 and = 2. y 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0
0th order
1st order
2nd order
3rd order
4th order
Exact
0.09932793 0.10383131 0.11774980 0.14234550 0.17984866 0.23365997 0.30865887 0.41164604 0.55196004 0.74232415 1.00000000
0.24856397 0.25678080 0.28139952 0.32225212 0.37885407 0.45023267 0.53495164 0.63153972 0.73944053 0.86039951 1.00000000
0.24997862 0.25747719 0.27997690 0.31748753 0.37001668 0.43755717 0.52008247 0.61756716 0.73002322 0.85749729 1.00000000
0.25000000 0.25750000 0.28000000 0.31750000 0.37000000 0.43750000 0.52000000 0.61750000 0.73000000 0.85750000 1.00000000
0.25000000 0.25750000 0.28000000 0.31750000 0.37000000 0.43750000 0.52000000 0.61750000 0.73000000 0.85750000 1.00000000
0.25000000 0.25750000 0.28000000 0.31750000 0.37000000 0.43750000 0.52000000 0.61750000 0.73000000 0.85750000 1.00000000
the entire fin were at the base temperature. Mathematically the fin efficiency ηe is given as (see for example Ref. 15), 1 + dθ ηe = . (92) β 2 dy y=1 Table 4 gives a comparison between the first four terms of the successive linearization approximation and the MATLAB numerical results for the fin efficiency for different values of β and in Example 2. The results from the table illustrate the Table 4. Comparison between the successive linearization results and MATLAB bvp4c for fin efficiency at different values of β and .
β
0th order
1st order
2nd order
3rd order
4th order
bvp4c
0.5
1 2 3 5 10
0.76159416 1.92805516 2.98516426 4.99954602 9.99999996
0.81994146 0.54943770 0.38245270 0.23102934 0.11552930
0.81939437 0.54898744 0.38223305 0.23091059 0.11547006
0.81939431 0.54898742 0.38223304 0.23091058 0.11547005
0.81939431 0.54898742 0.38223304 0.23091058 0.11547005
0.81939431 0.54898742 0.38223304 0.23091058 0.11547005
1.0
1 2 3 5 10
0.76159416 1.92805516 2.98516426 4.99954602 9.99999996
0.85838611 0.60596706 0.42745886 0.25882547 0.12943497
0.85593408 0.60326040 0.42616018 0.25815248 0.12909964
0.85593219 0.60325899 0.42615958 0.25815209 0.12909944
0.85593219 0.60325899 0.42615958 0.25815209 0.12909944
0.85593219 0.60325899 0.42615958 0.25815209 0.12909944
2.0
1 2 3 5 10
0.76159416 1.92805516 2.98516426 4.99954602 9.99999996
0.90613848 0.69774421 0.50674828 0.30859308 0.15434041
0.89834564 0.68437834 0.50001745 0.30540648 0.15275824
0.89831798 0.68431983 0.50000000 0.30539504 0.15275252
0.89831798 0.68431983 0.50000000 0.30539504 0.15275252
0.89831798 0.68431983 0.50000000 0.30539504 0.15275252
5.0
1 2 3 5 10
1.92805516 2.98516426 2.98516426 4.99954602 9.99999996
0.96728423 0.87946712 0.69768102 0.43440941 0.21737025
0.94683071 0.81286612 0.65035214 0.41584670 0.20838577
0.94657973 0.81041783 0.64943418 0.41540397 0.20816658
0.94657973 0.81041783 0.64943418 0.41540397 0.20816658
0.94657973 0.81041783 0.64943476 0.41540397 0.20816658
August 8, S1793962311000499
2011 16:37 WSPC/262-IJMSSC/S1793-9623
S. S. Motsa
370 1
1 0.9
0.95
1
0.8 0.7
5
0.9
0.6
0.85
0.5
2
2
0.4
0.8
0.3 0.2
0.75
0.5 0.7
0
0.2
0.1
0.4
0.6
0.8
1
0
5 0
0.2
0.4
0.6
0.8
1
Fig. 3. Comparison of the second-order successive linearization and bvp4c results when varying (with β = 1) and when varying β (with = 1).
fast convergence of the proposed method. Accuracy to 8 decimal places is achieved at only the third-order approximation. Figure 3 shows the temperature distribution in the fin when varying and β. It is noted from the figure that increasing the conductivity parameter results in an increase in temperature and increasing the convection parameter β results in a decrease in the temperature within the fin. This trend is consistent with the observation made in Refs. 1, 5, 9, 10, 13, 15 and 18 where a similar problem was studied using other methods of solution. The figure also demonstrates that there is excellent agreement between the second-order successive linearization results and the MATLAB bvp4c results. This shows the accuracy and rapid convergence of the linearization method. Table 5 shows the dimensionless temperature at θ(0) for different values of the radiation parameter β1 when β = 1 and = 1. It is noted from the table that the temperature decreases with an increase in radiation. It is only clear from the table that the successive linearization solutions converge to the numerical value. The effect of varying the radiation parameter on the temperature distribution is illustrated graphically in Fig. 4 for β = 1 and = 1 in Example 3. Again, it is noted that there is good agreement between the successive linearization results and the MATLAB numerical results. Table 5.
The dimensionless temperature at y = 0 for different values of β1 with = 1 and β = 1.
β1
0th order
2nd order
3rd order
4th order
6th order
Numerical
1 5 10 20 50 100
0.76159416 0.76159416 0.76159416 0.76159416 0.76159416 0.76159416
0.7063239 0.5718853 0.4986237 0.4329974 0.3824607 0.3699939
0.7063211 0.5718175 0.4968401 0.4206927 0.3336657 0.2942057
0.7063211 0.5718175 0.4968356 0.4204084 0.3263561 0.2668969
0.7063211 0.5718175 0.4968356 0.4204083 0.3262285 0.2643166
0.7063211 0.5718173 0.4968356 0.4204083 0.3262285 0.2643166
August 8, S1793962311000499
2011 16:37 WSPC/262-IJMSSC/S1793-9623
New Algorithm for Solving Nonlinear BVPs in Heat Transfer
371
1 0.9 0.8
1
0.7 0.6
10
0.5 0.4 0.3
100 0.2
0
0.2
0.4
0.6
0.8
1
Fig. 4. The comparison of the numerical solution (solid line) and fourth-order successive linearization solution (squares) for θ(y) using = 1, β = 1 for different values of β1 .
4. Conclusion A new method for solving two-point boundary-value problems has been presented. Using examples of highly nonlinear problems arising in heat transfer, it was clearly illustrated that the proposed SLM is very accurate and converges very rapidly to the exact solutions, which are reported for the first time in this work. Some of the advantages of this method are that solutions are obtained iteratively by integrating linear differential equations. Another desirable feature of the proposed method is that it reduces the volumes of calculations involved because it requires only a few iterations to converge to the exact results. The proposed method is also very simple and straightforward to implement. Because of its efficiency, the proposed method has great potential to be used in solving nonlinear equations in place of traditional numerical methods such as finite differences, shooting methods, and collocation methods. Acknowledgment The authors would like to thank the anonymous reviewers whose comments contributed significantly towards the improvement of this work. References 1. Arslanturk C., A decomposition method for fin efficiency of convective straight fins with temperature-dependent thermal conductivity, Int. Commun. Heat Mass Trans. 32:831–841, 2005.
August 8, S1793962311000499
372
2011 16:37 WSPC/262-IJMSSC/S1793-9623
S. S. Motsa
2. Chang M.-H. A., A decomposition solution for fins with temperature dependent surface heat flux, Int. J. Heat Mass Trans. 48:673–682, 2005. 3. Chiu C. H., Chen C. K., A decomposition method for solving the convective longitudinal fins with variable thermal conductivity, Int. J. Heat Mass Trans. 45:2067–2075, 2002. 4. Lesnic D., Heggs P. J., A decomposition method for power-law fin type problems, Int. Commun. J. Heat Mass Trans. 31:673–682, 2004. 5. Joneidi A. A., Ganji D. D., Babaelahi M., Differential transformation method to determine fin efficiency of convective straight fins with temperature dependent thermal conductivity, Int. Commun. Heat Mass Trans. 36:757–762, 2009. 6. Abbasbandy S., The application of the homotopy analysis method to nonlinear equations arising in heat transfer, Phys. Lett. A 360:109–113, 2006. 7. Domairry G., Fazeli M., Homotopoy analysis method to determine the fin efficiency of convective straight fins with temperature-dependent thermal conductivity, Commun. Nonlinear Sci. Numer. Simulat. 14:489–499, 2009. 8. Domairry G., Bararnia H., An approximation of the analytic solutions of some nonlinear heat transfer equations: A survery using homotopy analysis method, Adv. Studies Theor. Phys. 2:507–518, 2008. 9. Khani F., Ahmadzadeh M., Hamedi-Nezhad S., A series solution of the fin problem with temperature-dependent thermal conductivity, Commun. Nonlinear Sci. Numer. Simulat. 14:3007–3017, 2009. 10. Mustafa I., Application of homotopy analysis method for fin efficiency of convective straight fins with temperature-dependent thermal conductivity, Math. Comput. Simulat. 79:189–200, 2008. 11. Sajid M., Hayat T., Comparison of HAM and HPM methods in nonlinear heat conduction and convection equations, Nonlinear Anal. Real World Appl. 9:2296–2301, 2008. 12. Hosseini M. J., Gorji M., Ghanbarpour M., Solution of temperature distribution in a radiating fin using homotopy perturbatiod method, Math. Problem Eng. DOI:10.1155/2009/831362, 2009. 13. Rajabi A., Homotopy perturbation method for fin efficiency of convective straight fins with temperature-dependent thermal conductivity, Phys. Lett. A 364:33–37, 2007. 14. Rajabi A., Ganji D. D., Taherian H., Application of homotopy perturbation method in nonlinear heat conduction and convection equations, Phys. Lett. A 360:570–573, 2007. 15. Kim S., Huang C. H., A series solution of the fin problem with a temperaturedependent thermal conductivity, J. Phys. D: Appl. Phys. 39:4894–4901, 2006. 16. Kim S., Huang C. H., A series solution of the nonlinear fin problem with a temperature-dependent thermal conductivity and heat transfer coeffient, J. Phys. D: Appl. Phys. 40:2979–2987, 2006. 17. Ganji D. D., Hosseini M. J., Shayegh J., Some nonlinear heat transfer equations solved by three approximate methods, Int. Commun. Heat Mass Trans. 34:1003–1016, 2007. 18. Ganji D. D., Afrouzi G. A., Talarposhti F., Application of variational iteration method and homotopy-perturbation method for nonlinear heat diffusion and heat transfer equations, Phys. Lett. A 368:450–457, 2007. 19. Miansari Mo., Ganji D. D., Miansari Me., Application of He’s variational iteration method to nonlinear heat transfer equations, Phys. Lett. A 372:779–785, 2008. 20. Tari H., Ganji D. D., Babazadeh H., The application of He’s variational iteration method to nonlinear equations arising in heat transfer, Phys. Lett. A 363:213–217, 2007.
August 8, S1793962311000499
2011 16:37 WSPC/262-IJMSSC/S1793-9623
New Algorithm for Solving Nonlinear BVPs in Heat Transfer
373
21. Kierzenka J., Shampine L. F., A BVP solver based on residual control and the MATLAB pse, ACM Trans. Math. Software 27(3):299–316, 2001. 22. Canuto C., Hussaini M. Y., Quarteroni A., Zang T. A., Spectral Methods in Fluid Dynamics, Springer-Verlag, Berlin, 1988. 23. Don W. S., Solomonoff A., Accuracy and speed in computing the Chebyshev collocation derivative SIAM, J. Sci. Comput. 16:1253–1268, 1995. 24. Trefethen L. N., Spectral Methods in MATLAB, SIAM, 2000.
2011 16:37 WSPC/262-IJMSSC/S1793-9623
International Journal of Modeling, Simulation, and Scientific Computing Vol. 2, No. 3 (2011) 355–373 c World Scientific Publishing Company DOI: 10.1142/S1793962311000499
A NEW ALGORITHM FOR SOLVING NONLINEAR BOUNDARY VALUE PROBLEMS ARISING IN HEAT TRANSFER
S. S. MOTSA Mathematics Department University of Swaziland, Private Bag 4 Kwaluseni, M201, Swaziland [email protected] Received 27 May 2010 Revised 18 April 2011 Accepted 30 April 2011
In this paper, a very efficient and easy-to-use successive linearization approach for solving nonlinear differential equations is proposed. The implementation of the method is demonstrated by solving three nonlinear differential equations of different complexities arising in heat transfer. New closed form explicit analytic solutions of some of the governing nonlinear equations are obtained and compared with the results of the proposed method and with numerical solutions from the MATLAB in-built routine bvp4c. Keywords: Successive collocation.
linearization
method;
convective–radiative
fin;
spectral
1. Introduction Rectangular fins are often used to enhance heat exchange between heated surfaces and the surrounding environment. In recent years, the problem of heat transfer through extended fins has attracted the interest of many researchers in view of its numerous applications on heat exchangers, electronic components, semiconductors, solar collectors, and many other engineering applications. For most practical applications, the mathematical model of the underlying heat transfer problem consists of nonlinear differential equations whose closed form analytical solutions are either very difficult or impossible to obtain. As a result many researchers have recently devoted attention toward developing new analytical or numerical approximation algorithms or improving existing methods for solving these problems. In the last decade, several methods for solving nonlinear equations, including the Adomian decomposition method,1–4 differential transform method,5 homotopy analysis method,6–11 and homotopy perturbation method,12–14 Taylor series15,16 355
August 8, S1793962311000499
356
2011 16:37 WSPC/262-IJMSSC/S1793-9623
S. S. Motsa
and variational iteration method,17–20 have been applied to solve nonlinear heat transfer equations. In this work, a new numerical approach is proposed that is based on successive linearization of the governing equations. The applicability of the method is tested on three different nonlinear heat transfer equations of different complexity. One of the test problems is that of steady conduction in a slab with variable thermal conductivity that was previously discussed in Ref. 14 using the homotopy perturbation method (HPM) and later studied in Ref. 11 using the homotopy analysis method (HAM). It was proved in Ref. 11 that the HPM approach used in Ref. 14 is only valid for small thermal conductivity parameters (0 ≤ < 1). The second problem considered in this work is the nonlinear convective straight fin problem with temperature dependent thermal conductivity that was investigated in Refs. 5, 9, 10, 13, 15, 16 and 18. The new linearization approach was also applied on the nonlinear convective– radiative conduction problem with temperature dependent thermal conductivity that was investigated in Ref. 19 using the variational iteration method. Here, the nonlinear problems are solved to obtain the distribution of temperature using the new linearization approach and compared with results obtained using the Matlab ODE solver bvp4c, which is a finite difference approach that implements the 3-stage Lobatto III formula (see for example Ref. 21). Simple closed form analytic solutions of some of the governing nonlinear equations, which are subsequently used to validate the results of the proposed linearization method are also reported. 2. Successive Linearization Method Solution In this section, the basic idea of the proposed method of successive linearization is described and it is implemented in three test problems arising in nonlinear heat transfer. It is begun by describing the basic idea behind this method of solution. Consider a nonlinear boundary value problem of the form L[y(x), y (x), y (x)] + N [y(x), y (x), y (x)] = 0,
(1)
where y(x) is an unknown function, x is an independent variable and the primes denote ordinary differentiation with respect to x. The functions L and N represent the linear and nonlinear components of the governing equation, respectively. For illustrative purposes, Eq. (1) is to be solved for x ∈ [a, b] subjected to the boundary conditions y(a) = ya ,
y(b) = yb ,
(2)
where ya and yb are given constants. For the initial guess of the solution of Eq. (1), a function that satisfies the boundary conditions (2) is chosen. It is convenient to consider polynomial functions. Thus, a suitable initial approximation, denoted by y0 (x), that satisfies Eq. (2) would be a straight line.
August 8, S1793962311000499
2011 16:37 WSPC/262-IJMSSC/S1793-9623
New Algorithm for Solving Nonlinear BVPs in Heat Transfer
357
A function Y1 (x) is defined to represent the vertical difference between y(x) and the initial guess y0 (x) (see Fig. 1), that is Y1 (x) = y(x) − y0 (x),
or y(x) = y0 (x) + Y1 (x).
(3)
Substituting Eq. (3) in (1) gives L[Y1 , Y1 , Y1 ] + N [y0 + Y1 , y0 + Y1 , y0 + Y1 ] = −L[y0 , y0 , y0 ].
(4)
Since y0 (x) is given, solving Eq. (4) would yield an exact solution for Y1 (x). However, since the equation is nonlinear, it may not be possible to find an exact solution. Therefore, an approximate solution is needed, which is obtained by solving the linear part of the equation assuming that Y1 and its derivatives are small. This assumption enables us to use the Taylor series method to linearize the equation. If Y1 (x) is the solution of the full Eq. (4), let y1 (x) denote the solution of the linearized version of (4). Expanding (4) using Taylor series and neglecting higher order terms gives ∂N ∂N ∂N + y1 + y1 L[y1 , y1 , y1 ] + y1 ∂y1 (y0 ,y ,y ) ∂y1 (y0 ,y ,y ) ∂y1 (y0 ,y ,y ) 0
=
−L[y0 , y0 , y0 ]
−
0
0
0
N [y0 , y0 , y0 ].
0
0
(5)
Equation (5) can be written in compact form as L[y1 , y1 , y1 ] + a0,0 y1 + a1,0 y1 + a2,0 y1 = r0 (x),
(6)
subject to the boundary conditions y1 (a) = 0,
y1 (b) = 0,
(7)
where ∂N a0,0 (x) = ∂y1 (y0 ,y ,y ) 0 0 ∂N a1,0 (x) = ∂y1 (y0 ,y ,y ) 0 0 ∂N a2,0 (x) = ∂y1
(8)
(9)
(10)
(y0 ,y0 ,y0 )
r0 (x) = −L[y0 , y0 , y0 ] − N [y0 , y0 , y0 ].
(11)
Since the right-hand side of Eq. (6) is known and the left-hand side is linear, the equation can be solved for y1 (x). Assuming that the solution of the linear part (6)
August 8, S1793962311000499
358
2011 16:37 WSPC/262-IJMSSC/S1793-9623
S. S. Motsa
is close to the solution of Eq. (4), that is, Y1 (x) ≈ y1 (x), the current estimate (first-order) of the solution y(x) is y(x) ≈ y0 (x) + y1 (x).
(12)
To improve on this solution, a slack function, Y2 (x) is defined that when added to y1 (x) gives Y1 (x) (see Fig. 1), that is Y1 (x) = Y2 (x) + y1 (x).
(13)
Since y1 (x) is now known (as a solution of Eq. 6), Eq. (13) is substituted in Eq. (4) to obtain L[Y2 , Y2 , Y2 ] + N [y0 + y1 + Y2 , y0 + y1 + Y2 , y0 + y1 + Y2 ] = −L[y0 + y1 , y0 + y1 , y0 + y1 ].
(14)
Solving Eq. (14) would result in an exact solution for Y2 (x). But, since the equation is nonlinear, it may not be possible to find an exact solution. Therefore, the equation is linearized using Taylor series expansion and the resulting linear equation is solved. The solution of the linear version of Eq. (14) is denoted by y2 (x), such that Y2 (x) ≈ y2 (x). Setting Y2 (x) = y2 (x) and expanding Eq. (14), for small y2 (x) and its derivatives gives L[y2 , y2 , y2 ] + a0,1 y2 + a1,1 y2 + a2,1 y2 = r1 (x),
(15)
subject to the boundary conditions y2 (a) = 0,
y2 (b) = 0,
(16)
where ∂N a0,1 (x) = ∂y2 (y0 +y1 ,y +y ,y +y ) 0 1 0 1 ∂N a1,1 (x) = ∂y2 (y0 +y1 ,y +y ,y +y ) 0 1 0 1 ∂N a2,1 (x) = ∂y2 (y0 +y1 ,y +y ,y +y ) 0
1
0
(17)
(18)
(19)
1
r1 (x) = −L[y0 + y1 , y0 + y1 , y0 + y1 ] − N [y0 + y1 , y0 + y1 , y0 + y1 ].
(20)
After solving (15), the current (second-order) estimate of the solution y(x) is y(x) ≈ y0 (x) + y1 (x) + y2 (x).
(21)
August 8, S1793962311000499
2011 16:37 WSPC/262-IJMSSC/S1793-9623
New Algorithm for Solving Nonlinear BVPs in Heat Transfer
359
Next, Y3 (x) is defined (see Fig. 1) such that Y2 (x) = Y3 (x) + y2 (x).
(22)
Equation (22) is substituted in the nonlinear equation (14) and the linearization process described above is repeated. This process is repeated for m = 3, 4, 5, . . . , i. In general, one have Yi (x) = Yi+1 (x) + yi (x).
(23)
Thus y(x) is obtained as y(x) = Y1 (x) + y0 (x),
(24)
= Y2 (x) + y1 (x) + y0 (x),
(25)
= Y3 (x) + y2 (x) + y1 (x) + y0 (x), .. .
(26)
= Yi+1 (x) + yi (x) + · · · + y3 (x) + y2 (x) + y1 (x) + y0 (x),
(27)
= Yi+1 (x) +
i
ym (x).
m=0
The procedure for obtaining each Yi (x) is illustrated in Fig. 1 for i = 1, 2, 3.
Fig. 1.
Geometric representation of the functions Yi (i = 1, 2, 3).
(28)
August 8, S1793962311000499
360
2011 16:37 WSPC/262-IJMSSC/S1793-9623
S. S. Motsa
It is noted that when i becomes large, Yi+1 becomes increasingly smaller. Thus, for large i, the ith order solution of y(x) can be approximated by i
y(x) =
ym (x) = yi (x) +
m=0
i−1
ym (x).
(29)
m=0
Starting from a known initial guess y0 (x), the solutions for yi (x) can be obtained by successively linearizing the governing equation (1) and solving the resulting linear equation for yi (x) given that the previous guess yi−1 (x) is known. The general form of the linearized equation that is successively solved for yi (x) is given by L[yi , yi , yi ] + a0,i−1 yi + a1,i−1 yi + a2,i−1 yi = ri−1 (x),
(30)
subject to the boundary conditions yi (a) = 0,
yi (b) = 0,
(31)
where
∂N a0,i−1 (x) = ∂yi (Pi−1
m=0
ym ,
∂N a1,i−1 (x) = ∂yi (Pi−1
m=0
ym ,
∂N a2,i−1 (x) = ∂yi (Pi−1
m=0
ym ,
Pi−1
m=0
Pi−1
m=0
Pi−1
m=0
, ym
, ym
, ym
Pi−1
m=0
Pi−1
m=0
Pi−1
m=0
(32) ym )
(33) ym )
(34) ym )
i−1 i−1 i−1 i−1 i−1 i−1 ym , ym , ym − N ym , ym , ym . ri−1 (x) = −L m=0
m=0
m=0
m=0
m=0
(35)
m=0
2.1. Example 1: Nonlinear heat conduction in a slab with variable thermal conductivity The problem of one-dimensional heat conduction in a slab of finite thickness given in Refs. 11 and 14 is considered as d dθ (1 + θ) = 0, (36) dy dy θ(0) = 1,
θ(1) = 0.
The exact solution of the problem is obtained by integrating Eq. (36) once separating variables and integrating again to obtain θ + θ 2 = c1 y + c2 , 2 where c1 , c2 are integration constants. Applying the boundary conditions gives 2 2 c2 = −1 − , c3 = 1 + .
(37) then (38) (37)
(39)
August 8, S1793962311000499
2011 16:37 WSPC/262-IJMSSC/S1793-9623
New Algorithm for Solving Nonlinear BVPs in Heat Transfer
361
Thus, the exact solution of (36) is given as
θ(y) =
−1 +
1 − 2 y − 2y + 2 + 2 .
(40)
To investigate the successive linearization method (SLM) application on (36), the equation is written as L(θ, θ , θ ) + N (θ, θ , θ ) = 0,
(41)
where L(θ, θ , θ ) = θ ,
N (θ, θ , θ ) = θθ + θ2 .
(42)
Using Eqs. (30)–(35), the governing equation to be solved is found as θi + a0,i−1 θi + a1,i−1 θi + a2,i−1 θi = ri−1 (x),
(43)
subject to the boundary conditions θi (0) = 0,
θi (1) = 0,
(44)
where ∂N a0,i−1 (x) = ∂θi (Pi−1
m=0
θm ,
∂N a1,i−1 (x) = ∂θi (Pi−1
m=0 θm ,
∂N a2,i−1 (x) = ∂θi (Pi−1 ri−1 (x) = −
m=0 θm ,
i−1
n=0
θn +
Pi−1
m=0
, θm
Pi−1
m=0 θm ,
Pi−1
m=0 θm ,
i−1 n=0
θn
Pi−1
m=0
θm )
=
m=0 θm )
i−1 n=0
m=0 θm )
θn +
θm ,
(45)
m=0
= 2
Pi−1
Pi−1
i−1
i−1
θm ,
(46)
m=0
=
i−1
θm ,
(47)
m=0
i−1
2
θn
.
(48)
n=0
The initial approximation θ0 (y) must be chosen in such a way that it satisfies the boundary conditions (37). An appropriate initial guess is θ0 (y) = 1 − y.
(49)
August 8, S1793962311000499
362
2011 16:37 WSPC/262-IJMSSC/S1793-9623
S. S. Motsa
Note that Eq. (43) can be written as b1,i−1 θi + b2,i−1 θi + b3,i−1 θi = ri−1 ,
(50)
where b1,i−1 = 1 +
i−1
θn ,
b2,i−1 = 2
n=0
i−1
θn ,
b3,i−1 =
n=0
i−1
θn .
(51)
n=0
Starting from the initial approximation, θ0 , the subsequent solutions for θi , i ≥ 1 are obtained by iteratively by solving Eq. (51) subject the boundary conditions (44). Once each solution for θi (i ≥ 1) has been obtained, the approximate solutions for θ(y) are obtained as θ(y) ≈
i
θn (y),
(52)
n=0
where i is the order of SLM approximation. The coefficient parameters and the righthand side of Eq. (50) for i = 1, 2, 3, . . . , are known (from previous iterations). Thus, Eq. (50) can easily be solved using numerical methods such as finite differences, finite elements, Runge–Kutta based shooting methods or collocation methods. For this particular example, the solution for θi can be found exactly by integrating Eq. (50) subject to the boundary conditions (44) to obtain, θi =
1+
1 i−1
n=0 θm
i−1
1 + (1 − y) − θm − 2 2 n=0
i−1
2 θm
.
(53)
n=0
In this work, other equations similar to (50) are solved using the Chebyshev spectral collocation method (see for example Refs. 22–24). This method is based on approximating the unknown functions by the Chebyshev interpolating polynomials in such a way that they are collocated at the Gauss–Lobatto points defined as ξj = cos
πj , N
j = 0, 1, . . . , N,
(54)
where the variable ξ is defined as ξ = 2y − 1,
(55)
and is used to transform the physical domain of the problem [0, 1] to the space [−1, 1] on which the Chebyshev spectral collocation method can be implemented.
August 8, S1793962311000499
2011 16:37 WSPC/262-IJMSSC/S1793-9623
New Algorithm for Solving Nonlinear BVPs in Heat Transfer
363
The Chebyshev spectral collocation method involves approximating the unknown function θ(y) as a truncated series of Chebyshev polynomials of the form θi ≈
N
θi (ξk )Tk (ξj ),
j = 0, 1, . . . , N,
(56)
k=0
where Tk is the kth Chebyshev polynomial and N + 1 is the number of collocation points. Derivatives of the functions θi (y) at the collocation points are represented as N
ds θi s = Dkj θk (ξj ), s dy
(57)
k=0
where s is the order of differentiation and D = 2D, with D being the Chebyshev spectral differentiation matrix whose entries are defined as cj (−1)j+k ck z j − z k zk =− 2(1 − zk2 )
Djk =
j = k; j, k = 0, 1, . . . , N,
Dkk
k = 1, 2, . . . , N − 1,
D00 =
(58)
2N 2 + 1 = −DN N . 6
Applying the Chebyshev spectral method, by substituting Eqs. (56)–(58) in (50) gives Ai−1 Gi = Ri−1 ,
(59)
θi (ξ0 ) = θi (ξN ) = 0,
(60)
subject to
where Ai−1 = a1,i−1 D2 + a2,i−1 D + a3,i−1 , T
Gi = [θi (ξ0 ), θi (ξ1 ), . . . , θi (ξN )] , Ri−1 = [ri−1 (ξ0 ), ri−1 (ξ1 ), . . . , ri−1 (ξN )]T .
(61) (62) (63)
In the above definitions T stands for transpose and ap,i−1 (p = 1, 2, 3) denotes a diagonal matrix of size (N + 1) × (N + 1). The boundary conditions (60) are implemented by deleting the first and last rows and columns of Ai−1 , and deleting the first and last rows of Gi−1 and Ri−1 . The solutions for θi (ξ1 ), θi (ξ2 ), . . . , θi (ξN −1 ) are then obtained from Gi = A−1 i−1 Ri−1 .
(64)
August 8, S1793962311000499
364
2011 16:37 WSPC/262-IJMSSC/S1793-9623
S. S. Motsa
2.2. Example 2: Nonlinear convective–conduction fin problem The problem of one-dimensional heat conduction in convective straight fin with temperature dependent conductivity that was discussed in Refs. 1, 3, 5, 7, 9, 10, 13, 15 and 18, is considered and is given by dθ d (65) (1 + θ) − β 2 θ = 0, dy dy subject to θ (0) = 1,
θ(1) = 1.
(66)
Before applying the SLM, an exact solution for Eq. (65) is found. The transformation is introduced d2 θ du du =u , = dy 2 dy dθ
dθ = u, dy
(67)
and substituted in (65) to obtain (u2 − β 2 θ)dθ + (u + θu)du = 0.
(68)
Multiplying (68) by the integrating factor ρ = 1 + θ gives (u2 − β 2 θ)(1 + θ)dθ + u(1 + θ)2 du = 0, which is an exact differential equation with the following solution: 1 θ 1 2 θ 2 2 2 2 + + u θ + −β θ = c1 , 2 2 2 3 where c1 is an integration constant. From (67) and (70) it is found that 2θ 2 2 c2 + β θ 1 + 3 dθ = , dy 1 + θ
(69)
(70)
(71)
where c2 = 2c1 . The above equation is difficult to integrate, but it can be shown, using long division of the right-hand side of (71), that for the special case when and β are carefully selected in such a way that c2 = −
β2 , 32
(72)
the problem then reduces to (β 2 /32 )(2θ − 1)(1 + θ)2 (2θ − 1) dθ = = . dy 1 + θ 3β 2
(73)
It can be shown that the solution of the above equation is θ(y) =
β3y2 + 3 , 6
whenever =
β2 + 3 . 6
(74)
August 8, S1793962311000499
2011 16:37 WSPC/262-IJMSSC/S1793-9623
New Algorithm for Solving Nonlinear BVPs in Heat Transfer
365
For example, for the special case when = 2 and β = 3, the exact solution of (65) is θ(y) =
3y 2 + 1 . 4
(75)
When the SLM is applied on (65), using the approach described in the previous subsection, one obtains b1,i−1 θi + b2,i−1 θi + b3,i−1 θi = ri−1 ,
(76)
where b1,i−1 = 1 +
i−1
θn ,
b2,i−1 = 2
n=0
ri−1 = −
i−1
θn ,
b3,i−1 =
n=0
θn +
n=0
i−1
i−1 n=0
θn
i−1
θn +
n=0
i−1
θn − β 2 ,
(77)
n=0
i−1
2 θn
i−1
− β2
n=0
θn .
(78)
n=0
One chooses θ0 (y) =
cosh βy , cosh β
(79)
as an initial approximation. Thus, starting from the initial approximation, θ0 , the subsequent solutions for θi , i ≥ 1 are obtained by iteratively solving Eq. (76) subject to the boundary conditions θi (0) = 0,
θi (1) = 0.
(80)
Applying the Chebyshev spectral method in (76) and (80) gives Ai−1 Gi = Ri−1 ,
(81)
subject to N
DN k θi (ξk ) = 0, θi (ξ0 ) = 0,
(82)
k=0
where T
Ri−1 = [ri−1 (ξ0 ), ri−1 (ξ1 ), . . . , ri−1 (ξN )] ,
(83)
and Ai−1 and Gi−1 have the same form as Eqs. (61) and (62). After incorporating the boundary conditions (82) in (81), the solution Gi can be obtained from Gi = A−1 i−1 Ri−1 .
(84)
August 8, S1793962311000499
366
2011 16:37 WSPC/262-IJMSSC/S1793-9623
S. S. Motsa
2.3. Example 3: Nonlinear convective–radiative–conduction fin equation In this section, the method of successive linearization is applied to the nonlinear convective–radiative–conduction problem with temperature-dependent thermal conductivity that was investigated in Ref. 19. The problem extends the convection problem described in Example 2 to include a radiation term (β1 ) and the governing equations are given by d dθ (85) (1 + θ) − β 2 θ − β1 θ4 = 0, dy dy subject to θ (0) = 1,
θ(1) = 1.
(86)
When the SLM is applied on (85), using the approach described in the previous subsection, one obtains b1,i−1 θi + b2,i−1 θi + b3,i−1 θi = ri−1 ,
(87)
where b1,i−1 = 1 +
i−1
θn ,
b2,i−1 = 2
n=0
b3,i−1 =
i−1
θn − β 2 − 4β1
ri−1 = −
i−1
n=0
θn ,
(88)
n=0
i−1
n=0
i−1
3 θn
,
(89)
n=0
θn +
i−1
θn
n=0
i−1 n=0
θn +
i−1 n=0
2 θn
− β2
i−1 n=0
θ n − β1
i−1
4 θn .
n=0
(90) The rest of the analysis used to obtain θ(y) is exactly the same as the analysis given in the previous section. 3. Results and Discussion In this section, the results showing the temperature distribution and rate of heat transfer at the boundaries for the problems discussed in Examples 1 and 2 are presented. To assess the accuracy of the proposed SLM, comparison is made with some of the exact analytical solutions and the numerical solutions obtained using the MATLAB bvp4c routine that is used for solving boundary value problems. In generating the results presented in this study, N + 1 = 60 collocation points were used in the proposed SLM. Table 1 gives a comparison between the successive linearization results and the exact results for θ(y) at selected values of y obtained using (40) for different
August 8, S1793962311000499
2011 16:37 WSPC/262-IJMSSC/S1793-9623
New Algorithm for Solving Nonlinear BVPs in Heat Transfer
367
Table 1. Comparison between the successive linearization results and the exact results for θ(y) at selected values of y.
y
0th order
1st order
2nd order
3rd order
4th order
Exact
0.5
0.2 0.4 0.6 0.8
0.8000000 0.6000000 0.4000000 0.2000000
0.8285714 0.6461538 0.4500000 0.2363636
0.8284271 0.6457513 0.4494898 0.2360680
0.8284271 0.6457513 0.4494897 0.2360680
0.8284271 0.6457513 0.4494897 0.2360680
0.8284271 0.6457513 0.4494897 0.2360680
2
0.2 0.4 0.6 0.8
0.8000000 0.6000000 0.4000000 0.2000000
0.8615385 0.7090909 0.5333333 0.3142857
0.8601478 0.7041695 0.5247312 0.3062657
0.8601471 0.7041595 0.5246951 0.3062258
0.8601471 0.7041595 0.5246951 0.3062258
0.8601471 0.7041595 0.5246951 0.3062258
5
0.2 0.4 0.6 0.8
0.8000000 0.6000000 0.4000000 0.2000000
0.8800000 0.7500000 0.6000000 0.4000000
0.8770370 0.7381579 0.5750000 0.3666667
0.8770330 0.7380832 0.5745968 0.3656863
0.8770330 0.7380832 0.5745967 0.3656854
0.8770330 0.7380832 0.5745967 0.3656854
values of . It is noted from the table that the results from the two methods are in good agreement with accuracy up to 7 decimal places being achieved at third-order of approximation. This shows that the proposed linearization method converges rapidly and is very accurate. The temperature profile for nonlinear conduction problem described in Example 1 is depicted in Fig. 2 for different values of . The figure gives a comparison between the exact solution and the second-order successive linearization solution. It can be
1 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2
5 2 0.5
0.1 0
0
0.2
0.4
0.6
0.8
1
Fig. 2. The comparison of the exact solution (solid line) and second-order successive linearization solution (squares) for θ(y) using = 0.5, = 2, and = 5.
August 8, S1793962311000499
368
2011 16:37 WSPC/262-IJMSSC/S1793-9623
S. S. Motsa
Table 2. Comparison of the successive linearization results for θ (0) against the HAM results of Ref. 11 and the exact solution.
Ref. 11 (HAM)
Present results θ (0)
Exact θ (0)
Order of approximation
θ (0)
Order of approximation
0.5
1 5 10 15 20
−0.75 −0.828125 −0.833496 −0.833328 −0.833333
1 2 3 4 5
−0.83333333 −0.83333333 −0.83333333 −0.83333333 −0.83333333
−0.83333333
2
1 5 10 15 20
−0.75 −0.666992 −0.666668 −0.666667 −0.666667
1 2 3 4 5
−0.66666667 −0.66666667 −0.66666667 −0.66666667 −0.66666667
−0.66666667
5
1 5 10 15 20
−0.75 −0.5876 −0.583377 −0.583333 −0.583333
1 2 3 4 5
−0.58333333 −0.58333333 −0.58333333 −0.58333333 −0.58333333
−0.58333333
seen from the figure that there is good agreement between the two results. This shows the accuracy of the linearization method for different values of . Table 2 gives a comparison of the results for the rate of heat transfer θ (0), for the conduction problem from the current linearization approach versus the HAM of Ref. 11 and the exact solution. It can be seen from the table that the exact solution matches the current linearization results even at first-order approximation. A close approximation to the exact solution that is accurate to 6 decimal places is obtained at the 15th order using the HAM approach of Ref. 11. Consequently, it can be concluded that the current method is efficient and converges more rapidly than the standard HAM. When the solution (53) is differentiated and evaluated at y = 0, one obtains θ (0) = −
1 + (/2) . 1+
(91)
The present results given in Table 2 match the results obtained using Eq. (91). Table 3 gives a comparison between the successive linearization results and the exact results, given by Eq. (75) for θ(y) when β = 3 and = 2 in Example 2. It can be seen from the table that there is very good agreement between the exact results and the results from the proposed method. Accuracy to eight decimal places is achieved at only the third-order of approximation. This shows that the linearization method converges rapidly and is very accurate. The performance of a fin is characterized by fin efficiency. This is defined as the ratio of the fin heat transfer rate to the heat transfer rate of the entire surface if
August 8, S1793962311000499
2011 16:37 WSPC/262-IJMSSC/S1793-9623
New Algorithm for Solving Nonlinear BVPs in Heat Transfer
369
Table 3. Comparison between the successive linearization results and the exact results for θ(y), when β = 3 and = 2. y 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0
0th order
1st order
2nd order
3rd order
4th order
Exact
0.09932793 0.10383131 0.11774980 0.14234550 0.17984866 0.23365997 0.30865887 0.41164604 0.55196004 0.74232415 1.00000000
0.24856397 0.25678080 0.28139952 0.32225212 0.37885407 0.45023267 0.53495164 0.63153972 0.73944053 0.86039951 1.00000000
0.24997862 0.25747719 0.27997690 0.31748753 0.37001668 0.43755717 0.52008247 0.61756716 0.73002322 0.85749729 1.00000000
0.25000000 0.25750000 0.28000000 0.31750000 0.37000000 0.43750000 0.52000000 0.61750000 0.73000000 0.85750000 1.00000000
0.25000000 0.25750000 0.28000000 0.31750000 0.37000000 0.43750000 0.52000000 0.61750000 0.73000000 0.85750000 1.00000000
0.25000000 0.25750000 0.28000000 0.31750000 0.37000000 0.43750000 0.52000000 0.61750000 0.73000000 0.85750000 1.00000000
the entire fin were at the base temperature. Mathematically the fin efficiency ηe is given as (see for example Ref. 15), 1 + dθ ηe = . (92) β 2 dy y=1 Table 4 gives a comparison between the first four terms of the successive linearization approximation and the MATLAB numerical results for the fin efficiency for different values of β and in Example 2. The results from the table illustrate the Table 4. Comparison between the successive linearization results and MATLAB bvp4c for fin efficiency at different values of β and .
β
0th order
1st order
2nd order
3rd order
4th order
bvp4c
0.5
1 2 3 5 10
0.76159416 1.92805516 2.98516426 4.99954602 9.99999996
0.81994146 0.54943770 0.38245270 0.23102934 0.11552930
0.81939437 0.54898744 0.38223305 0.23091059 0.11547006
0.81939431 0.54898742 0.38223304 0.23091058 0.11547005
0.81939431 0.54898742 0.38223304 0.23091058 0.11547005
0.81939431 0.54898742 0.38223304 0.23091058 0.11547005
1.0
1 2 3 5 10
0.76159416 1.92805516 2.98516426 4.99954602 9.99999996
0.85838611 0.60596706 0.42745886 0.25882547 0.12943497
0.85593408 0.60326040 0.42616018 0.25815248 0.12909964
0.85593219 0.60325899 0.42615958 0.25815209 0.12909944
0.85593219 0.60325899 0.42615958 0.25815209 0.12909944
0.85593219 0.60325899 0.42615958 0.25815209 0.12909944
2.0
1 2 3 5 10
0.76159416 1.92805516 2.98516426 4.99954602 9.99999996
0.90613848 0.69774421 0.50674828 0.30859308 0.15434041
0.89834564 0.68437834 0.50001745 0.30540648 0.15275824
0.89831798 0.68431983 0.50000000 0.30539504 0.15275252
0.89831798 0.68431983 0.50000000 0.30539504 0.15275252
0.89831798 0.68431983 0.50000000 0.30539504 0.15275252
5.0
1 2 3 5 10
1.92805516 2.98516426 2.98516426 4.99954602 9.99999996
0.96728423 0.87946712 0.69768102 0.43440941 0.21737025
0.94683071 0.81286612 0.65035214 0.41584670 0.20838577
0.94657973 0.81041783 0.64943418 0.41540397 0.20816658
0.94657973 0.81041783 0.64943418 0.41540397 0.20816658
0.94657973 0.81041783 0.64943476 0.41540397 0.20816658
August 8, S1793962311000499
2011 16:37 WSPC/262-IJMSSC/S1793-9623
S. S. Motsa
370 1
1 0.9
0.95
1
0.8 0.7
5
0.9
0.6
0.85
0.5
2
2
0.4
0.8
0.3 0.2
0.75
0.5 0.7
0
0.2
0.1
0.4
0.6
0.8
1
0
5 0
0.2
0.4
0.6
0.8
1
Fig. 3. Comparison of the second-order successive linearization and bvp4c results when varying (with β = 1) and when varying β (with = 1).
fast convergence of the proposed method. Accuracy to 8 decimal places is achieved at only the third-order approximation. Figure 3 shows the temperature distribution in the fin when varying and β. It is noted from the figure that increasing the conductivity parameter results in an increase in temperature and increasing the convection parameter β results in a decrease in the temperature within the fin. This trend is consistent with the observation made in Refs. 1, 5, 9, 10, 13, 15 and 18 where a similar problem was studied using other methods of solution. The figure also demonstrates that there is excellent agreement between the second-order successive linearization results and the MATLAB bvp4c results. This shows the accuracy and rapid convergence of the linearization method. Table 5 shows the dimensionless temperature at θ(0) for different values of the radiation parameter β1 when β = 1 and = 1. It is noted from the table that the temperature decreases with an increase in radiation. It is only clear from the table that the successive linearization solutions converge to the numerical value. The effect of varying the radiation parameter on the temperature distribution is illustrated graphically in Fig. 4 for β = 1 and = 1 in Example 3. Again, it is noted that there is good agreement between the successive linearization results and the MATLAB numerical results. Table 5.
The dimensionless temperature at y = 0 for different values of β1 with = 1 and β = 1.
β1
0th order
2nd order
3rd order
4th order
6th order
Numerical
1 5 10 20 50 100
0.76159416 0.76159416 0.76159416 0.76159416 0.76159416 0.76159416
0.7063239 0.5718853 0.4986237 0.4329974 0.3824607 0.3699939
0.7063211 0.5718175 0.4968401 0.4206927 0.3336657 0.2942057
0.7063211 0.5718175 0.4968356 0.4204084 0.3263561 0.2668969
0.7063211 0.5718175 0.4968356 0.4204083 0.3262285 0.2643166
0.7063211 0.5718173 0.4968356 0.4204083 0.3262285 0.2643166
August 8, S1793962311000499
2011 16:37 WSPC/262-IJMSSC/S1793-9623
New Algorithm for Solving Nonlinear BVPs in Heat Transfer
371
1 0.9 0.8
1
0.7 0.6
10
0.5 0.4 0.3
100 0.2
0
0.2
0.4
0.6
0.8
1
Fig. 4. The comparison of the numerical solution (solid line) and fourth-order successive linearization solution (squares) for θ(y) using = 1, β = 1 for different values of β1 .
4. Conclusion A new method for solving two-point boundary-value problems has been presented. Using examples of highly nonlinear problems arising in heat transfer, it was clearly illustrated that the proposed SLM is very accurate and converges very rapidly to the exact solutions, which are reported for the first time in this work. Some of the advantages of this method are that solutions are obtained iteratively by integrating linear differential equations. Another desirable feature of the proposed method is that it reduces the volumes of calculations involved because it requires only a few iterations to converge to the exact results. The proposed method is also very simple and straightforward to implement. Because of its efficiency, the proposed method has great potential to be used in solving nonlinear equations in place of traditional numerical methods such as finite differences, shooting methods, and collocation methods. Acknowledgment The authors would like to thank the anonymous reviewers whose comments contributed significantly towards the improvement of this work. References 1. Arslanturk C., A decomposition method for fin efficiency of convective straight fins with temperature-dependent thermal conductivity, Int. Commun. Heat Mass Trans. 32:831–841, 2005.
August 8, S1793962311000499
372
2011 16:37 WSPC/262-IJMSSC/S1793-9623
S. S. Motsa
2. Chang M.-H. A., A decomposition solution for fins with temperature dependent surface heat flux, Int. J. Heat Mass Trans. 48:673–682, 2005. 3. Chiu C. H., Chen C. K., A decomposition method for solving the convective longitudinal fins with variable thermal conductivity, Int. J. Heat Mass Trans. 45:2067–2075, 2002. 4. Lesnic D., Heggs P. J., A decomposition method for power-law fin type problems, Int. Commun. J. Heat Mass Trans. 31:673–682, 2004. 5. Joneidi A. A., Ganji D. D., Babaelahi M., Differential transformation method to determine fin efficiency of convective straight fins with temperature dependent thermal conductivity, Int. Commun. Heat Mass Trans. 36:757–762, 2009. 6. Abbasbandy S., The application of the homotopy analysis method to nonlinear equations arising in heat transfer, Phys. Lett. A 360:109–113, 2006. 7. Domairry G., Fazeli M., Homotopoy analysis method to determine the fin efficiency of convective straight fins with temperature-dependent thermal conductivity, Commun. Nonlinear Sci. Numer. Simulat. 14:489–499, 2009. 8. Domairry G., Bararnia H., An approximation of the analytic solutions of some nonlinear heat transfer equations: A survery using homotopy analysis method, Adv. Studies Theor. Phys. 2:507–518, 2008. 9. Khani F., Ahmadzadeh M., Hamedi-Nezhad S., A series solution of the fin problem with temperature-dependent thermal conductivity, Commun. Nonlinear Sci. Numer. Simulat. 14:3007–3017, 2009. 10. Mustafa I., Application of homotopy analysis method for fin efficiency of convective straight fins with temperature-dependent thermal conductivity, Math. Comput. Simulat. 79:189–200, 2008. 11. Sajid M., Hayat T., Comparison of HAM and HPM methods in nonlinear heat conduction and convection equations, Nonlinear Anal. Real World Appl. 9:2296–2301, 2008. 12. Hosseini M. J., Gorji M., Ghanbarpour M., Solution of temperature distribution in a radiating fin using homotopy perturbatiod method, Math. Problem Eng. DOI:10.1155/2009/831362, 2009. 13. Rajabi A., Homotopy perturbation method for fin efficiency of convective straight fins with temperature-dependent thermal conductivity, Phys. Lett. A 364:33–37, 2007. 14. Rajabi A., Ganji D. D., Taherian H., Application of homotopy perturbation method in nonlinear heat conduction and convection equations, Phys. Lett. A 360:570–573, 2007. 15. Kim S., Huang C. H., A series solution of the fin problem with a temperaturedependent thermal conductivity, J. Phys. D: Appl. Phys. 39:4894–4901, 2006. 16. Kim S., Huang C. H., A series solution of the nonlinear fin problem with a temperature-dependent thermal conductivity and heat transfer coeffient, J. Phys. D: Appl. Phys. 40:2979–2987, 2006. 17. Ganji D. D., Hosseini M. J., Shayegh J., Some nonlinear heat transfer equations solved by three approximate methods, Int. Commun. Heat Mass Trans. 34:1003–1016, 2007. 18. Ganji D. D., Afrouzi G. A., Talarposhti F., Application of variational iteration method and homotopy-perturbation method for nonlinear heat diffusion and heat transfer equations, Phys. Lett. A 368:450–457, 2007. 19. Miansari Mo., Ganji D. D., Miansari Me., Application of He’s variational iteration method to nonlinear heat transfer equations, Phys. Lett. A 372:779–785, 2008. 20. Tari H., Ganji D. D., Babazadeh H., The application of He’s variational iteration method to nonlinear equations arising in heat transfer, Phys. Lett. A 363:213–217, 2007.
August 8, S1793962311000499
2011 16:37 WSPC/262-IJMSSC/S1793-9623
New Algorithm for Solving Nonlinear BVPs in Heat Transfer
373
21. Kierzenka J., Shampine L. F., A BVP solver based on residual control and the MATLAB pse, ACM Trans. Math. Software 27(3):299–316, 2001. 22. Canuto C., Hussaini M. Y., Quarteroni A., Zang T. A., Spectral Methods in Fluid Dynamics, Springer-Verlag, Berlin, 1988. 23. Don W. S., Solomonoff A., Accuracy and speed in computing the Chebyshev collocation derivative SIAM, J. Sci. Comput. 16:1253–1268, 1995. 24. Trefethen L. N., Spectral Methods in MATLAB, SIAM, 2000.
