# Successive Linearization Analysis of the Effects of Partial Slip, Thermal Diffusion, and Diffusion-Thermo on Steady MHD Convective Flow due to a Rotating Disk

Authors: S. S. Motsa; S. Shateyi
Journal: Mathematical Problems in Engineering
 (2012)

---

## Abstract

We proposed a general formulation of the successive linearization method for solving highly nonlinear boundary value problem arising in rotating disk flow. The problem was studied under the effects of partial slip, thermal diffusion, and diffusion-thermo. The governing fundamental conservation equations of mass, momentum, angular momentum, energy, and concentration are transformed into a system of ordinary differential equations by means of similarity transformations. A parametric study illustrating the influence of the magnetic field, slip factor, Eckert number, Dufour and Soret numbers was carried out.


---

## Body

## 1. Introduction

The growing concern of global warming effects of convectional fossil fuels forces energy research into the direction of the renewable energy sources, like in the design of turbines and turbomachines. The energy, chemical, and automobile industries extensively use disk-shaped bodies which are often encountered in many engineering applications and heat transfer problem of natural convection boundary layer flow over a rotating disk, which occurs in rotating heat exchangers, rotating disks for biofuels production and turbines.In his pioneering work, von Karman [1] considered the case of an infinite disk and gave a formulation of the hydrodynamic problem. Since then, a lot of work has been done in this field of study (see, e.g., Cochran [2]; Benton [3]; Attia [4]; Miklavčič and Wang [5]; Arikoglu and Ozkol [6], among others).The relations between the fluxes and driving potentials are more intricate in nature, when heat and mass transfer occurs simultaneously in a moving fluid. It is now a known fact that an energy flux can be generated not only by temperature gradients but also by composition gradients. The energy flux caused by composition gradient is called the diffusion thermal or Dufour effect. Temperature gradients can also create mass fluxes. This phenomenon is called the thermal diffusion or Soret effect.Generally, Soret and Dufour effects are of smaller-order magnitude than the effects prescribed by Fourier's or Fick's laws and are often neglected in heat and mass transfer processes. However, there are exceptions. The Soret effect plays an important role in the operation of solar ponds, biological systems, and the microstructure of the world oceans. On the other hand Dufour effect was found to be of order of considerable magnitude in mixture between gasses with very low molecular weight (He, H2) and of medium molecular weight (N2, air), (Kafoussias and Williams [7]).Siginer [8] observed that the Soret effect plays an important role in the concentration distribution of different components in hydrocarbon mixtures as these are mainly driven by phase separation and diffusion.Osalusi et al. [9] investigated Soret and Dufour effects on combined heat and mass transfer of steady hydromagnetic convective and slip flow due to a rotating disk in the presence of viscous dissipation and ohmic heating. Rashidi et al. [10] solved the same problem analytically using the Homotopy analysis method. Turkyilmazoglu [11] derived analytical expressions for the solution of steady, laminar, incompressible, viscous fluid of the boundary layer flow due to a rotating disk in the presence of a uniform suction/injection. Shateyi et al. [12] investigated the influence of a magnetic field on heat and mass transfer by mixed convection from vertical surfaces in the presence of Hall, radiation, Soret and Dufour effects. Alam and Ahammad [13] investigated the effects of variable chemical reaction and variable electric conductivity on free convection flow with heat and mass transfer over an inclined permeable stretching sheet under the influence of Dufour and Soret effects with variable heat and mass fluxes.Srinivasacharya and Kaladhar [14] employed the Homotopy analysis method to investigate the Hall and ion-slip effects on steady free convective heat transfer flow between two cylinders in couple stress fluid flow. Pal and Mondal [15] carried out an analysis of the combined effects of Soret and Dufour on unsteady MHD non-Darcy mixed convection over a stretching sheet embedded in a saturated porous medium in the presence of thermal radiation, viscous dissipation, and first-order chemical reaction.More improved models which more accurately predict the concentrations of the different components of the fluid physics in crude oil reservoirs are necessary. Through the advancement of the symbolic computation software such as MATHEMATICA, MAPLE, and MATLAB, approximate analytic methods for nonlinear problems have been adopted by many researchers. Among these are homotopy perturbation method, Dehghan et al. [16], homotopy analysis method, Liao [17], Dehghan and Salehi [18], spectral homotopy analysis method, Motsa et al. [19], improved spectral homotopy analysis method, Motsa et al. [20], and successive linearization method (SLM), Motsa and Shateyi [21].In this work we present a general formulation of the successive linearization method that can be used to solve any two-point nonlinear boundary value problem. We test the validity of the method of solution on the model of steady MHD convective flow due to a rotating disk with partial slip in the presence of thermal diffusion and diffusion-thermo effects. The objective of this study is to demonstrate the ease of use of the SLM approach and its accuracy when solving nonlinear BVPs arising in rotating disk flow.


## 2. Problem Formulation

The problem of steady hydromagnetic convective and slip flow due to a rotating disk in the presence of viscous dissipation and ohmic heating is considered in this study. The effects of thermal diffusion, diffusion-thermo, and heat and mass transfer are also considered in the present study. The disk rotates atz=0 with constant angular velocity Ω. We have z being the vertical axis in the cylindrical coordinates system with r and ϕ as the radial and tangential axes, respectively. The surface of the rotating disk is maintained at a uniform temperature Tw and uniform concentration Cw. Far away from the surface, the free stream is kept at a constant temperature T∞, concentration C∞, and pressure P∞. The viscous fluid is an electrically conducting one. A uniform magnetic field is applied normally to the surface of the disk and has a constant magnetic flux density B0 which is assumed unchanged by taking small magnetic Reynolds number. In this study we also assume that there is no applied electric field (Rashidi et al. [10]) and that the Hall effects are negligible. The equations governing the motion of MHD laminar flow take the following form:
(2.1)∂u∂r+ur+∂w∂z=0,(2.2)u∂u∂r-v2r+w∂u∂z+1ρ∂P∂r=ν(∂2u∂r2+1r∂u∂r-ur2+∂2u∂z2)-σB02uρ,(2.3)u∂v∂r+uvr+w∂v∂z=ν(∂2v∂r2+1r∂v∂r-vr2+∂2v∂z2)-σB02vρ,(2.4)u∂w∂r+w∂w∂z+1ρ∂P∂z=ν(∂2w∂r2+1r∂w∂r+∂2w∂z2),(2.5)u∂T∂r+w∂T∂z=&κρcp(∂2T∂r2+1r∂T∂r+∂2T∂z2)+σB02ρcp(u2+v2)+μρcp[(∂u∂z)2+(∂v∂z)2]+DkTcscp(∂2C∂r2+1r∂C∂r+∂2C∂z2),(2.6)u∂C∂r+w∂C∂z=D(∂2C∂r2+1r∂C∂r+∂2C∂z2)+DkTTm(∂2T∂r2+1r∂T∂r+∂2T∂z2),
where u,v,w are the velocity components of the flow in the directions of increasing (r,ϕ,z), respectively, P is the pressure, ρ is the fluid density. T and C are the fluid temperature and concentration, respectively. ν=μ/ρ is the kinematic viscosity of the ambient fluid, κ is the thermal conductivity, σ is the electrical conductivity, cp is the specific heat at constant pressure, D is the molecular diffusion coefficient, kT is the thermal-diffusion ratio, cs is the concentration susceptibility, and Tm is the mean fluid temperature.The boundary conditions for the model can be expressed as(2.7)u=2-ξξλ∂u∂z,v=rΩ+2-ξξλ∂v∂z,w=0,T=Tw,C=Cwatz=0,u⟶0,v⟶0,T⟶T∞,C⟶C∞,P⟶P∞asz⟶∞.Introducing a dimensionless normal distance from the disk,η=zν-1/2 along with the following von Karman transformations, transforms the governing equations into a system of nonlinear ordinary differential equations:
(2.8)u=rf(η),v=Rg(η),w=ν1/2h(η),T=θ(η),C=ϑ,
where f,g,h,θ, and ϑ are nondimensional functions of modified dimensionless vertical coordinate η. The dimensionless temperature and concentration are θ=(T-T∞)/(Tw-T∞) and ϑ=(C-C∞)/(Cw-C∞).We now have the following nonlinear ordinary differential equations:(2.9)f′′-(f)2+g2-f′h-Mf=0,(2.10)g′′-2fg-hg′-Mg=0,(2.11)h′+2f=0,(2.12)1Prθ′′-hθ′+MEc(f2+g2)+Ec(f′2+g′2)+Duϑ′′=0,(2.13)1Scϑ′′-hϑ′+Srθ′′=0,
where M=σB02/ρΩ is the magnetic parameter, Pr=ρνcp/κ is the Prandtl number, Ec=L2Ω2R2/cp(Tw-T∞) is the Eckert number, Du=D(Cw-C∞)kT/cscpν(Tw-T∞) is the Dufour number, Sc=ν/D is the Schmidt number, Sr=D(Tw-T∞)kT/ν(Cw-C∞) is the Soret number. The continuity equation (2.1) is transformed into (2.11). Also the pressure term can be obtained from (2.4) as P/ρ=ν(∂w/∂z)-(w2/2)+constant, and this expression is then used to eliminate all the pressure terms [22]. The prime denotes derivative with respect to η. The transformed boundary conditions now become
(2.14)f(0)=γf′(0),g(0)=1+γg′(0),h(0)=0,θ(0)=1,ϑ(0)=1,(2.15)f(∞)=0,g(∞)=0,θ(∞)=0,ϑ(0)=0,
with γ=[(2-ξ)λ]/ξν1/2 being the slip factor. The boundary conditions show that the radial component f, the tangential component g, temperature and concentration vanish sufficiently away from the rotating disk, whereas the axial velocity component h is anticipated to approach an unknown asymptotic limit for large η-values.


## 3. Numerical Method of Solution

In this section we present the numerical solution of the governing nonlinear systems of (2.9)–(2.13). The equation system will be solved using the successive linearization method (SLM). The SLM was recently introduced as an efficient and robust method for solving boundary value problems [23]. In its basic application the SLM seeks to linearize the governing nonlinear differential equations to form an iterative system of linear differential equations which, in most cases, cannot be solved analytically. The Chebyshev pseudospectral method (or any other collocation method or numerical scheme) is then used to transform the iterative sequence of linearized differential equations into a system of linear algebraic equations which are converted into a matrix system. In this paper we present a compact general version of the SLM which can be applied to any system of nonlinear boundary value problems. We observe that (2.9)–(2.11) can be solved independently of (2.12)-(2.13). To generate a more general solution method, we introduce the following notation:
(3.1)y1=f,y2=g,y3=h.Using the definition (3.1), the equation system (2.9)–(2.11) can be separated into its constituent linear and nonlinear components as
(3.2)arryr′′+∑k=13‍brkyk′+∑k=13‍crkyk+Fr(y1,y2,y3)=0
subject to
(3.3)y1(0)=γy1′(0),y2(0)=1+γy2′(0),y3(0)=0,(3.4)y1(∞)=0,y2(∞)=0,
where the coefficients arr, brk, crk (for r,k=1,2,3) represent the constant factors of the various derivatives and Fr(y1,y2,y3) denotes the nonlinear components of (2.9)–(2.11) for r=1,2,3.The SLM approach assumes that the solution of the system (3.2) can be expressed as
(3.5)yr(η)=yr,s(η)+∑m=0s-1‍yr,m(η).Starting from a suitable initial approximationyr,0(t), the solution for yr,s(t) can be obtained by successively linearising equation (3.2) and solving the resulting linear system. The general form of the linearized equations for the SLM algorithm corresponding to (3.2) is
(3.6)arryr′′+∑k=13‍brkyk′+∑k=13‍crkyk+Fr(y1,y2,y3)=0,(3.7)arryr,s′′+∑k=13‍(br,k+∂Fr*∂yk′)yk,s′+∑k=13‍(cr,k+∂Fr*∂yk)yk,s=Gr,s-1,
subject to
(3.8)γy1,s′(0)-y1,s(0)=-(γ∑m=0s-1‍y1,m′(0)-∑m=0s-1‍y1,m(0)),(3.9)γy2,s′(0)-y2,s(0)=-(γ∑m=0s-1‍y2,m′(0)-∑m=0s-1‍y2,m(0)+1),(3.10)y3,s(0)=0,y1,s(∞)=0,y2,s(∞)=0,
where
(3.11)∂Fr*∂yk=∂Fr∂yk(y-1,y-2,y-3),withy-r=∑m=0s-1‍yr,m,(3.12)Gr,s-1=-[arry-r′′+∑k=13‍brky-k′+∑k=13‍crky-k+Fr(y-1,y-2,y-3)].The SLM algorithm (3.7) requires an initial approximation yr,0. This initial guess is chosen as a function that satisfies the boundary conditions (3.3). Suitable functions are
(3.13)y1,0=0,y2,0=11+γe-η,y3,0=0.Note that these initial guesses are chosen in such a way that they satisfy the boundary conditions. Thus, starting from the initial approximation (3.13), the solution of the linearized equation system (3.7) is solved iteratively using spectral collocation methods (or any other numerical method) for yr,s(η) (for s=1,2,3,…). The approximate solution for each yr(η) is determined as the series solution
(3.14)yr(η)=yr,0(η)+yr,1(η)+yr,2(η)+⋯An SLM solution is said to be of orderi if the above series is truncated at s=i, that is, if
(3.15)yr(η)=∑m=0i‍yr,m(η).To solve the linearised system (3.7) we use the Chebyshev collocation spectral method in which the solution space is discretized using the Chebyshev-Gauss-Lobatto collocation point:
(3.16)ηj=cos⁡(πjN),j=0,1,…,N,
which are the extrema of the Nth order Chebyshev polynomial
(3.17)TN(η)=cos⁡(Ncos⁡-1η).Before applying the spectral method, it is convenient to transform the governing physical region[0,∞] for the problem to the interval [−1,1] on which the spectral method is defined. This can be achieved by using the linear transformation τ=L(η+1)/2. Here, L is chosen to be sufficiently large enough to numerically approximate infinity. The Chebyshev spectral collocation method (see, e.g., [24–26]) is based on the idea of introducing a differentiation matrix D which is used to approximate the derivatives of the unknown variables yr,s at the collocation points as the matrix vector product:
(3.18)dyr,sdη=∑k=0N‍Djkyr,s(τk)=DYr,s,j=0,1,…,N,
where D=2D/L and Yr,s=[yr,s(τ0),yr,s(τ1),…,yr,s(τN)]T is the vector function at the collocation points τj. The entries of D can be computed in different ways. In this work we use the method proposed by Trefethen [25] in the cheb.m MATLAB m-file. If we denote the entries of the derivative matrix D by Djk, we can apply the spectral collocation method, with derivative matrices on the linear boundary value system (3.7) and boundary conditions (3.8)–(3.10); we obtain the following linear matrix system:
(3.19)As-1Ys=Gs-1,
with the boundary conditions
(3.20)γ∑k=0N‍DNky1,s(τk)-y1,s(τN)=ps-1,(3.21)γ∑k=0N‍DNky2,s(τk)-y2,s(τN)=qs-1,(3.22)y3,s(τN)=0,y1,s(τ0)=0,y2,s(τ0)=0,
where ps-1 and qs-1 denote the right hand side of (3.8) and (3.9), respectively, when evaluated at the collocation points τj(j=0,2,…,N). The vectors Ys and Gs-1 have dimension 3(N+1)×1 and are defined as
(3.23)Ys=[Y1,s;Y2,s;Y3,s]T,Gs-1=[G1,s-1;G2,s-1;G3,s-1]T,
with Yr,s and Gr,s-1 being the vectors yr,s1 and Gr,s-1(r=1,2,3), respectively, evaluated at the collocation points. The matrix As-1 is an 3(N+1)×3(N+1) matrix that is defined as
(3.24)As-1=[A11A12A13A21A22A23A31A32A33],
with
(3.25)Apq={appD2+(bppI+∂Fp*∂yp′)D+cppI+∂Fp*∂yp,p=q,(bpqI+∂Fp*∂yp′)D+cpqI+∂Fp*∂yqI,p≠q,
where I is an identity matrix of order N+1. After imposing the boundary conditions (3.20)–(3.22) on the matrix system (3.19), the solutions for yr,s(η)(s=1,2,3…) can be obtained by solving the iterative matrix system:
(3.26)Ys=As-1-1Gs-1.The solutions forf(η), g(η), and h(η) obtained from (3.26) are substituted in (2.12) and (2.13) which now become a linear coupled system for the unknown variables θ and ϑ. Applying the Chebyshev spectral collocation method on the resulting linear system, we obtain the following matrix system:
(3.27)[1PrD2-hDDuDSrD21ScD2-hD][θ~ϑ~]=[KO],
subject to the boundary conditions
(3.28)θ(τN)=ϑ(τN)=1,θ(τ0)=ϑ(τ0)=0,
where
(3.29)K=-MEc(f2+g2)-Ec(f′2+g′2)
and O is a (N+1)×1 vector of zeros, f,g,h are diagonal matrices corresponding to the solutions f,g, and h, respectively, when evaluated at the collocation points τj(j=0,1,2,…,N), and θ~, ϑ~ are the approximations of θ(η) and ϑ(η) at the collocation points. After the boundary conditions (3.28) have been imposed on the matrix system (3.27), the resulting system can easily be inverted and solved as
(3.30)[θ~ϑ~]=[1PrD2-hDDuDSrD21ScD2-hD]-1[K1K2].


## 4. Numerical Results

The semianalytical results are obtained by solving (2.10)–(2.15) using the method elucidated in the previous section for various values of physical parameters to describe the physics of the problem. We remark that, unless otherwise specified, the SLM results presented in this analysis were obtained using N=200 collocation points, and L=40 was used as a numerical approximation infinity. To test the accuracy of the results we consider the residuals which are obtained by substituting the SLM approximate solutions and checking in the governing equations. The residuals corresponding to f and g, for instance, are given by
(4.1)Res(f);f′′-(f)2+g2-f′h-Mf=0,Res(g);g′′-2fg-hg′-Mg=0.In Figure1 we give the absolute value for the residual functions of f and g, respectively, for different iterations. We observe that the residual becomes increasingly smaller as the number of iterations increases. We remark that the SLM algorithm was implemented using MATLAB which has a machine epsilon of 2.2204e-16 which means that numbers are stored with about 15-16 digits of precision. We observe that after only 5 iterations, the residual error has almost converged to the machine epsilon of the computational software. This clearly shows the accuracy of the proposed method. We remark then that residual graphs of the same problem obtained using the homotopy analysis method (HAM) in [10] produced residuals values whose magnitude is significantly larger than the magnitude of the residuals depicted in Figure 1. This is the case even when the (HAM) iterations are increased to order 30. This indicates that the proposed SLM method is much more accurate than the HAM for the type of problems discussed in this paper. In Figure 2 we show the norm of the residual against the iterations. Again, it can be seen from this graph that as the iterations increase, the norm of the residual becomes smaller.Resf and Res g for different order of approximation, when γ=0.2,M=0.2,Pr=0.71,Sc=0.2,Sr=1,Du=0.06,Ec=0.2.
(a)
(b)Norm of the residuals against iterations, whenγ=0.2,M=0.2,Pr=0.71,Sc=0.2,Sr=1,Du=0.06,Ec=0.2.
(a)
(b)In Figures3–5 we show the graphical representations of the numerical results when varying some of the governing parameters.Sixth order SLM results for the variation ofγ on the flow properties when M=0.2,Pr=0.71,Sc=0.2,Sr=1,Du=0.06,Ec=0.2.
(a)
(b)
(c)
(d)
(e)Sixth order SLM results for the variation ofM on the flow properties when γ=0.2,Pr=0.71,Sc=0.2,Sr=1,Du=0.06,Ec=0.2.
(a)
(b)
(c)
(d)
(e)Effects ofEc and Du on temperature and also effects of Sc on concentration.
(a)
(b)
(c)The effect of the slip factorγ on the radial velocity f(η), tangential velocity, g(η), axial velocity h(η), temperature θ(η), and concentration ϑ(η) profiles is shown in Figure 3. In this investigation, we took Prandtt number of air (Pr=0.71) and Schmidt to be 0.2. We observe in this figure that the radial velocity is reduced as the notating surface becomes more slippery. Also the shear-driven flow (tangential velocity) is reduced by the increasing values of the slip factor. However, towards the rotation disk, the axial velocity is increased as the slip factor increases. The centrifugal force associated with this circular motion causes the reductions of both the radial flow f(η) and the shear-driven flow g(η). In Figure 3, we also observe that the radial outflow f(η) is compensated by an axial inflow h(η) towards the rotating disk as γ increases, in accordance with (2.11). We also observe that the magnitude of the temperature and concentration profiles slightly increase with an increase in a slip parameter.Figure4 depicts the effect of the magnetic field parameter M on the radial velocity f(η) profiles, tangential velocity profiles g(η) profiles, axial velocity h(η) profiles, temperature θ(η) profiles, and concentration ϑ(η).It is observed in this figure that both the radial and tangential velocity profiles decrease with the increase of the magnetic field parameter values. This is because the presence of a magnetic field in an electrically conducting fluid introduces a force called Lorentz force which creates a drag-like force to slow down the flow along the disk surface. However, the Lorentz force accelerates flow in the axial direction causing the axial velocity profiles to increase asM increases, (Alam and Ahammad [13], Motsa and Sibanda [27], and Shateyi et al. [28]). The temperature and concentration profiles increase as M increases. This is because the drag-like forces caused by the pressure of a magnetic field cause both the thermal and solutal boundary layer to increase, hence increasing the temperature and concentration of the slowed down flow.Figure5 shows the effect of the Eckert number Ec on the temperature profiles and the effect on the Soret number Sr on the concentration profiles. The effect of increasing the values of the Eckert number is to enhance the temperature of the fluid at any point. This is expected according to (2.12). This is also because the heat energy is stored in the liquid due to frictional heating. The increased values of Ec lead to a strong viscous dissipation which significantly increases the temperature profiles, (Pal and Mondal [29]). The temperature of the fluid is also increased as the values of the Dufour parameter increase. Pal and Mondal [30] observed that the effect of increasing Du is to reduce the Nusselt number thereby causing the temperature profiles to increase as depicted in Figure 5 of the current investigation. We observe also in this figure that the concentration profiles increase with increasing Soret number values. This is because the increase of the values of Sr is caused by the increase of the temperature gradient which in turn leads to the increase of the concentration distributions with the fluid flow.


## 5. Conclusion

Numerical analysis has been carried out to investigate the MHD convective flow due to a rotating disk in the presence of viscous dissipation and ohmic heating, thermal-diffusion, and diffusion-thermo effects. The partial differential equations which describe the problem are transformed by using a suitable similarity transformation. The resultant nonlinear ordinary differential equations are then solved numerically using the SLM with the Chebyshev spectral collocation method. It was observed that the slip factor has significant effects on the fluid properties. The magnetic field helps to decelerate flow in the radial and tangential directions but accelerates flow in the axial direction. It also increases the temperature and concentration distributions. Both thermal diffusion and diffusion-thermo have some significant effect on the concentration and temperature profiles.



---

*Source: 397637-2012-06-27.xml*