# The Effects of Chemical Reaction,  Hall, and Ion-Slip Currents on MHD Micropolar Fluid Flow with Thermal Diffusivity Using a Novel Numerical Technique

Authors: S. S. Motsa; S. Shateyi
Journal: Journal of Applied Mathematics
 (2012)

---

## Abstract

The problem of magnetomicropolar fluid flow, heat, and mass transfer with suction through a porous medium is numerically analyzed. The problem was studied under the effects of chemical reaction, Hall, ion-slip currents, and variable thermal diffusivity. The governing fundamental conservation equations of mass, momentum, angular momentum, energy, and concentration are converted into a system of nonlinear ordinary differential equations by means of similarity transformation. The resulting system of coupled nonlinear ordinary differential equations is the then solved using a fairly new technique known as the successive linearization method together with the Chebyshev collocation method. A parametric study illustrating the influence of the magnetic strength, Hall and ion-slip currents, Eckert number, chemical reaction and
permeability on the Nusselt and Sherwood numbers, skin friction coefficients, velocities, temperature, and concentration was carried out.


---

## Body

## 1. Introduction

Eringen [1] proposed the theory of micropolar fluids, which shows microrotation effects as well as microinertia,  as these flow properties cannot be described by the classical Navier-Stokes theory. Since the pioneering work by Eringen,  the theory of micropolar fluid has generated a lot of interest. Extension has been done, to include studies of magneto-micropolar fluid with Hall current and ion-slip currents with heat transfer due to vast possible engineering applications in areas like power generators, MHD accelerators, refrigeration coils, electric transformers, and heating elements. MHD flows of a viscous and incompressible fluid have been extensively studied with the effect of Hall current by Chamkha [2], Seddeek [3], Takhar et al. [4], Shateyi et al. [5, 6], Salem and Abd El-aziz [7], among others.The momentum, heat, and mass transport on stretching sheet has several applications in polymer processing as well as in electrochemistry. The heat transfer problem associated with the boundary layer micropolar fluid under different physical conditions has been studied by several authors. Takhar et al. [8] considered diffusion of a chemical reactive species from a stretching sheet. Muthucumaraswamy and Ganesan [9] studied diffusion and first-order chemical reaction on impulsively started infinite vertical plate with variable temperature. Shateyi [10] investigated thermal radiation and buoyancy effects on heat and mass transfer over a semi-infinite stretching sheet with suction and blowing. Shateyi and Motsa [11] numerically investigated the unsteady heat, mass, and fluid transfer over a horizontal stretching sheet.The porous media heat and mass transfer problems have several practical engineering applications such as geothermal systems, crude oil extraction, and ground-water pollution. The study of chemical reaction heat transfer in porous medium has important applications such as in tabular reactors, oxidation of solid materials, and synthesis of ceramic materials. Elgazery [12] numerically analyzed the effects of chemical reaction, Hall, ion-slip currents, variable viscosity and variable thermal diffusivity on the problem of magnetomicropolar fluid flow, heat and mass transfer with suction and blowing through a porous medium. Srinivasacharya and RamReddy [13] presented natural convection heat and mass transfer along a vertical plate embedded in a doubly stratified micropolar fluid saturated non-Darcy porous medium. Ishak et al. [14] investigated the unsteady boundary layer flow over a stretching permeable surface.Elshehawey et al. [15] applied the Chebyshev finite-difference method to investigate the effects of Hall and ion-slip currents on magneto-hydrodynamic flow with variable thermal conductivity. Seddeek and Salama [16]. studied effects of chemical reaction and variable viscosity on hydromagnetic mixed convection heat and mass transfer for Hiemenz flow through porous media with radiation. Shateyi and Motsa [17] analyzed numerically the problem of magnetohydrodynamic flow and heat transfer of a viscous, incompressible and electrically conducting fluid past a semi-infinite unsteady stretching sheet with Hall currents, variable viscosity and thermal diffusivity. Mahmoud and Waheed [18] performed a theoretical analysis to study heat transfer characteristics of magnetohydrodynamic mixed convection flow of a micropolar fluid past a stretching surface with slip velocity at the surface and heat generation/absorption.The aim of this work is to analyze the effects of chemical reaction, Hall and ion-slip currents on the MHD flow of a micropolar fluid through a porous medium using the successive linearization method (SLM). The SLM is based on a novel idea of iteratively linearizing the underlying governing nonlinear equations, which are written in similarity form, and solving the resulting equations using spectral methods. This approach has been successfully applied to different flow problem (see, e.g., Makukula et al. [19–22]; Shateyi and Motsa [17]; Motsa and Shateyi [23, 24]). The influences of the governing parameters on the flow characteristics are illustrated graphically and using tables.


## 2. Problem Formulation

We consider a steady, incompressible, magneto-micropolar and electrically conducting viscous fluid flowing over a horizontal plate in thex-direction under the influence of mass transfer with chemical reaction through a porous medium. We let the reaction of a species, say A with B be the first order homogeneous chemical reaction with a constant rate, κ. We also assume that the concentration of dissolved A is small enough and physical properties ρ and D are virtually constant throughout the fluid. The flow under consideration is also subjected to a strong transverse magnetic field B0 with a constant intensity along the y-axis (see Figure 1).Figure 1
The physical configuration and coordinate system for the boundary layer flow.Generally, an electrically conducting fluid is affected by Hall and ion-slip currents in the presence of a magnetic field. The effect of Hall current gives rise to a force in thez-direction, which induces a cross-flow in the z-direction and hence the flow becomes three-dimensional.Following Elgazery [12], the generalized Ohms law including Hall currents is given by(2.1)J=σ1+(ω/νe)2(E+(V+B)-1ene),where σ is the electrical conductivity, J=(Jx,Jy,Jz) is the electric current vector, V is the velocity vector,   E is the intensity of the electric field, B is the magnetic induction vector, 1/ene is the Hall factor, ω/νe is the Hall parameter, ne is the density number of the electrons, e is the charge of electron, νe is the electron-atom collision frequency and ω exceeds the collision frequency. The electric magnetic field can force both the ions and electrons to produce a relative drift between them and neutral particles, when the ratio, ω/νe becomes very large. The drift is called “ion-slip” and is however negligible for highly ionized gases.We assume that all the fluid properties are isotropic and constant. The velocity componentu on a stretching sheet is proportional to its distance from the slit. Under the usual boundary layer and Boussinesq approximations, the governing equation in (x,y,z)-coordinate for the problem under consideration can be written as follows:(2.2)∂u∂x+∂v∂y=0,(2.3)u∂u∂x+v∂u∂y=ν∂2u∂y2+Kρ∂N∂y-σB02ρ(αe2+βe2)(αeu+βew)-μρk*u,(2.4)u∂w∂x+v∂w∂y=ν∂2w∂y2+σB02ρ(αe2+βe2)(βeu-αew)-μρk*w,(2.5)G1K∂2N∂y2-2N-∂u∂y=0,(2.6)u∂T∂x+v∂T∂y=∂∂y(α∂T∂y)+μρcp[(∂u∂y)2+(∂w∂y)2]+1ρcpσB02(αe2+βe2)(u2+w2),(2.7)u∂C∂x+v∂C∂y=D∂2C∂y2-k(C-C∞).
The boundary conditions are
(2.8)u(x,0)=Us=bx,v(x,0)=-Vw,w(x,0)=0,T(x,0)=Tw,C(x,0)=Cw,u(x,∞)=0,w(x,∞)=0,T(x,∞)=T∞,C(x,∞)=C∞.
Here (u,v,w) are the fluid velocity components in the (x,y,z)-directions, respectively. ρ is the fluid density, μ is the dynamic viscosity, cp the specific heat at constant pressure, βt the coefficient of thermal expansion, βc the coefficient of concentration expansion, U∞, T∞ and C∞ the fluid free stream velocity, temperature, and concentration, respectively, N is the component of microrotation of the fluid, K is the vortex viscosity, G1 is the spin gradient viscosity. and k* permeability of the porous medium. βe and βi are respectively, the Hall parameter and the ion-slip parameter, where αe=1+βiβe. b is a constant with dimension (time)-1. Us is the surface velocity and Vw is suction/injection velocity, ν is the kinematic viscosity, D is the chemical molecular diffusivity, α is the thermal diffusivity, k is the chemical reaction, and B0 is the magnetic field of constant strength.The governing differential equations (2.2)–(2.7) together with the boundary conditions (2.8) are non-dimensionalised using the following similarity transformations:
(2.9)η=(bν)1/2y,u=bxf′(η),v=-bνf(η),w=bνg(η),N=b3νxh(η),θ=T-T∞Tw-T∞,ϕ=C-C∞Cw-C∞,
where f(η),g(η),h(η),θ and ϕ(η) are the dimensionless stream function and microrotation functions, temperature, and concentration distributions functions, respectively. η is the similarity variable. The variation of thermal diffusivity with the dimensionless temperature is written in the form:
(2.10)α=α0(1+β2θ),
where β2 is parameter depending on the nature of the fluid, and α0 is the thermal diffusivity of the fluid at temperature Tw. Then introducing the relations (2.9) into the equations (2.3)–(2.7) respectively, we obtain the following local similarity equations:
(2.11)f′′′+ff′′-(f′)2+N1h′-Mαe2+βe2(αef′+βeReg)-1kpf′=0,(2.12)g′′+fg′+Mαe2+βe2(βeRef′-αeg)-1kpg=0,(2.13)Gh′′-2h-f′′=0,(2.14)θ′′+β2(1+β2θ)(θ′)2+PrEc(1+β2θ)[(f′′)2+(g′)2Re+Mαe2+βe2((f′)2+g2Re)]+Pr(1+β2θ)fθ′=0,(2.15)1Scϕ′′+fϕ′-γReϕ=0,
where the primes denote differentiation with respect to η.M=σB02/ρb is the magnetic parameter, N1=K/ρν is the coupling constant parameter, G=G1b/Kν is the microrotation parameter, Pr=ν/α is the Prandtl number, Sc=ν/D is the Schmidt number, Ec=Us2/cp(Tw-T∞) is the Eckert number, kp=k*b/μ0 is the permeability parameter, Re=xUs/ν0, is the Reynolds number, γ=κν/Us2 is the nondimensional chemical reaction parameter. In view of the similarity transformations, the boundary conditions transform into:
(2.16)f(0)=fw,f′(0)=1,g(0)=0,h(0)=0,θ(0)=1,ϕ(0)=1,f′(∞)=0,g(∞)=0,h(∞)=0,T(∞)=0,C(∞)=0,
where fw=Vw/bν is the mass transfer coefficient such that fw>0 indicates suction and fw<0 indicates blowing at the surface.For practical applications, quantities of interest include the velocity componentsu and w, temperature and concentration, the local skin friction coefficients Cfx=τw/ρbxbν and Cfz=τz/ρbxbν, the local Nusselt number Nu=-(∂T/∂y)y=0/(Tw-T∞)b/ν, and the local Sherwood number Sh=(∂C/∂y)y=0/(Cw-C∞)b/ν, where τw=[(μ+K)(∂u/∂y)+KN]y=0 and τz=[(μ+K)(∂w/∂y)]y=0 are the shear stress at the wall. By using (2.9), these quantities can be expressed as:
(2.17)Cfx=(1+N1)f′′(0),Cfz=(1+N1)g′(0),Nu=-θ′(0),Sh=-ϕ′(0).


## 3. Generalization of the Successive Linearisation Method (SLM)

In this section we describe the basic idea behind the proposed method of successive linearisation method (SLM) and we derive formulas that can be used to implement the generalized SLM in any system of nonlinear boundary value problems. We consider a generaln-order nonlinear system of ordinary differential equations which is represented by the nonlinear boundary value problem of the form
(3.1)L[Y(x),Y′(x),Y′′(x),…,Y(n)(x)]+N[Y(x),Y′(x),Y′′(x),…,Y(n)(x)]=0,
where Y(x) is a vector of unknown functions, x is an independent variable, and the primes denote ordinary differentiation with respect to x. The functions L and N are vector functions which represent the linear and nonlinear components of the governing system of equations, respectively, defined by
(3.2)L=[L1(y1,y2,…,yk;y1′,y2′,…,yk′;…;y1(n),y2(n),…,yk(n))L2(y1,y2,…,yk;y1′,y2′,…,yk′;…;y1(n),y2(n),…,yk(n))⋮Lk(y1,y2,…,yk;y1′,y2′,…,yk′;…;y1(n),y2(n),…,yk(n))],N=[N1(y1,y2,…,yk;y1′,y2′,…,yk′;…;y1(n),y2(n),…,yk(n))N2(y1,y2,…,yk;y1′,y2′,…,yk′;…;y1(n),y2(n),…,yk(n))⋮Nk(y1,y2,…,yk;y1′,y2′,…,yk′;…;y1(n),y2(n),…,yk(n))],Y(x)=[y1(x)y2(x)⋮yk(x)],
where y1,y2,…,yk are the unknown functions. We define an initial guess Y0(x) of the solution of (3.1) as
(3.3)Y0(x)=[y1,0(x)y2,0(x)⋮yk,0(x)].
For illustrative purposes, we assume that (3.1) is to be solved for x∈[a,b] subject to the boundary conditions
(3.4)Y(a)=Ya,Y(b)=Yb,
where Ya and Yb are given constants. In previous implementation of the SLM (see, e.g.,  [19–25]) an appropriate initial guess was considered to be functions that satisfy the governing boundary conditions of (3.1). The general approach proposed in this study assumes that the initial approximation is a solution of the equation
(3.5)L[Y0(x),Y0′(x),Y0′′(x),…,Y0(n)(x)]=0,
which is solved subject to the underlying problem’s boundary conditions.Define a functionZ1(x) to represent the vertical difference between Y(x) and the initial guess Y0(x), that is,
(3.6)Z1(x)=Y(x)-Y0(x),orY(x)=Y0(x)+Z1(x).
For example, the vertical displacement between the variable y1(x) and its corresponding initial guess y1,0(x) is z1,1=y1(x)-y1,0(x). This is shown in Figure 2.Figure 2
Geometric representation of the functionz1,1(x).Substituting (3.6) in (3.1) gives
(3.7)L[Z1,Z1′,Z1′′,…,Z1(n)]+N[Y0+Z1,Y0′+Z1′,Y0′′+Z1′′,…,Y0(n)+Z1(n)]=-L[Y0,Y0′,Y0′′,…,Y0(n)].
Since Y0(x) is an known function, solving (3.7) would yield an exact solution for Z1(x). However, since the equation is nonlinear, it may not be possible to find an exact solution. We therefore look for an approximate solution which is obtained by solving the linear part of the equation assuming that Z1 and its derivatives are small. This assumption enables us to use the Taylor series method to linearise the equation. If Z1(x) is the solution of the full equation (3.7) we let Y1(x) denote the solution of the linearised version of (3.7). Expanding (3.7) using Taylor series (assuming Z1(x)≈Y1(x)) and neglecting higher order terms gives
(3.8)L[Y1,Y1′,Y1′′,…,Y1(n)]+[∂N∂Y1](Y0,Y0′,Y0′′,…,Y0(n))Y1+[∂N∂Y1′](Y0,Y0′,Y0′′,…,Y0(n))Y1′+[∂N∂Y1′′](Y0,Y0′,Y0′′,…,Y0(n))Y1′′+⋯+[∂N∂Y1(n)](Y0,Y0′,Y0′′,…,Y0(n))Y1(n)=-L[Y0,Y0′,Y0′′,…,Y0(n)]-N[Y0,Y0′,Y0′′,…,Y0(n)].
The partial derivatives inside square brackets in (3.8) represent Jacobian matrices of size k×k, defined as
(3.9)[∂N∂Yi]=[∂N1∂y1,i∂N1∂y2,i⋯∂N1∂yk,i∂N2∂y1,i∂N2∂y2,i⋯∂N2∂yk,i⋮⋮⋮∂Nk∂y1,i∂Nk∂y2,i⋯∂Nk∂yk,i],[∂N∂Yi(p)]=[∂N1∂y1,i(p)∂N1∂y2,i(p)⋯∂N1∂yk,i(p)∂N2∂y1,i(p)∂N2∂y2,i(p)⋯∂N2∂yk,i(p)⋮⋮⋮∂Nk∂y1,i(p)∂Nk∂y2,i(p)⋯∂Nk∂yk,i(p)],
where i=1 and p is the order of the derivatives.Since the right hand side of (3.8) 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 (3.8) is close to the solution of (3.7), that is, Z1(x)≈Y1(x), the current estimate (1st order) of the solution Y(x) is
(3.10)Y(x)≈Y0(x)+Y1(x).
To improve on this solution, we define a slack function, Z2(x), which when added to Y1(x) gives Z1(x) (see, e.g., Figure 3), that is
(3.11)Z1(x)=Z2(x)+Y1(x).
Since Y1(x) is now known (as a solution of (3.8) ), we substitute (3.11) in (3.7) to obtain
(3.12)L[Z2,Z2′,Z2′′,…,Z2(n)]+N[Y0+Y1+Z2+Y0′+Y1′+Z2′,…,Y0(n)+Y1(n)+Z2(n)]=-L[Y0+Y1,Y0′+Y1′,Y0′′+Y1′′,…,Y0(n)+Y1(n)].
Solving (3.12) would result in an exact solution for Z2(x). But since the equation is nonlinear, it may not be possible to find an exact solution. We therefore linearise the equation using Taylor series expansion and solve the resulting linear equation. We denote the solution of the linear version of (3.12) by Y2(x), such that Z2(x)≈Y2(x). Setting Z2(x)=Y2(x) and expanding (3.12), for small Y2(x) and its derivatives gives
(3.13)L[Y2,Y2′,…,Y2(n)]+[∂N∂Y2](Y0+Y1,Y0′+Y1′,…,Y0(n)+Y1(n))Y2+[∂N∂Y2′](Y0+Y1,Y0′+Y1′,…,Y0(n)+Y1(n))Y2′+[∂N∂Y2′′](Y0+Y1,Y0′+Y1′,…,Y0(n)+Y1(n))Y2′′+⋯+[∂N∂Y2(n)](Y0+Y1,Y0′+Y1′,…,Y0(n)+Y1(n))Y1(n)=-L[Y0+Y1,Y0′+Y1′,…,Y0(n)+Y1(n)]-N[Y0+Y1,Y0′+Y1′,…,Y0(n)+Y1(n)],
where the partial derivatives inside square brackets in (3.13) represent Jacobian matrices defined as in (3.9) with i=2.Figure 3
Geometric representation of the functionsz2,1.After solving (3.13), the current (2nd-order) estimate of the solution Y(x) is
(3.14)Y(x)≈Y0(x)+Y1(x)+Y2(x).
Next we define Z3(x) (see Figure 4) such that(3.15)Z2(x)=Z3(x)+Y2(x).
Equation (3.15) is substituted in the nonlinear (3.12) and the linearisation process described above is repeated. This process is repeated for m=3,4,5,…,i. In general, we have
(3.16)Zi(x)=Zi+1(x)+Yi(x).
Thus, Y(x) is obtained as
(3.17)Y(x)=Z1(x)+Y0(x),=Z2(x)+Y1(x)+Y0(x),=Z3(x)+Y2(x)+Y1(x)+Y0(x),⋮=Zi+1(x)+Yi(x)+⋯+Y3(x)+Y2(x)+Y1(x)+Y0(x),=Zi+1(x)+∑m=0iYm(x).
The procedure for obtaining each Zi(x) is illustrated in Figures 2, 3, and 4, respectively, for i=1,2,3.Figure 4
Geometric representation of the functionsz3,1.We note that wheni becomes large, Zi+1 becomes increasingly smaller. Thus, for large i, we can approximate the ith-order solution of Y(x) by
(3.18)Y(x)=∑m=0iYm(x)=Yi(x)+∑m=0i-1Ym(x).
Starting from a known initial guess Y0(x), the solutions for Yi(x) can be obtained by successively linearising the governing equation (3.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 linearised equation that is successively solved for Yi(x) is given by
(3.19)L[Yi,Yi′,Yi′′,…,Yi(n)]+a0,i-1Yi(n)+a1,i-1Yi(n-1)+⋯+an-1,i-1Yi′+an,i-1Yi=Ri-1(x),
where
(3.20)a0,i-1(x)=[∂N∂Yi(n)](∑m=0i-1Ym,∑m=0i-1Ym′,∑m=0i-1Ym′′,…,∑m=0i-1Ym(n)),a1,i-1(x)=[∂N∂Yi(n-1)](∑m=0i-1Ym,∑m=0i-1Ym′,∑m=0i-1Ym′′,…,∑m=0i-1Ym(n)),an-1,i-1(x)=[∂N∂Yi′](∑m=0i-1Ym,∑m=0i-1Ym′,∑m=0i-1Ym′′,…,∑m=0i-1Ym(n)),an,i-1(x)=[∂N∂Yi](∑m=0i-1Ym,∑m=0i-1Ym′,∑m=0i-1Ym′′,…,∑m=0i-1Ym(n)),Ri-1(x)=-L(∑m=0i-1Ym,∑m=0i-1Ym′,∑m=0i-1Ym′′,…,∑m=0i-1Ym(n))-N(∑m=0i-1Ym,∑m=0i-1Ym′,∑m=0i-1Ym′′,…,∑m=0i-1Ym(n)).


## 4. Numerical Solution

In this section we solve the governing (2.11)–(2.15) using the SLM method described in the last section. We note that (2.11)–(2.13) can be solved independently of equations (2.14)–(2.15). We begin by using the SLM approach to solve for f,g and h. We write the governing equations (2.11)–(2.13) as a sum of the linear and nonlinear components as
(4.1)-L[f,f′,f′′,f′′′,g,g′,g′′,h,h′,h′′]+N[f,f′,f′′,f′′′,g,g′,g′′,h,h′,h′′]=0,
where the primes denote differentiation with respect to η and
(4.2)L=[L1L2L3]=[f′′′+N1h′-Mαe2+βe2(αef′+βeReg)-1kpf′g′′+Mαe2+βe2(βeRef′-αeg)-1kpgGh′′-2h-f′′],N=[N̅1N̅2N̅3]=[ff′′-(f′)2fg′0].
Using (3.19), the general equation to be solved for Yi, where
(4.3)Yi=[fgh],
is
(4.4)L[Yi,Yi′,Yi′′,Yi′′′]+a0,i-1Yi′′′+a1,i-1Yi′′+a2,i-1Yi′+a3,i-1Yi=Ri-1(η),
subject to the boundary conditions
(4.5)fi(0)=fi′(0)=gi(0)=hi(0)=fi′(∞)=gi(∞)=hi(∞)=0,
where
(4.6)a0,i-1=[∂N̅1∂f′′′∂N̅1∂g′′′∂N̅1∂h′′′∂N̅2∂f′′′∂N̅2∂g′′′∂N̅2∂h′′′∂N̅3∂f′′′∂N̅3∂g′′′∂N̅3∂h′′′]=[000000000],a1,i-1=[∂N̅1∂f′′∂N̅1∂g′′∂N̅1∂h′′∂N̅2∂f′′∂N̅2∂g′′∂N̅2∂h′′∂N̅3∂f′′∂N̅3∂g′′∂N̅3∂h′′]=[∑fm00000000],a2,i-1=[∂N̅1∂f′∂N̅1∂g′∂N̅1∂h′∂N̅2∂f′∂N̅2∂g′∂N̅2∂h′∂N̅3∂f′∂N̅3∂g′∂N̅3∂h′]=[-2∑fm′000∑fm0000],a3,i-1=[∂N̅1∂f∂N̅1∂g∂N̅1∂h∂N̅2∂f∂N̅2∂g∂N̅2∂h∂N̅3∂f∂N̅3∂g∂N̅3∂h]=[∑fm′′00∑gm′00000],Ri-1=[r1,i-1r2,i-1r3,i-1],r1,i-1=-[∑fm′′′+∑fm∑fm′′+N1∑hm′-Mαe2+βe2(αe∑fm′+βeRe∑gm)-1kp∑fm′],r2,i-1=-[∑gm′′+∑fm∑gm′+Mαe2+βe2(βeRe∑fm′-αe∑gm)-1kp∑gm]r3,i-1=-[G∑hm′′-2∑hm-∑fm′′]
and the sums in (4.6) denote ∑=∑m=0i-1. Once each solution for fi, gi, hi (i≥1) has been found from iteratively solving equations (4.4)–(4.5), the approximate solutions for f(η), g(η), and h(η) are obtained as
(4.7)f(η)≈∑m=0ifm(η),g(η)≈∑m=0igm(η),h(η)≈∑m=0ihm(η),
where i is the order of SLM approximation. Since the coefficient parameters and the right hand side of (4.4)), for i=1,2,3,…, are known (from previous iterations). The equation system is solved using the Chebyshev spectral collocation method. 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
(4.8)ξj=cos⁡πjN,j=0,1,…,N,
where N+1 is the number of collocation points used (see, e.g., [26–28]). In order to implement the method, the physical region [0,∞) is transformed into the region [-1,1] using the domain truncation technique in which the problem is solved on the interval [0,L] instead of [0,∞). This leads to the mapping
(4.9)ηL=ξ+12,-1≤ξ≤1,
where L is the scaling parameter used to invoke the boundary condition at infinity. The unknown functions fi, θi and ϕi are approximated at the collocation points by
(4.10)fi(ξ)≈∑k=0Nfi(ξk)Tk(ξj),gi(ξ)≈∑k=0Ngi(ξk)Tk(ξj),hi(ξ)≈∑k=0Nhi(ξk)Tk(ξj),j=0,1,…,N,
where Tk is the kth Chebyshev polynomial defined as
(4.11)Tk(ξ)=cos⁡[kcos⁡-1(ξ)].
The derivatives of the variables at the collocation points are represented as
(4.12)dafidηa=∑k=0NDkjafi(ξk),dagidηa=∑k=0NDkjagi(ξk),dahidηa=∑k=0NDkjahi(ξk),j=0,1,…,N,
where a is the order of differentiation and D=(2/L)D with D being the Chebyshev spectral differentiation matrix (see, e.g., [26, 28]). Substituting (4.9)–(4.12) in (4.4))–((4.5) leads to the matrix equation given as
(4.13)Mi-1Yi=Ri-1,
subject to the boundary conditions
(4.14)fi(ξN)=∑k=0ND0kfi(ξk)=∑k=0NDNkfi(ξk)=gi(ξN)=gi(ξ0)=hi(ξN)=hi(ξ0)=0,
where
(4.15)Mi-1=A+a1,i-1D̅2+a2,i-1D̅+a3,i-1,(4.16)A=[D3-(Mαeαe2+βe2+1kp)D-M(αe2+βe2)ReIN1DMβeReαe2+βe2DD2-(Mαeαe2+βe2+1kp)O-D2OGD2-2I],(4.17)D̅=[DOOODOOOD],
and Yi and Ri-1 are (3N+1)×1 column vectors defined by
(4.18)Yi=[FiGiHi],Ri-1=[r1,i-1r2,i-1r3,i-1],
where
(4.19)Fi=[fi(ξ0),fi(ξ1),…,fi(ξN-1),fi(ξN)]T,Gi=[gi(ξ0),gi(ξ1),…,gi(ξN-1),gi(ξN)]T,Hi=[hi(ξ0),hi(ξ1),…,hi(ξN-1),hi(ξN)]T,r1,i-1=[r1,i-1(ξ0),r1,i-1(ξ1),…,r1,i-1(ξN-1),r1,i-1(ξN)]T,r2,i-1=[r2,i-1(ξ0),r2,i-1(ξ1),…,r2,i-1(ξN-1),r2,i-1(ξN)]T,r3,i-1=[r3,i-1(ξ0),r3,i-1(ξ1),…,r3,i-1(ξN-1),r3,i-1(ξN)]T.
In the above definitions, ak,i-1, (k=1,2,3) are now diagonal matrices of size 3(N+1)×3(N+1) and the superscript T is the transpose.To impose the boundary conditions (4.14) on (4.13) we begin by splitting the matrix M in (4.13) into 9 blocks each of size (N+1)×(N+1) in such a way that M takes the form
(4.20)Mi-1=[M11M12M13M21M22M23M31M32M33].
We then modify the first and last rows of Mmn (m,n=1,2,3) and rm,i-1 and the N-1th row of M1,1,M1,2,M1,3 in such a way that the modified matrices Mi-1 and Ri-1 take the form
(4.21)Mi-1=(D0,0D0,1⋯D0,N-1D0,N00⋯0000⋯00M11M12M13DN,0DN,1⋯DN,N-1DN,N00⋯0000⋯0000⋯0100⋯0000⋯0000⋯0010⋯0000⋯00M21M22M2300⋯0000⋯0100⋯0000⋯0000⋯0010⋯00M31M32M3300⋯0000⋯0000⋯01),Ri-1=(0r1,i-1(ξ1)⋮r1,i-1(ξN-2)000r2,i-1(ξ1)⋮r2,i-1(ξN-2)r2,i-1(ξN-1)00r3,i-1(ξ1)⋮r3,i-1(ξN-2)r3,i-1(ξN-1)0).
After modifying the matrix system (4.13) to incorporate boundary conditions, the solution is obtained as
(4.22)Yi=Mi-1-1Ri-1.
We use (3.5) to solve solve for the initial approximation Y0. If we use the Chebyshev spectral method to solve for Y0, we arrive at the following:
(4.23)AY0=0,
subject to the boundary conditions
(4.24)∑k=0NDNkf0(ξk)=1,f0(ξN)=∑k=0ND0kf0(ξk)=g0(ξN)=g0(ξ0)=h0(ξN)=hi(ξ0)=0,
where A is as defined in (4.16) and Y0 are (3N+1)×1 column vectors defined by
(4.25)Y0=[F0G0H0],
where
(4.26)F0=[f0(ξ0),f0(ξ1),…,f0(ξN-1),f0(ξN)]T,G0=[g0(ξ0),g0(ξ1),…,g0(ξN-1),g0(ξN)]T,H0=[h0(ξ0),h0(ξ1),…,h0(ξN-1),h0(ξN)]T.
To impose the boundary conditions (4.24) on the system (4.23) we begin by splitting the matrix A into 9 blocks each of size (N+1)×(N+1) in such a way that A takes the form
(4.27)A=[A11A12A13A21A22A23A31A32A33].
We then modify the first and last rows of Amn (m,n=1,2,3) and the N-1th row of A1,1,A1,2,A1,3 in such a way that the modified matrix A takes the form;
(4.28)A=(D0,0D0,1⋯D0,N-1D0,N00⋯0000⋯00A11A12A13DN,0DN,1⋯DN,N-1DN,N00⋯0000⋯0000⋯0100⋯0000⋯0000⋯0010⋯0000⋯00A21A22A2300⋯0000⋯0100⋯0000⋯0000⋯0010⋯00A31A32A3300⋯0000⋯0000⋯01).
We introduce a matrix R0 on the right hand side of (4.23) such that
(4.29)AY0=R0,
where
(4.30)R0=(00⋯010|00⋯000|00⋯000)T.
The solution Y0 is obtained as
(4.31)Y0=A-1R0.
Thus, starting from Y0, the solutions for Yi are then obtained iteratively from solving equation (4.22). Once the solution for f and g and h have been obtained, the solutions for θ and ϕ can be found using the same procedure as described above.


## 5. Results and Discussion

In this section, we present the results obtained using both the successive linearization method and Chebyshev spectral collocation numerical method. The number of collocation points employed in this investigation isN=50. The results are also validated against those obtained using the Matlab built-in solver bvp4c as well as reported by Elgazery [12]. Table 1 represents comparison between the SLM, Chebyshev pseudospectral method [12], and the bvp4c results of the local Nusselt number -θ′(0). We note from this table that the SLM and the bvp4c results match exactly but the results reported by Elgazery [12] are slightly different from the current results. This is an interesting observation, since Elgazery [12] also applied the pseudospectral method, this substantiate the claim that the SLM technique improves the accuracy of the Chebyshev method. We also note in this table that the absolute values of the Nusselt number decrease with increasing values of magnetic parameter and the ion-slip parameter. However, it increases with increasing values of the Hall parameter.Table 1
Comparison of the present results (4th-order) for the local Nusselt number-θ′(0) against the bvp4c results and the results of [12] for various values of βe, βi and M at Pr=0.72,β2=0.5,fw=0.1, kp=2, Re=1, G=2, N1=0.2, and Ec=0.02.
Mβeβi[12]Presentbvp4c050.40.3274980.327495280.327495280.50.3248440.324841640.3248416410.3219120.321909120.321909121.50.3187720.318768090.318768090.300.3112340.311228850.3112288510.3192970.319293920.3192939220.3230830.323080610.323080610.3500.3266960.326693280.326693280.20.3262250.326222740.326222740.50.3258650.325862310.32586231Table2 gives the results of the wall stresses (Cfx) for various values of kp,N1,βe, and βi when other parameters are fixed. For fixed values of N1, we observe that the absolute values of the local-skin friction decrease with the increasing values of the permeability parameter (kp). The local shear stress (Cfx) increases in absolute values as the coupling parameter increases. The local wall shear stress (Cfx) increases in absolute values when the Hall parameter or ion-slip parameter increases.Table 2
The wall shear stressesCfx=(1+N1)f′′(0) for various values of kp, N1, βe, and βi at M=0.3, Pr=0.72,fw=0.1, Re=1, G=2, and Ec=0.02.
kpN1βeβi2nd-Order3rd-Order4th-Order6th-Order0.1150.4−6.65334091−6.65334091−6.65334091−6.653340910.5150.4−3.46988916−3.46989224−3.46989224−3.469892241150.4−2.83779413−2.83787814−2.83787815−2.837878152150.4−2.46510143−2.46622986−2.4662317−2.466231720.150.4−1.40914015−1.40961763−1.40961811−1.4096181120.550.4−1.88961638−1.89034144−1.89034233−1.890342332150.4−2.46510143−2.46622986−2.4662317−2.46623172250.4−3.52816005−3.53067494−3.53068356−3.5306835620.500.4−2.04791788−2.04807626−2.04807629−2.0480762920.520.4−1.91920286−1.91972357−1.91972398−1.9197239820.540.4−1.89445825−1.89514419−1.89514497−1.8951449720.560.4−1.88649265−1.88724415−1.88724512−1.8872451220.550−1.88077116−1.88155788−1.88155895−1.8815589520.552−1.88698782−1.88774243−1.88774341−1.8877434120.554−1.8814535−1.88225343−1.88225455−1.8822545520.556−1.8789626−1.87978367−1.87978486−1.87978486The values of the wall shear stresses(Cfx,Cfz), the local Nusselt number -θ′(0) and the local Sherwood number -ϕ′(0) for different Hartman number M, Hall parameter βe, ion-slip parameter βi, and coupling parameter N1 are presented in Table 3. It is observed in this table that the local skin-friction coefficients f′′(0) and g′(0) increase with increasing values of M and N1 but g′(0) decreases as values of βe and βi increase. Both the local Nusselt number and the Sherwood number are found to be decreasing with increasing values of the Hartman number M and they also slightly decrease with increasing values of the ion-slip parameter and the coupling parameter. The local Nusselt number and the Sherwood number increases as Hall currents increase.Table 3
The local Nusselt number-θ′(0), Sherwood number -ϕ′(0), and the wall shear stresses Cfx=(1+N1)f′′(0) and Cfz=(1+N1)g′(0) for various values of βe, βi, M and N1 at Pr=0.72,fw=0.1, kp=2, Re=1, G=2 and Ec=0.02.
MβeβiN1CfxCfz-θ′(0)-ϕ′(0)150.40.2−1.563229140.085991360.321909120.48724548250.40.2−1.610774020.164225720.315478470.48232094350.40.2−1.659723360.235543360.308643360.47715106450.40.2−1.709225640.30089270.301698440.47195166550.40.2−1.758728950.36115650.294828950.46685360.300.40.2−1.6581241900.311228850.479411840.320.40.2−1.555081220.048859960.323080610.488183230.340.40.2−1.535302140.031718990.325472470.490013020.360.40.2−1.528942060.022919910.326242050.490608330.3500.2−1.524377450.035104970.326693280.490989710.350.50.2−1.532225190.024303180.325862310.490308870.3510.2−1.532737380.014833090.325839130.490279230.351.50.2−1.531103550.009319870.326037780.490432660.350.50−1.287381390.020385920.327578640.491380140.350.50.2−1.532225190.024303180.325862310.490308870.350.50.5−1.891325860.030062730.323075380.488586110.350.51−2.467535290.039313680.317755920.48535366Graphical representation of the numerical results are illustrated in Figure5 through Figure 13 to depict the influence of different parameters on the flow characteristics. Figure 5 depicts the effects of the Hall current parameter βe on the velocity components f′(η) and g(η). We see that the stream velocity profiles f′(η) increases as βe increases. We also observe in this figure that velocity distribution across the stretching sheet g(η) increases with increasing values of βe when βe≤1 but decreases with increasing values of βe greater than unity.(a)(b)Figure6 depicts the effects of the Hall current parameter βe on the angular velocity h(η), temperature and concentration distributions. We observe that the angular velocity steeply rises up to maximum peaks as the Hall parameter increases. In the same figure we see that both the temperature and concentration profiles approach their classical values when the Hall parameter βe becomes large. They both decrease with increasing values of βe.(a)(b)(c)The influence of the Hartman number on the stream velocity and the velocity across the plate is depicted in Figure7. It is observed that the stream velocity of the fluid decreases with the increase of the magnetic field parameter values. The decrease in this velocity component as the Hartman number (M) increases is because the presence of a magnetic field in an electrically conducting fluid introduces a force called the Lorentz force, which acts against the flow if the magnetic field is applied in the normal direction, as in the present study. This resistive force slows down the fluid velocity component as shown in Figure 7. In Figure 7 we have the influence of the magnetic field parameter on the lateral velocity. It can be seen that as the values of this parameter increase, the lateral velocity increases.(a)(b)Figure8 shows the effects of the magnetic parameter M on the angular velocity h(η), temperature and concentration profiles. As expected, the angular velocity steeply increases with every value of the magnetic value until attaining a peak and thereafter the angular velocity decreases monotonically with increasing values of M. In this figure, we also observe that both the temperature and concentration boundary layers become thick as values of the magnetic parameter increase. The effects of a transverse give rise to a resistive-type force called the Lorentz force. This force has the tendency to slow down the motion of the fluid and increase its thermal and concentration boundary layers hence increasing the temperature and concentration fields of the flow.(a)(b)(c)Figures9 and 10 give the effects of the permeability parameter kp on the velocity components, temperature and concentration distributions. As shown in these figures, all the velocity components are increasing with increasing values of the permeability parameter, whereas the temperature and concentration decrease as the permeability parameter kp increases. Physically, this means that the porous medium impact on the boundary layer growth is significant due to the increase in the thickness of the thermal and concentration boundary layers. It is expected that, an increase in the permeability of porous medium lead to a rise in the flow of the fluid through it, since when the holes of the porous medium become large, the resistance of the medium may be neglected. Figure 11 depicts the effects of coupling constant or material parameter on the velocity components. The stream velocity is slightly affected by the coupling constant. As we move away from the stretching sheet surface, η>1, the stream velocity decreases as N1 increases. We observe also in Figure 11 that the velocity across the plate monotonically decreases as N1 increases. We observe in Figure 11 that for values of η less than two, the angular velocity increases when N1 increases, and thereafter we have a cross-flow angular velocity when increasing values of N1 cause the angular velocity to decrease.(a)(b)(a)(b)(c)(a)(b)(c)Figure12 shows the effects of the Eckert number Ec and the variable thermal diffusivity parameter β2 on the temperature. The temperature distribution θ(η) increases as Ec and β2 increase as shown in Figure 12.(a)(b)(a)(b)Finally, Figure13 shows the effects of the ion-slip parameter βi on the velocity across the plate g(η) and those of the chemical reaction on the concentration distribution ϕ(η). We observe in this figure that the ion-slip parameter has a significant effect on the induced velocity in the z-direction. The velocity component decreases with increases in the parameter βi. We also observe in this study (not shown for brevity) that f′(η),h(η), and θ(η) profiles increase as the ion-slip parameter increases. From Figure 13 we also observe that the effects of chemical reaction parameter is to reduce the concentration distribution ϕ(η) when γ increases.


## 6. Conclusions

In this study, the effects of Hall and ion-slip currents, chemical reaction and variable thermal diffusivity on magnetomicropolar fluid flow through a porous medium past a stretching sheet in the presence of heat and mass transfer have been numerically analyzed using the successive linearization method together with Chebyshevcollocation method. Results for the wall stresses, local Nusselt and Sherwood numbers as well as the details of the velocities, temperature and concentration distributions are presented in tabular and/or graphical forms for the governing parameters. The successive linearization method was found to be a very efficient and accurate method which we hope to apply to different problems in fluid mechanics and related fields. The study observed that higher values of the coupling parameter results in lower stream velocity and velocity across the plate but have little effects on the angular velocity. The porous medium impact on the boundary layer growth is significant due to the increase in the thickness of the hydrodynamic boundary layer and the decrease in the thickness of the thermal and concentration boundary layers. All the velocity components increase with increasing values of the Hall current. The ion-slip current causes the induced velocity in thez-direction to decrease but has opposing effects on the stream velocity, angular velocity, as well as temperature profiles.



---

*Source: 689015-2012-02-16.xml*