Effects of radiation on MHD free convection of a Casson fluid from a horizontal circular cylinder with partial slip in non-Darcy porous medium with viscous dissipation

Abstract

In the present study, the effects of radiation on MHD free convection from a cylinder with partial slip in a Casson fluid in non-Darcy porous medium is investigated. The surface of the cylinder is heated under constant surface temperature with partial slip. Partial slip factors are considered on the surface for both velocity and temperature. The boundary layer equations are normalized into a system of non-similar partial differential equations and are then solved using a bi-variate quasilinearization method (BQLM). The boundary layer velocity and temperature profiles are computed for different values of the physical parameters. Increasing the Forchheimer parameter decreases the temperature profiles. The decrease of the velocity profiles with the increase in magnetic parameter is more enhanced in the presence of the velocity slip factor. Increasing the Eckert number increases the temperature profiles in both suction and blowing cases. This study considers the unique problem of the effect of transpiration in a Casson fluid in the presence of radiation, a magnetic field, and viscous dissipation. The results obtained in this study are compared with other numerical methods and were found to be in excellent agreement.

1 Introduction

The flow of non-Newtonian fluids is applied in many situations in industry such as processing of materials and chemical engineering. These fluids show different characteristics from the Newtonian fluids which cannot be fully represented by the Navier-Stokes equations. To represent these non-Newtonian fluids some modifications to the Navier-Stokes equations are used and these are seen in many research works which studied viscoelastic and micropolar fluids [1, 2]. These fluids are categorized as viscoelastic, thixotropic, and power-law fluids. The constitutive equations of such fluids cannot fully represent the actual behavior of these fluids. These fluids include contaminated lubricants, molten metal, synovial fluids, etc.

The study of radiation effects on MHD free convection of a Casson fluid in porous medium with partial slip is an important aspect due to its practical application in the design of automatic cooking machines and the design of internal engine parts in mechanical engineering. Other examples arise in petroleum production, multiphase mixtures, pharmaceutical formulations, coal in water, paints, lubricants, jams, sewage, soup, blood. There has been a significant improvement in the study of non-Newtonian fluids in which many different situations have been considered. Studies in a Casson fluid include work by among others Mukhopadhyay et al. [3] and Nadeem et al. [4].

Several research workers have studied radiation effects, these include Kameswaran et al. [5] who studied radiation effects on hydromagnetic Newtonian liquid flow due to an exponentially stretching sheet. They studied radiation effects in the presence of a magnetic field, advancing the studies of radiation effects in Newtonian fluids. Shateyi and Marewo [6] investigated numerical analysis of MHD stagnation point flow of a Casson fluid, they considered thermal radiation in their work. Chamkha et al. [7] studied thermal radiation effects on MHD forced convection flow adjacent to a non-isothermal wedge in the presence of a heat source or sink. Pramanik [8] studied Casson fluid flow and heat transfer past an exponentially porous stretching surface in the presence of thermal radiation.

The study of fluid flow in the presence of a magnetic fluid has also been performed by many authors, among others Ece [9], who investigated free convection flow about a cone under mixed thermal boundary conditions and a magnetic field. Narayana et al. [10] studied free magnetohydrodynamic flow and convection from a vertical spinning cone with cross diffusion effects. Nadeem et al. [11] studied numerically MHD boundary layer flow of a Maxwell fluid past a stretching sheet in the presence of nanoparticles. Chen [12] investigated the combined heat and mass transfer in MHD free convection from a vertical surface with ohmic heating and viscous dissipation.

The study of fluid flow past a cylindrical geometry was performed by among others Anwa et al. [13], who studied mixed convection boundary layer flow of a viscoelastic fluid over a horizontal cylinder. Deka et al. [14] investigated transient free convection flow past an accelerated vertical cylinder in a rotating cylinder. Ribeiro et al. [15] studied viscoelastic flow past a confined cylinder with three dimensional effects and stability. Patel and Chhabra [16] studied steady flow of Bingham plastic fluids past an elliptical cylinder.

Studies in porous media and viscous dissipation have been carried out by among others, Awad et al. [17] who studied natural convection of viscoelastic fluid from a cone embedded in a porous medium with viscous dissipation. Awad et al. [17] investigated convection from an inverted cone in a porous medium with cross diffusion effects. Hayat et al. [18] studied heat and mass transfer for Soret and Dufour effects on mixed convection boundary layer flow over a stretching vertical surface in a porous medium filled with a viscoelastic fluid. Cheng [19] studied Soret and Dufour effects on free convection boundary layer over a vertical cylinder in a saturated porous medium. Chamkha and Rashad [20] investigated natural convection from a vertical permeable cone in nanofluid saturated porous media for uniform heat and nanoparticles volume fraction fluxes.

In this study we investigate the effects of radiation in MHD free convection of a Casson fluid from a horizontal circular cylinder with partial slip in non-Darcy porous medium with viscous dissipation. The surface of the cylinder is perforated in which we have the effects of transpiration which acts transversely in the direction ξ. The force which causes transpiration is sometimes called the Forchheimer drag force term $-\xi\Lambda f^{\prime2}$, which appears in the momentum equation (12). This term is also associated with the geometry of the porous medium. In this work we extended the work of Ramachandra et al. [1] in which we introduced magnetic field, radiation effects and viscous dissipation. We also solved the system of the resulting partial differential equations by the bi-variate quasilinearization method (BQLM). The method is described in detail in Motsa et al. [21]. We study the effects of the magnetic field, radiation. and the Eckert number on velocity and temperature profiles for different values of the Casson parameter β, the local inertia (Forchheimer parameter) Λ, and the Darcian drag force coefficient $\Lambda_{1}$. The results of this work are validated by comparison with other methods, which are the successive linearization method (SLM) and MATLAB’s routine ‘bvp4c’.

2 Mathematical formulation

The steady, laminar, two dimensional MHD free convection of a Casson fluid flow from a horizontal circular cylinder with partial slip in a non-Darcy porous medium, with viscous dissipation and radiation effects, is considered. The fluid is maintained at a uniform temperature $T_{w}$ ($>T_{\infty}$), the transpiration velocity $V_{w}$ is in the direction ξ. x is the tangential coordinate of the circle and y is the radial coordinate of the circle. u and v are the velocity components in the x and y directions, respectively, as shown in Figure 1.

Figure 1

Schematic of the problem circle in a Casson fluid with transpiration.

The rheological equation of state for an isotropic and incompressible flow of a Casson fluid is given as in [1–3] by

$$ \tau_{ij}=\left \{ \begin{array}{l@{\quad}l} 2(\mu_{B} + \frac{P_{y}}{\sqrt{2\pi}})e_{ij}, & \pi> \pi _{c}, \\ 2(\mu_{B} + \frac{P_{y}}{\sqrt{2\pi_{c}}})e_{ij}, & \pi< \pi_{c}, \end{array} \right . $$

(1)

$\pi= e_{ij}e_{ij}$ and $e_{ij}$ is the $(i,j)$th component of the deformation rate, π is the product of the deformation rate with itself, $\pi_{c}$ is a critical value of this product based on the non-Newtonian model, $\mu_{B}$ is the plastic dynamic viscosity of the non-Newtonian fluid, $P_{y}$ is the yield stress of the fluid.

The governing equations in this flow are given as

$$\begin{aligned}& \frac{\partial}{\partial x}(ru)+\frac{\partial}{\partial y}(rv)=0, \end{aligned}$$

(2)

$$\begin{aligned}& u\frac{\partial u}{\partial x}+v\frac{\partial u}{\partial y} =\nu\biggl(1+ \frac{1}{\beta}\biggr)\frac{\partial ^{2}u}{\partial y^{2}}+ \bar{g}\beta_{T}(T-T_{\infty}) \sin\biggl(\frac{x}{a}\biggr) - \Gamma u^{2} -\frac{\nu}{\bar{K}}u - \frac{\sigma B_{0}^{2}}{\rho}u, \end{aligned}$$

(3)

$$\begin{aligned}& u\frac{\partial T}{\partial x}+v\frac{\partial T}{\partial y}=\alpha\frac{\partial^{2}T}{\partial y^{2}} - \frac {1}{\rho C_{p}}\frac{\partial q_{r}}{\partial y} +\frac{\nu}{\rho C_{p}} \biggl(1+ \frac{1}{\beta}\biggr) \biggl(\frac{\partial u}{\partial y} \biggr)^{2} + \frac{\sigma B_{0}^{2}u^{2}}{\rho C_{p}}, \end{aligned}$$

(4)

where u and v are the velocity components in the x and y directions, respectively, a is the radius of the cylinder, ν is the kinematic viscosity of the Casson fluid, $\beta=\mu _{B}\sqrt{2\pi_{c}}/P_{y}$ is the non-Newtonian Casson parameter, $\alpha= \kappa/\rho C_{p}$ is the thermal diffusivity, κ is thermal conductivity of the fluid, $q_{r}$ is the radiative heat flux, $C_{p}$ is the specific heat, $\bar {g}$ is the acceleration due to gravity, $\beta_{T}$ is the coefficient of thermal expansions, and T is the temperature of the fluid. $T_{\infty}$ is the free stream temperature, ρ is the fluid density, $\bar{K}$ and Γ are the permeability and inertia coefficients respectively. σ is the electric conductivity and $B_{0}$ is the magnetic flux density. The Rosseland approximation for radiation may be written as follows:

$$ q_{r} = -\frac{4\sigma^{*}}{3k^{*}}\frac{\partial T^{4}}{\partial y}, $$

(5)

where $\sigma^{*}$ is the Stefan-Boltzmann constant and $k^{*}$ is the absorption coefficient. If the temperature difference within the flow is such that $T^{4}$ may be expanded in a Taylor series about $T_{\infty}$ and neglecting higher powers, we obtain $T^{4}-4T^{3}_{\infty} - 3T^{4}_{\infty}$ and therefore (4) can be written as

$$ u\frac{\partial T}{\partial x}+v\frac{\partial T}{\partial y}=\alpha\frac{\partial^{2}T}{\partial y^{2}} + \frac {16\sigma^{*} T^{3}_{\infty} }{3\rho C_{p}k^{*}}\frac{\partial^{2} T}{\partial y^{2}} +\frac{\nu}{\rho C_{p}} \biggl(1+ \frac{1}{\beta}\biggr) \biggl(\frac{\partial u}{\partial y} \biggr)^{2} + \frac{\sigma B_{0}^{2}u^{2}}{\rho C_{p}}. $$

(6)

The boundary conditions are given as

$$\begin{aligned}& u = N_{0}\biggl(1+\frac{1}{\beta}\biggr) \frac{\partial u}{\partial y},\qquad v=-V_{w},\qquad T=T_{w}+K_{0} \frac{\partial T}{\partial Y}, \qquad y=0, \end{aligned}$$

(7)

$$\begin{aligned}& u\rightarrow0, \qquad T\rightarrow T_{\infty},\quad \text{as } y\rightarrow\infty, \end{aligned}$$

(8)

where $N_{0}$ is the velocity slip factor, $K_{0}$ is the thermal slip factor. $N_{0}=K_{0}=0$ corresponds to no-slip conditions. The subscripts w and ∞ refer to surface and ambient conditions, respectively.

We introduce the non-dimensional variables

$$ \left. \begin{aligned} &\xi= \frac{x}{a}, \qquad \eta=\frac{y}{a}\mathit{Gr}^{\frac{1}{4}},\qquad \mathit{Pr}=\frac {\nu}{\alpha}, \\ &\mathit{Gr}=\frac{g\beta_{T}(T_{w}-T_{\infty})a^{3}}{\nu^{2}},\qquad \beta _{c}=\mu_{B} \frac{\sqrt{2\pi_{c}}}{P_{y}} , \\ &\Lambda= \Gamma a,\qquad \Lambda_{1}= \frac{1}{\mathit{Da} \mathit{Gr}^{\frac {1}{2}}}, \qquad \mathit{Da} = \frac{\bar{K}}{a^{2}}, \qquad f_{w}=-\frac{V_{w}a}{\nu \mathit{Gr}^{\frac{1}{4}} }, \qquad U_{0}=\frac{\nu \mathit{Gr}^{\frac{1}{2}} }{a}. \end{aligned} \right\} $$

(9)

Introducing the stream function ψ and similarity variables

$$\begin{aligned}& u =\frac{\partial\psi}{\partial y}\quad \text{and}\quad v=-\frac {\partial\psi}{\partial x}, \end{aligned}$$

(10)

$$\begin{aligned}& f(\xi,\eta)=\frac{\psi}{\nu\xi \mathit{Gr}^{\frac{1}{4}}},\qquad \theta (\xi,\eta )= \frac{T-T_{\infty}}{T_{w}-T_{\infty}}, \end{aligned}$$

(11)

using the stream function defined in (10) and the similarity variable in (11), (2)-(6) together with the boundary conditions (7) and (8) reduce to the following system of partial differential equations:

$$\begin{aligned}& \biggl(1+\frac{1}{\beta} \biggr)f''' + ff'' - (1+\Lambda\xi) f^{\prime2} - \bigl(\Lambda_{1} + M^{2}\bigr)f' + \frac{\sin(\xi)}{\xi} \theta=\xi \biggl(f'\frac {\partial f'}{\partial\xi} - f''\frac{\partial f}{\partial\xi } \biggr), \end{aligned}$$

(12)

$$\begin{aligned}& \frac{1}{\mathit{Pr}}\biggl(1+\frac{4}{3}K\biggr) \theta'' + f\theta' + \xi^{2} \mathit{Ec} \biggl[ \biggl(1+\frac{1}{\beta} \biggr)f^{\prime\prime2} + M^{2}\bigl(f'\bigr)^{2} \biggr]=\xi \biggl(f'\frac{\partial\theta}{\partial\xi} - \theta'\frac{\partial f}{\partial\xi} \biggr) \end{aligned}$$

(13)

with boundary conditions

$$\begin{aligned}& \eta= 0, \qquad f = f_{w}, \qquad f'= \biggl(1+\frac{1}{\beta} \biggr)S_{f}f'', \qquad \theta= 1 + S_{T}\theta', \end{aligned}$$

(14)

$$\begin{aligned}& \eta\rightarrow\infty,\qquad f'\rightarrow0, \qquad \theta \rightarrow0, \end{aligned}$$

(15)

where β is the Casson parameter, Λ is the Forchheimer parameter, $\Lambda_{1}$ is the Darcian drag coefficient, M is the magnetic field parameter, K is the radiation parameter, Pr is the Prandtl number, Ec is the Eckert number, $f_{w}>0$ corresponds to suction, and $f_{w}<0$ corresponds to blowing, $S_{f} = N_{0}\mathit{Gr}^{\frac {1}{4}}/L$ is the velocity slip factor, and $S_{T} = k\mathit{Gr}^{\frac {1}{4}}/L$ is the thermal slip factor. In the above equations the primes refer to the derivative with respect to η. The engineering parameters of interest are the local skin friction and the local Nusselt numbers, which are defined as follows.

The shear stress at the surface of the cone is given by

$$ \tau_{w} = \mu \biggl[ \biggl(1+\frac{1}{\beta} \biggr)\frac {\partial u}{\partial y} \biggr]_{y=0}=\frac{\mu (1+\frac{1}{\beta } )\nu \xi \mathit{Gr}^{\frac{3}{4}} }{a^{2}} f''(0), $$

(16)

where μ is the coefficient of viscosity. The skin friction coefficient is given by

$$ C_{f} = \frac{\tau_{w}}{\frac{1}{2}\rho U_{0}^{2}}. $$

(17)

Using (16) and (17) together with (10) and (11) give

$$ C_{f}\mathit{Gr}^{-\frac{3}{4}}=\biggl(1+\frac{1}{\beta}\biggr)\xi f''(0). $$

(18)

The heat transfer from the surface of the circle into the fluid is given by

$$ q_{w}=-k \biggl[\frac{\partial T}{\partial y} \biggr]_{y=0}=\frac {-k(T_{w}-T_{\infty})}{a\mathit{Gr}^{-\frac{1}{4}}}X\theta'(0), $$

(19)

k is the thermal conductivity of the fluid. The Nusselt number under LST is given by

$$ \mathit{Nu} = \frac{a}{k}\frac{q_{w}}{T_{w}-T_{\infty}}. $$

(20)

Equations (19) and (20) together with (10) and (11) give

$$ \mathit{Nu}\mathit{Gr}^{-\frac{1}{4}} =-\theta'(0). $$

(21)

3 Solution method

In this section we describe the implementation of the bi-variate quasilinearization method (BQLM) which is based on the quasilinearization method (QLM) which is described in detail in Motsa et al. [21]. We apply the quasilinearization method (QLM) first proposed by Bellman and Kalaba [22] to (12)-(13), which is based on the Taylor series expansion with the assumption that the differences $(f_{r+1}-f_{r})$, $(\theta _{r+1}-\theta_{r})$, and all the derivatives are small. We obtain the following equations:

$$\begin{aligned}& \biggl(1+\frac{1}{\beta} \biggr) f'''_{r+1}+a_{1,r}( \eta,\xi) f''_{r+1} + a_{2,r}(\eta, \xi) f'_{r+1}+a_{3,r}(\eta,\xi) f_{r+1} \\& \quad {}+a_{4,r}(\eta,\xi)\frac{\partial f'_{r+1}}{\partial\xi} + a_{5,r}(\eta, \xi)\frac{\partial f_{r+1}}{\partial\xi} =a_{6,r}(\eta,\xi), \end{aligned}$$

(22)

$$\begin{aligned}& \frac{1}{\mathit{Pr}} \biggl(1+\frac{4}{3}K \biggr)\theta ''_{r+1}+b_{1,r}(\eta,\xi ) \theta'_{r+1} + b_{2,r}(\eta,\xi) \frac{\partial\theta _{r+1}}{\partial \xi}=b_{3,r}(\eta,\xi), \end{aligned}$$

(23)

where

$$\begin{aligned}& a_{1,r}=f_{r} + \xi\frac{\partial f_{r}}{\partial\xi}, \end{aligned}$$

(24)

$$\begin{aligned}& a_{2,r}=- \biggl[2(1+\Lambda\xi)f'_{r} + \Lambda_{1} + M^{2} + \xi \frac {\partial f'_{r}}{\partial\xi} \biggr], \end{aligned}$$

(25)

$$\begin{aligned}& a_{3,r}=f''_{r} , \end{aligned}$$

(26)

$$\begin{aligned}& a_{4,r}=-\xi f'_{r} , \end{aligned}$$

(27)

$$\begin{aligned}& a_{5,r}=\xi f''_{r} , \end{aligned}$$

(28)

$$\begin{aligned}& a_{6,r}=f_{r}f''_{r} -(1+\Lambda\xi) \bigl(f'_{r}\bigr)^{2} - \frac{\sin\xi }{\xi }\theta_{r} - \xi \biggl(\frac{\partial f'_{r}}{\partial\xi }-f''_{r} \frac {\partial f_{r}}{\partial\xi} \biggr), \end{aligned}$$

(29)

$$\begin{aligned}& b_{1,r}=f_{r} + \xi\frac{\partial f_{r}}{\partial\xi}, \end{aligned}$$

(30)

$$\begin{aligned}& b_{2,r}=-\xi f', \end{aligned}$$

(31)

$$\begin{aligned}& b_{3,r}=-\mathit{Ec} \xi^{2} \biggl[\biggl(1+ \frac{1}{\beta}\biggr) \bigl(f''_{r} \bigr)^{2} + M^{2}\bigl(f'_{r} \bigr)^{2} \biggr]. \end{aligned}$$

(32)

The solution for the now linear partial differential equations (22)-(23) is obtained by approximating the exact solutions of $f(\eta,\xi)$ and $\theta(\eta,\xi)$ by the Lagrange form of polynomial $F(\eta,\xi)$ and $\Theta(\eta,\xi)$ at the selected collocation points,

$$0=\xi_{0}< \xi_{1}<\xi_{2}<\cdots< \xi_{N_{\xi}}=1. $$

The approximation for $f(\eta,\xi)$ and $\theta(\eta,\xi)$ has the form

$$\begin{aligned}& f(\eta,\xi) \approx\sum_{j=0}^{N_{\xi}}F( \eta,\xi)L_{j}(\xi )=\sum_{j=0}^{N_{\xi}}F_{j}( \eta)L_{j}(\xi), \end{aligned}$$

(33)

$$\begin{aligned}& \theta(\eta,\xi) \approx\sum_{j=0}^{N_{\xi}} \Theta(\eta,\xi )L_{j}(\xi)=\sum_{j=0}^{N_{\xi}} \Theta_{j}(\eta)L_{j}(\xi), \end{aligned}$$

(34)

where $F_{j}(\eta)=F(\eta,\xi)$ and $\Theta_{j}(\eta)=\Theta(\eta ,\xi )$, $L_{j}$ is the characteristic Lagrange cardinal polynomial defined as

$$ L_{j}(\xi)=\prod_{k=0,k\neq j}^{M} \frac{\xi-\xi_{k}}{\xi _{j}-\xi_{k}}, $$

(35)

which obey the Kronecker delta equation

$$ L_{j}(\xi_{k})=\delta_{jk}=\left \{ \begin{array}{l@{\quad}l} 0 & \text{if } j\neq k, \\ 1 & \text{if } j=k. \end{array} \right . $$

(36)

The equations for the solution of $F_{j}(\eta)$ and $\Theta_{j}(\eta)$ are obtained by substituting (33)-(34) into (22)-(23) and letting the equations be satisfied at the points $\xi_{i}$, $i=0,1,2,\ldots,N_{\xi}$. To compute the derivatives of the Lagrange polynomial analytically we transform $\xi\in[0,L_{\xi}]$ to $\zeta\in[-1,1]$, and then we choose Chebyshev-Gauss-Lobatto points $\zeta_{i}=\cos\frac{i\pi}{N_{\xi}}$. After using the linear transformation $\xi= L_{\xi} (\zeta+ 1)/2$, the derivatives of $f'$ with respect to the collocation points $\zeta_{j}$ are computed as

$$ {\biggl.\frac{\partial f'}{\partial\xi}\biggr|_{\xi=\xi_{i}}}=2\sum_{j=0}^{N_{\xi }}F'_{j}( \eta)\frac{dL_{j}}{d\zeta}(\zeta_{i})=\sum_{j=0}^{N_{\xi }} \mathbf{d}_{i,j}F'_{j}(\eta),\quad i=0,1,2, \ldots,N_{\xi}, $$

(37)

where $\mathbf{d}_{i,j}=\frac{dL_{j}}{d\zeta}(\zeta_{i})$ ($i=0,1,\ldots ,N_{\xi}$) are entries of the standard Chebyshev differentiation matrix, $\mathbf{d}=\frac{2}{L_{\xi}}d$. We now apply the collocation $(\eta,\xi_{i})$ in (22)-(23) we obtain

$$\begin{aligned}& \biggl(1+\frac{1}{\beta} \biggr) F'''_{r+1,i}( \eta) + a_{1,r}^{(i)} F''_{r+1,i}( \eta) + a_{2,r}^{(i)} F'_{r+1,i}(\eta) + a_{3,r}^{(i)} F_{r+1,i}(\eta) \\& \quad {}+a_{4,r}^{(i)}\sum_{j=0}^{N_{\xi}} \mathbf {d}_{i,j}F'_{r+1,i}(\eta)+ a_{5,r}^{(i)}\sum_{j=0}^{N_{\xi}} \mathbf {d}_{i,j}F_{r+1,i}(\eta) =a_{6,r}^{(i)}, \end{aligned}$$

(38)

$$\begin{aligned}& \frac{1}{\mathit{Pr}} \biggl(1+\frac{4}{3}K \biggr) \Theta''_{r+1,i}(\eta) + b_{1,r}^{(i)} \Theta'_{r+1,i}(\eta) + b_{2,r}^{(i)}\sum _{j=0}^{N_{\xi }}\mathbf{d}_{i,j} \Theta_{r+1,i}(\eta)=b_{3,r}^{(i)}, \end{aligned}$$

(39)

where $a_{k,r}^{(i)} = a_{k,r}(\eta,\xi_{i})$ ($k=1,2,3,4,5,6$) and $b_{k,r}^{(i)} = b_{k,r}(\eta,\xi_{i})$ ($k=1,2,3$). Since the solution at $\xi=0$ ($\zeta=\zeta_{N_{\xi}}$) is known, we evaluate (38)-(39) for $i=0,1,\ldots,N_{\xi}-1$ and the system becomes

$$\begin{aligned}& \biggl(1+\frac{1}{\beta} \biggr) F'''_{r+1,i} + a_{1,r}^{(i)} F''_{r+1,i} + a_{2,r}^{(i)} F'_{r+1,i} + a_{3,r}^{(i)} F_{r+1,i}+a_{4,r}^{(i)} \sum_{j=0}^{N_{\xi}-1}\mathbf {d}_{i,j}F'_{r+1,i} \\& \quad{}+ a_{5,r}^{(i)}\sum_{j=0}^{N_{\xi}-1} \mathbf {d}_{i,j}F_{r+1,i} =a_{6,r}^{(i)} - a_{4,r}^{(i)}\mathbf{d}_{i,N_{\xi }}F'_{r+1,N_{\xi}} - a_{5,r}^{(i)}\mathbf{d}_{i,N_{\xi }}F_{r+1,N_{\xi }}, \end{aligned}$$

(40)

$$\begin{aligned}& \frac{1}{\mathit{Pr}} \biggl(1+\frac{4}{3}K \biggr) \Theta''_{r+1,i} + b_{1,r}^{(i)} \Theta'_{r+1,i} + b_{2,r}^{(i)}\sum _{j=0}^{N_{\xi }-1}\mathbf{d}_{i,j} \Theta_{r+1,i} \\& \quad =b_{3,r}^{(i)} - b_{2,r}^{(i)} \mathbf {d}_{i,N_{\xi}}\Theta_{r+1,N_{\xi}}. \end{aligned}$$

(41)

For each $\xi_{i}$, (40)-(41) forms a system of linear ordinary differential equations with variable coefficients. In this system we apply the Chebyshev spectral collocation independently in the η direction by choosing $N_{\eta}+1$ Chebyshev-Gauss-Lobatto points $0=\eta_{0}<\eta_{1}<\eta_{2}<\cdots<\eta_{N_{\eta}}=\eta_{e}$, where $\eta_{e}$ is a finite value that is chosen to be adequately large to approximate the conditions at ∞. We now implement the collocation in the interval $[0,\eta_{e}]$ on the η-axis, which is then transformed into the interval $[-1,1]$ using a linear transformation $\eta= \eta_{e}(\tau+1)/2$. The collocation points are chosen as $\tau_{j}=\cos\frac{j\pi}{N_{\eta}}$. The derivatives with respect to η are defined in terms of the Chebyshev differentiation matrix as

$$ {\biggl.\frac{d^{p} F'_{r+1,i}}{d \eta^{p}} \biggr\vert _{\eta=\eta _{j}}}= \biggl( \frac{2}{\eta_{e}} \biggr)^{p}\sum_{k=0}^{N_{\eta }}D^{p}_{j,k}F_{r+1,i}( \tau_{k})= \bigl[\mathbf{D}^{p}\mathbf {F}_{r+1,i} \bigr]_{j}, $$

(42)

where p is the order of the derivative, $\mathbf{D}=\frac{2}{\eta _{e}}D$ ($j,k=0,1,2,\ldots,N_{\eta}$) with D being an $(N_{\eta }+1)\times(N_{\eta}+1)$ Chebyshev derivative matrix, and the vector $\mathbf{F}_{r+1,i}$ is defined as

$$\begin{aligned}& \mathbf{F}_{r+1,i}=\bigl[F_{r+1,1}(\tau_{0}), F_{r+1,1}(\tau_{1}),\ldots ,F_{r+1,1}( \tau_{N_{\eta}})\bigr]^{T}, \end{aligned}$$

(43)

$$\begin{aligned}& \boldsymbol{\Theta}_{r+1,i} = \bigl[\Theta_{r+1,i}( \tau_{0}), \Theta _{r+1,i}(\tau_{1}),\ldots, \Theta_{r+1,i}(\tau_{N_{\eta}})\bigr]^{T}, \end{aligned}$$

(44)

substituting (42) into (40) we get

$$\begin{aligned}& \mathbf{A}^{(i)}\mathbf{F}_{r+1,i} + \mathbf{a}_{4,r}^{(i)}\sum_{j=0}^{N_{\xi}-1}d_{i,j} \mathbf{D}\mathbf{F}_{r+1,j}+ \mathbf {a}_{5,r}^{(i)}\sum _{j=0}^{N_{\xi}-1}d_{i,j} \mathbf{F}_{r+1,j} = R_{1}^{(i)}, \end{aligned}$$

(45)

$$\begin{aligned}& \mathbf{A}^{(i)} = \biggl(1+\frac{1}{\beta}\biggr) \mathbf{D}^{3} + \mathbf {a}_{1,r}^{(i)} \mathbf{D}^{2} + a_{2,r}^{(i)}\mathbf{D} + \mathbf {a}_{2,r}^{(i)}, \end{aligned}$$

(46)

$$\begin{aligned}& \mathbf{R}_{1}^{(i)}=\mathbf{a}_{6,r}^{(i)} - \mathbf {a}_{4,r}^{(i)}d_{i,N_{\xi}}\mathbf{D} \mathbf{F}_{r+1,N_{\xi}} - \mathbf {a}_{5,r}^{(i)}d_{i,N_{\xi}} \mathbf{F}_{r+1,N_{\xi}} \end{aligned}$$

(47)

$$\begin{aligned}& \mathbf{B}^{(i)}\boldsymbol{\Theta}_{r+1,i} + \mathbf {b}_{2,r}^{(i)}\sum_{j=0}^{N_{\xi}-1}d_{i,j} \boldsymbol{\Theta}_{r+1,j} = R_{2}^{(i)}, \end{aligned}$$

(48)

$$\begin{aligned}& \mathbf{B}^{(i)} = \frac{1}{\mathit{Pr}}\biggl(1+\frac{4}{3}K \biggr)\mathbf{D}^{2} + \mathbf{b}_{1,r}^{(i)} \mathbf{D}, \end{aligned}$$

(49)

$$\begin{aligned}& \mathbf{R}_{2}^{(i)}=\mathbf{b}_{3,r}^{(i)} - \mathbf {b}_{2,r}^{(i)}d_{i,N_{\xi}}\boldsymbol{ \Theta}_{r+1,N_{\xi}} . \end{aligned}$$

(50)

$\mathbf{a}_{k,r}$ ($k=1,2,3,4,5,6$), $\mathbf{b}_{k,r}$ ($k=1,2,3$) are the diagonal matrices with vectors $[a_{k,r}(\tau_{0}), a_{k,r}(\tau _{1}), \ldots ,a_{k,r}(\tau_{N_{x}})]^{T}$ and $[b_{k,r}(\tau_{0}), b_{k,r}(\tau_{1}),\ldots,b_{k,r}(\tau _{N_{x}})]^{T}$. We then impose boundary conditions and a matrix system is formed as follows:

$$\begin{aligned}& \mathbf{AA}=\left ( \begin{array}{@{}c|c@{}} \begin{array}{@{}c@{\quad}c@{\quad}c@{\quad}c@{}} A_{0,0} & A_{0,1}&\cdots & A_{0,N_{\xi}-1} \\ A_{1,0} & A_{1,1}&\cdots & A_{1,N_{\xi}-1} \\ \vdots &\vdots& \ddots& \vdots \\ A_{N_{\xi}-1,0}& A_{N_{\xi}-1,1}&\cdots& A_{N_{\xi}-1,N_{\xi}-1} \end{array}& \begin{array}{@{}c@{\quad}c@{\quad}c@{\quad}c@{}} 0 &0&\cdots& 0 \\ 0& 0&\cdots& 0 \\ \vdots &\vdots &\ddots &\vdots \\ \vdots& \vdots& \ddots& \vdots \end{array} \\ \hline \begin{array}{@{}c@{\quad}c@{\quad}c@{\quad}c@{}} 0 &0&\cdots&0 \\ 0 &0&\cdots &0 \\ \vdots &\vdots& \ddots& \vdots \\ 0& 0&\cdots& 0 \end{array} & \begin{array}{@{}c@{\quad}c@{\quad}c@{\quad}c@{}} B_{0,0} & B_{0,1}&\cdots& B_{0,N_{\xi}-1} \\ B_{1,0}& B_{1,1}&\cdots& B_{1,N_{\xi}-1} \\ \vdots& \vdots& \ddots& \vdots \\ B_{N_{\xi}-1,0}& B_{N_{\xi}-1,1}&\cdots& B_{N_{\xi}-1,N_{\xi}-1} \end{array} \end{array} \right ), \\& \mathbf{FF}=\left ( \begin{array}{@{}c@{}} \mathbf{F}_{r+1,0} \\ \mathbf{F}_{r+1,1} \\ \vdots \\ \mathbf{F}_{r+1,N_{\xi}-1} \\ \hline \boldsymbol{\Theta}_{r+1,0} \\ \boldsymbol{\Theta}_{r+1,1} \\ \vdots \\ \boldsymbol{\Theta}_{r+1,N_{\xi}-1} \end{array} \right ),\qquad \mathbf{RR}=\left ( \begin{array}{@{}c@{}} R_{1}^{(0)} \\ R_{1}^{(1)} \\ \vdots \\ R_{1}^{(N_{\xi}-1)} \\ \hline R_{2}^{(0)} \\ R_{2}^{(1)} \\ \vdots \\ R_{2}^{(N_{\xi}-1)} \end{array} \right ). \end{aligned}$$

The above system is then solved as

$$ \mathbf{FF} = \mathbf{AA}^{-1}\mathbf{RR}. $$

(51)

We have

$$\begin{aligned}& \begin{aligned} &A_{i,i}=\mathbf{A}^{(i)} + a_{4,r}^{(i)}d_{i,i} \mathbf{D} + a_{5,r}^{(i)}d_{i,i}, \\ &B_{i,i}=\mathbf{B}^{(i)} + b_{2,r}^{(i)}d_{i,i}, \quad i=0,1,\ldots ,N_{\xi}-1, \end{aligned} \end{aligned}$$

(52)

$$\begin{aligned}& \begin{aligned} &A_{i,j}= a_{4,r}^{(i)}d_{i,j} \mathbf{D} + a_{5,r}^{(i)}d_{i,j},\quad i\neq j, \\ &B_{i,j}= b_{2,r}^{(i)}d_{i,j} ,\quad i\neq j. \end{aligned} \end{aligned}$$

(53)

4 Results and discussion

In Table 1, we compare the results of the bi-variate quasilinearization method (BQLM) to the successive linearization method (SLM) and MATLAB’s ‘bvp4c’ method. The results of the BQLM are in excellent agreement with these two methods to seven decimal places showing the accuracy of this method.

Table 1 Comparison of the values of the skin friction coefficient and heat transfer coefficient obtained by BQLM against the SLM and bvp4c for $\pmb{\beta \rightarrow\infty}$ , $\pmb{\Lambda=\Lambda_{1}=M=\xi =K=f_{w}= S_{f}= S_{T}=0}$

To understand the behavior of the velocity and temperature profiles the illustrations for the numerical solution obtained are depicted as graphs in Figures 2-8.

Figure 2

Effects of Casson parameter on (a) velocity and (b) temperature profiles for different Λ and $\pmb{S_{T}}$ , at $\pmb{\Lambda_{1}=0.5}$ , $\pmb{K=0.5}$ , $\pmb{\mathit{Pr}=1}$ , $\pmb{\mathit{Ec} =1}$ , $\pmb{M=1}$ , $\pmb{f_{w}=0.5}$ , $\pmb{S_{f}=0.5}$ .

Figure 2 shows the influence of Casson parameter β on the velocity and the temperature profiles at different values of the Forchheimer parameter Λ in the case of velocity profiles (see Figure 2(a)), and thermal slip factor $S_{T}$ in the case of temperature profiles (see Figure 2(b)). Increasing the Casson parameter increases the velocity profiles close to the boundary, as $\beta\rightarrow\infty$, $\frac{1}{\beta}\rightarrow0$, the fluid becomes Newtonian. The same result is noted in Ramachandra et al. [1]. The difference in these results from those of Ramachandra et al. [1] is that there is a reverse effect that is noted further away from the boundary. This is due to the presence of the velocity slip factor which tends to assist flow at the boundary. Increasing the Forchheimer parameter reduce velocity profiles, this is caused by the transpiration effect taking place at the surface of the circle. Increasing the Casson parameter β reduces temperature profiles (see [1–3]). In the case of no thermal slip factor, higher temperature profiles are noticed and lower temperature profiles are noticed in the presence of the thermal slip factor, as shown in Figure 2(b).

Figure 3 depicts the effect of the magnetic parameter M on the velocity and temperature profiles at different values of the velocity slip factor $S_{f}$ in the case of velocity profiles (see Figure 3(a)), and thermal slip factor in the case of temperature profiles (see Figure 3(b)). M suppress velocity profiles. A magnetic field incident across a flowing fluid retards the velocity of the flow as noted in [4, 12]. The velocity profiles in the case of no-slip conditions resemble those of Narayana et al. [10]. It is noted that the velocity profiles are increased with the increase in the velocity slip factor as expected on a lubricated surface. In Figure 3(b) the temperature profiles are increased with increasing magnetic parameter and are decreased by the increase in the Forchheimer parameter Λ. Increasing the inertia causes more friction on the surface changing surface temperature, thereby causing the observed thermal slip.

Figure 3

Effects of magnetic parameter on (a) velocity and (b) temperature profiles for different $\pmb{S_{f}}$ and $\pmb{S_{T}}$ , at $\pmb{\Lambda_{1}=0.5}$ , $\pmb{K=0.5}$ , $\pmb{\mathit{Pr}=1}$ , $\pmb{\mathit{Ec} =1}$ , $\pmb{f_{w}=0.5}$ , $\pmb{\Lambda =0.5}$ , $\pmb{\beta=2}$ .

Figure 4 shows the effect of the radiation parameter K on the velocity and temperature profiles at different values of the velocity slip factor $S_{f}$ in the case of velocity profiles (see Figure 4(a)), and thermal slip factor in the case of temperature profiles (see Figure 4(b)). The same observation is noted in Figure 4(b); increasing K increases the temperature profiles.

Figure 4

Effects of radiation parameter on (a) velocity and (b) temperature profiles for different $\pmb{S_{f}}$ and $\pmb{S_{T}}$ , at $\pmb{\beta=2}$ , $\pmb{\Lambda_{1}=0.5}$ , $\pmb{K=0.5}$ , $\pmb{\mathit{Pr}=1}$ , $\pmb{\mathit{Ec} =1}$ , $\pmb{M=1}$ , $\pmb{f_{w}=0.5}$ , $\pmb{\Lambda=0.5}$ .

Figure 5 illustrates the velocity response to the change in the Darcian drag force term $\Lambda_{1}$ at different values of the suction/injection parameter $f_{w}$ (see Figure 5(a)), and the effect of the Eckert number Ec at different values of the suction/injection parameter $f_{w}$ (see Figure 5(b)). In this case the $\Lambda_{1}$ is inversely proportional to the Darcy number, Da (also inversely proportional to the permeability of the medium). Increasing $\Lambda_{1}$ decreases the velocity profiles. It is noted that the case $f_{w}<0$ (blowing) assists the flow, while the case $f_{w}>0$ (suction) retards the velocity profiles. The differences of these results to those of Ramachandra et al. [1] are shown in the presence of a magnetic field which generally tend to suppress velocity profiles. The presence of viscous dissipation is also noted in the increase in the temperature profiles. In Figure 5(b) increasing the Eckert number results in the increase in the temperature profiles.

Figure 5

Effects of $\pmb{\Lambda_{1}}$ and Ec on (a) velocity and (b) temperature profiles for different $\pmb{f_{w}}$ , at $\pmb{\beta =2}$ , $\pmb{\Lambda =0.5}$ , $\pmb{K=0.5}$ , $\pmb{M=1}$ , $\pmb{K=0.5}$ , $\pmb{\mathit{Pr}=1}$ , $\pmb{S_{T}=0.5}$ .

Figure 6 shows the plot of the skin friction coefficient versus ξ for different values of both the velocity slip factor $S_{f}$ (see Figure 6(a)), and the thermal slip factor $S_{T}$ (see Figure 6(b)). The skin friction coefficient decreases with increasing velocity slip factor and also decreases with increasing ξ. In Figure 6(b) the skin friction coefficient decreases with increasing thermal slip factor.

Figure 6

Skin friction coefficient vs. ξ for different $\pmb{S_{f}}$ (a) and skin friction coefficient vs. ξ for different $\pmb{S_{T}}$ (b) at $\pmb{\beta =2}$ , $\pmb{\Lambda _{1}=0.5}$ , $\pmb{K=0.5}$ , $\pmb{M=1}$ , $\pmb{\mathit{Pr}=1}$ , $\pmb{\mathit{Ec} =1}$ , $\pmb{f_{w}=0.5}$ , $\pmb{\Lambda=0.5}$ .

Figure 7 shows the plot of the heat transfer coefficient versus ξ for different values of both the velocity slip factor $S_{f}$ (see Figure 7(a)), and the thermal slip factor $S_{T}$ (see Figure 7(b)). It is noted that, at the surface of the circle, there is a more significant change in heat transfer than further away from the surface of the circle. In Figure 7(b) increasing the velocity slip factor does not change the heat transfer at the surface. It is noticed that there is a large difference in the heat transfer further away from the circle surface.

Figure 7

Heat transfer coefficient vs. ξ for different $\pmb{S_{f}}$ (a) and skin friction coefficient vs. ξ for different $\pmb{S_{T}}$ (b) at $\pmb{\beta =2}$ , $\pmb{\Lambda _{1}=0.5}$ , $\pmb{K=0.5}$ , $\pmb{M=1}$ , $\pmb{\mathit{Pr}=1}$ , $\pmb{\mathit{Ec} =1}$ , $\pmb{f_{w}=0.5}$ , $\pmb{\Lambda=0.5}$ .

Figure 8 shows the plot of the skin friction coefficient versus K for different values of both the velocity slip factor $S_{f}$ (see Figure 8(a)), and the Nusselt number versus thermal slip factor $S_{T}$ (see Figure 8(b)). In Figure 8(a) the skin friction coefficient increases with increasing thermal radiation K and decreases with increasing velocity slip factor, which is expected on a lubricated surface. There is a noticeable difference in the skin friction coefficient at different values of the velocity slip factor. In Figure 8(b) increasing the thermal radiation decreases the heat transfer coefficient. Increasing the thermal slip factor increases the heat transfer coefficient.

Figure 8

Skin friction coefficient vs. K for different $\pmb{S_{f}}$ (a) and skin friction coefficient vs. ξ for different $\pmb{S_{T}}$ (b) at $\pmb{\beta=2}$ , $\pmb{\Lambda _{1}=0.5}$ , $\pmb{M=1}$ , $\pmb{\mathit{Pr}=1}$ , $\pmb{\mathit{Ec} =1}$ , $\pmb{f_{w}=0.5}$ , $\pmb{\Lambda=0.5}$ .

5 Conclusion

The study presented in this analysis of effects of radiation on MHD free convection of a Casson fluid from a horizontal circular cylinder with partial slip in non-Darcy porous medium with viscous dissipation provides numerical solutions for the boundary layer flow and heat transfer. The coupled nonlinear governing partial differential equations were solved using the bi-variate quasilinearization method (BQLM) and validated by the successive linearization method (SLM) and MATLAB’s ‘bvp4c’. This paper also describes the BQLM numerical method which uses collocation methods in both directions η (direction of increasing boundary layer thickness) and ξ (radial transpiration direction). The most important results are those reflected in the presence of a magnetic field and viscous dissipation in a Casson fluid, which were never reported before.