Flow Dynamics of the Homogeneous and Heterogeneous Reactions with an Internal Heat Generation and Thermal Radiation between Two Squeezing Plates

Abstract

This paper explores the time dependent squeezing flow of a viscous fluid between parallel plates with internal heat generation and homogeneous/heterogeneous reactions. The motive of the present effort is to upgrade the heat transformation rate for engineering and industrial purpose with the rate of chemical reaction. For this purpose the equations for the conservation of mass, momentum, energy and homogeneous/heterogeneous reactions are transformed to a system of coupled equations using the similarity transformation. According to HAM, with the proper starting assumptions and other factors, a similarity solution may be found. On the way to verifying the validity and correctness of HAM findings, we compare the HAM solution with numerical solver programme BVP4c to see whether it matches up. The results of a parametric inquiry are summarized and presented with the use of graphs.

1. Introduction

It is easy to modify and regulate the convergence area of an infinite series solution using homotopy analysis by introducing an additional parameter. Tests are conducted on a few test cases to verify the method’s validity, efficiency and applicability. As a consequence of the results, the to be highly efficient, clear, and simple. HAM has several advantages over many traditional analytical methods. It is a series expansion method that does not straightforwardly rely upon small or large physical parameters. Accordingly, it is pertinent for weakly as well as powerfully nonlinear issues. It gives fantastic adaptability to the declaration of the solution and how the solution is specifically gotten, and gives incredible opportunity in picking the base elements of the ideal solution the corresponding auxiliary auxiliary linear operator of homotopy. HAM provides a simple way to ensure the convergence of the solution series and is able to combine with other techniques employed in nonlinear differential equations such as spectral methods, Padé approximates, and so on. Squeezing flow of the viscous fluid between parallel plates is one of the most important type of the flow. The importance of the squeezing flow could be observed in many industries such as injection, moulding liquid-metal, lubricated bearings, compression, squeezed films in power transmission, cooling and heat water. Heterogeneous and homogeneous reactions in chemically reactive systems are known. In food processing, ceramics production, and crop loss owing to freezing and wetness over agricultural fields and orchards, the interplay between these two reactions have uses. The researchers are taking great interest because of its vast choices of useful application for the investigation of the velocity field and heat transfer between squeezing plates. Stefan [1] did the initial work on the squeezing flow and investigated lubrication system in their research. Stefan opened a new field of research for the scientists on squeezing flow. In recent years the study of the squeezing flow has been pushed forward by different researchers. Ran [2] investigated the quasi-steady axisymmetric Newtonian fluid squeezed between two parallel plates and found an analytical solution by homotopy analysis method (HAM). Mustafa [3] studied an unsteady viscous and squeezing flow between squeezing plates with heat and mass transfer. Sheikholislami [4] studied a nanofluid of a squeezed duct by neural network. He investigated this problem for the effect of growth ratio, nanopartical absorption, power law indexed and Reynolds number.

Hayat [5] explored the effects of the chemical reaction and convective conditions on squeezing flow and mostly awareness has been focused on the coverlet flow. Ahmad [6] has recently examined the impact on the transportation features of the squeezed fluid on velocity, thermal and solvent effects on the transport properties. The fluid velocity near both plates decays with larger Darcy number, the flow have obtained more resistance. Hayat [7] showed in a theoretical analysis in a Darcy–Forchheimer-permeable medium with thermal radiation in the gripping flow of water-based carbon nanotubes. In another paper, Ahmad [8] has proposed to investigate in theoretical terms the combination of thermal-radiation effect and chemical reaction in squeezing fluency of the Sutterby fluid in the squeezed tube. It was recently studied by Cortell [9] Nonlinear thermal radiation with the effect of heat transfer on the flow of an extended sheet and increased Prandtl number, $P_r$, lead to thermal boundary layer thickness decrease. Recently Shehzad [10] have been addressed the nonlinear thermal effects of Jeffrey nanofluid on the 3D flow and obtained the increased functions of the radiation parametric quantity of temperature and nanoparticle absorption. Mustafa [11] study the effect of heat generation/absorption and radiation on the movement of the second-grade nanofluid on a three-dimensional boundary layer and its transformation into a fluid operating in nanoparticles and the combination of solar incident radiation. More natural convection of nano-fluids is important for next-generation nuclear reactor equipment, heat exchangers, solar film collectors, etc. Hayat [12] has investigated a chemical reaction between homogeneous and heterogeneous reactions with nano-fluid and obtained that the strength of the heterogeneous parameter increases the profiles of attention. Ravikiran and Radhakrishnamacharya [13] investigated the effect of homogeneous and heterogeneous effects on peristaltic flow of the wall through the porous medium has shown the effective dispersion coefficient to be increased with Darcy and the slippy parameter in both homogenates, mixed and homogeneous chemical reactions. Hayat [14] discussed Homogeneous and heterogeneous chemical reactions with effects of heat transfer on Carreau and hall stream flows and showed the similar deeds are experimental on temperature profile and heat transfer coefficient in submissive wall parameters.

Sheikholeslami [15] discussed a nano-fluid flow between two magnetic field impact parallel panels of time-dependent and achieved at the surface decreases with increasing values of the strong points of homogeneous and heterogeneous chemical reaction parameters. Sheikholeslami [16] studied the effect of radiation on the movement of nanofluids using a two-stage model and MHD. The recession by the radiation parameter would be the boundary layer thickness of the temperature. Hayat [17] studied the effects of MHD on the bidirectional flow of a homogenous and heterogeneous nanofluid, considering slip velocity and velocity profiles in the second-order, decreases the functions of the slip velocity and volume fraction of nanoparticles. It was suggested by Merkin [18], a mathematical pattern for homogenous and heterogeneous isothermal reactions with boundary layer flow and has been shown that the primary production of the surface reaction is dominated only by the results of the homogeneous reaction as the flow from the leading surface progresses. Raju [19] discussed the chemical reactions and impacts of radiation on thermophoretic MHD fluidity over the associated isothermal surface and numerical observations show that the magnetic field surface parameter velocity profiles are more effectively affecting. Toghraie [20] studied the numerical thermal analysis of water’s boiling heat transfer based on a turbulent jet impingement on heated surface. The researchers in [21] investigated the Impact of variable fluid properties on forced convection of Fe3O4/CNT/water hybrid nanofluid in a double-pipe mini-channel heat ex-changer.

According to the current survey of the literature review, the earlier researchers have not tried to examine the study of Homogeneous and Heterogeneous reactions with an internal heat generation/absorption and thermal radiation between squeezing plates. An analysis is made on the impact of the homogeneous and heterogeneous chemical reaction on the flow between the parallel squeezing plates. In this paper, the effects of various parameters are discussed on velocity field, heat distribution, homogeneous and heterogeneous reactions, coefficient of skin friction, heat and mass fluxes.

2. Mathematical Formulation

Consider a laminar, axisymmetric, incompressible and viscous fluid between horizontally parallel and squeezing plates. The two plates are separated by a distance $h\left(t\right)=l\sqrt{1-\alpha t}$, where l is the representative length equivalent to the plates separation at $t=0$ as shown in Figure 1. For $\alpha >0$, the two plates are squeezed until they touch at $t=\frac{1}{\alpha }$ and for $\alpha <0$, the two plates are separated. The velocity field is also under the influence of a uniform magnetic field $B\left(t\right)=\frac{B_0}{\sqrt{1-\alpha t}}$ applied in the $y-direction$. The upper and lower plates are kept at constant temperature $T_u$ and $T_l$, respectively. An analysis of heat transfer is made in the presence of internal heat generation/absorption and thermal radiation. The upper and lower plates are maintained at a constant temperature $T_u$ and $T_l$, respectively. The Chaudhary [22] model is adopted for the investigation of homogeneous and heterogeneous reactions. The homogeneous reaction in cubic autocatalysis looks like this:

$$A1+2A2\to 3A2,$$

(1)

and the rate is $K_{c}a{b}^2$, whereas on the catalyst surface, heterogeneous reaction is

$$A1\to A2,$$

(2)

The concentration of the chemical species $A_1$ and $A_2$ are presented by a and b, respectively, and the rate constants are denoted by $k_c$. The reaction rate vanishes beyond the boundary layer edge in external flow, as shown by the equations above. The key time-dependent conservation of mass, momentum, thermal energy, homogeneous and heterogeneous reactions in the Cartesian coordinates $(x,y,z)$ taken at the center of the lower plate are [23,24,25].

Conservation of Mass Equation:

$$\frac{\partial u^{\ast }}{\partial x}+\frac{\partial v^{\ast }}{\partial y}=0$$

(3)

Conservation of Momentum Equation:

$$\frac{\partial u^{\ast }}{\partial t}+u^{\ast }\frac{\partial u^{\ast }}{\partial x}+v^{\ast }\frac{\partial u^{\ast }}{\partial y}=-\frac{1}{\rho }\frac{\partial p}{\partial x}+\frac{\mu }{\rho }\left(\frac{{\partial }^2{u}^{\ast }}{\partial x^2}+\frac{{\partial }^2{u}^{\ast }}{\partial y^2}\right)-\sigma {\varrho }^2{u}^{\ast }$$

(4)

$$\frac{\partial v^{\ast }}{\partial t}+u^{\ast }\frac{\partial v^{\ast }}{\partial x}+v^{\ast }\frac{\partial v^{\ast }}{\partial y}=-\frac{1}{\rho }\frac{\partial p}{\partial y}+\frac{\mu }{\rho }\left(\frac{{\partial }^2{v}^{\ast }}{\partial x^2}+\frac{{\partial }^2{v}^{\ast }}{\partial y^2}\right)$$

(5)

Conservation of Energy Equation:

$$\frac{\partial T}{\partial t}+u^{\ast }\frac{\partial T}{\partial x}+v^{\ast }\frac{\partial T}{\partial y}=\frac{k}{\left(\rho c_p\right)}\left(\frac{{\partial }^{2}T}{\partial x^2}+\frac{{\partial }^{2}T}{\partial y^2}\right)+\frac{1}{\left(\rho c_p\right)}\frac{16{\sigma }^{\ast }{T}_0^3}{3{\kappa }^{\ast }}\frac{{\partial }^{2}T}{\partial y^2}+\frac{Q^{\ast }}{\left(\rho c_p\right)}\left(T-T_0\right)$$

(6)

Homogeneous and Heterogeneous Equations:

$$\frac{\partial a}{\partial t}+u^{\ast }\frac{\partial a}{\partial x}+v^{\ast }\frac{\partial a}{\partial y}=D_A\frac{{\partial }^{2}a}{\partial y^2}-K_{c}a{b}^2$$

(7)

$$\frac{\partial b}{\partial t}+u^{\ast }\frac{\partial b}{\partial x}+v^{\ast }\frac{\partial b}{\partial y}=D_B\frac{{\partial }^{2}b}{\partial y^2}+K_{c}a{b}^2$$

(8)

In the following equations, $u^{\ast }$, $v^{\ast }$ represent the horizontal and vertical velocity components, respectively, where T is the temperature distribution, p is the fluid pressure and a, b are the homogeneous and heterogeneous reaction variables, $\rho$ is the fluid density, $\left(\rho c_p\right)$ is the heat capacity and $\sigma$ is electrical conductivity of the fluid, $D_A$, $D_B$ display chemical species of the diffusion coefficients A and B, respectively, $K^{\ast }$ is the permeability, k is the thermal conductivity and $Q^{\ast }$ is the heat generation.

3. Boundary Conditions

The boundary conditions are defined as:

$$\begin{array}{cc}\hfill & u^{\ast }=0,v^{\ast }=0,T=T_l,D_A\frac{\partial a}{\partial y}=k_{3}a,D_B\frac{\partial b}{\partial y}=-k_{3}aaty=0\hfill \\ & u^{\ast }=0,v^{\ast }=\frac{-\alpha D}{2\sqrt{1-\alpha t}},T=T_{u^{\ast }},a=a_{0}b=0aty=h\left(t\right)\hfill \end{array}$$

(9)

The differential equations system was converted into the system of ordinary differential equations using the following similarity transformations [26].

$$u^{\ast }=\frac{\alpha x{f}^{\prime }\left(\eta \right)}{2(1-\alpha t)},v^{\ast }=\frac{-\alpha lf\left(\eta \right)}{(1-\alpha t)},T=\theta \left(\eta \right)T_l,\eta =\frac{y}{l\sqrt{1-\alpha t}}$$

and

$$a=a_{0}H\left(\eta \right),b=a_{0}G\left(\eta \right),B\left(t\right)=\frac{B_0}{\sqrt{1-at}},\theta =\frac{T-T_u}{T_l-T_u},$$

Here f, $f^{\prime }$, $\theta$, G and H are the velocity components, fluid temperature, homogeneous and heterogeneous reactions, respectively.

The continuity Equation (3) is identically satisfied and the momentum, energy, homogeneous and heterogeneous equations takes the following form,

$$f^{⁗}-S\left(f^{″}{f}^{\prime }+2f^{″}-f{f}^{‴}+\eta f^{‴}\right)-H{a}^{2}L{f}^{″}=0$$

(10)

$${\theta }^{″}\left(1+\frac{4}{3}R\right)+P_{r}S\left(f{\theta }^{\prime }-2\eta \theta +Q\theta \right)=0$$

(11)

$$H^{″}-S_c{K}_{1}H{G}^2-S_{c}S\left(\eta H^{\prime }-fH\right)=0$$

(12)

$$G^{″}\delta -S_c{K}_{1}H{G}^2-S_{c}S\left(\eta G^{\prime }-f{G}^{\prime }\right)=0.$$

(13)

The boundary conditions are transformed as:

$$\begin{array}{cc}\hfill & f\left(0\right)=0,f^{\prime }\left(0\right)=0,\theta \left(0\right)=1,H^{\prime }\left(0\right)=K_{2}H\left(0\right),\delta G^{\prime }\left(0\right)=-K_{2}H\left(0\right)\hfill \\ & f\left(1\right)=0.5,f^{\prime }\left(1\right)=0,\theta \left(1\right)=0,H\left(1\right)=1G\left(1\right)=0,\hfill \end{array}$$

(14)

Here $S=\frac{\alpha l^2}{2v_f}$ is the squeezed Reynolds number, $Ha=l{B}_0\sqrt{\frac{{\sigma }_f}{{\mu }_f}}$, is the Hartman number, $L=\frac{\alpha {\sigma }_{nf}}{x}$ is the length parameter, $R=\frac{4{\sigma }^{\ast }{T}_0^3}{K{K}^{\ast }}$ is the radiation parameter, $P_r=\frac{\mu c_p}{K}$ is the Prandtl number, $Q=\frac{2Q^{\ast }(1-\alpha t)}{\rho c_p\alpha }$ is the heat generation parameter, $S_c=\frac{v_f}{D_A}$ is the Schmidt number, $K_1=\frac{8k_c{a}_0^2(1-\alpha t)}{{\alpha }^3}$ is the homogeneous reaction strength, $K_2=\frac{k_s}{D_A}$ is the heterogeneous reaction strength, $\delta =\frac{D_A}{D_B}$ is the ratio of the diffusion coefficients.

Here, it is consider that A and B diffusion coefficients of chemical species to be comparable size. The other hypothesis is, $D_A$ and $D_B$ are equivalent, so $\delta =\frac{D_A}{D_B}=1$, also $G\left(\eta \right)+H\left(\eta \right)=1$ [18],

Coefficients of Interest

The local Nusselt number ($Nu$), Sherwood number ($Sh$) and skin-friction coefficient ($Cf$), are some of the coefficient of interest in engineering. $C_f^{\ast }=\frac{S}{h}{C}_f=f^{″}\left(0\right)$, $-{\theta }^{\prime }\left(0\right)=Nu$, $-G^{\prime }\left(0\right)=-H^{\prime }\left(0\right)=Sh$,

4. Approximate Analytical Solution

To solve system of Equations (10)–(13) is used the analytic method HAM. Due to HAM, the functions $f\left(\varrho \right)$, $\theta \left(\varrho \right)$, $H\left(\varrho \right)$ and $G\left(\varrho \right)$ can be stated, by a set of base functions ${\varrho }^c,c\ge 0$ as “Rashidi (2015)”[27]:

$$f_{\mathsf{\Psi }}\left(\varrho \right)=\sum _{\xi =0}^{\infty }{a}_{\xi }{\varrho }^{\xi },$$

(15)

$${\theta }_{\mathsf{\Psi }}\left(\varrho \right)=\sum _{\xi =0}^{\infty }{b}_{\xi }{\varrho }^{\xi },$$

(16)

$$H_{\mathsf{\Psi }}\left(\varrho \right)=\sum _{\xi =0}^{\infty }{c}_{\xi }{\varrho }^{\xi },$$

(17)

$$G_{\mathsf{\Psi }}\left(\varrho \right)=\sum _{\xi =0}^{\infty }{d}_{\xi }{\varrho }^{\xi },$$

(18)

where $a_{\xi }$, $b_{\xi }$, $c_{\xi }$ and $d_{\xi }$ are the constant coefficients to be determined. Initial approximations are chosen follows:

$$f_0\left(\varrho \right)=\frac{3}{2}{\varrho }^2-{\varrho }^3,$$

(19)

$${\theta }_0\left(\varrho \right)=(S_r-1)\varrho +1,$$

(20)

$$H_0\left(\varrho \right)=1-\varrho ,$$

(21)

$$G_0\left(\varrho \right)=1-\varrho .$$

(22)

The auxiliary operators are chosen as

$${\ell }_f=\frac{{\partial }^4}{\partial {\varrho }^4},{\ell }_{\theta }=\frac{{\partial }^2}{\partial {\varrho }^2},{\ell }_H=\frac{{\partial }^2}{\partial {\varrho }^2},{\ell }_G=\frac{{\partial }^2}{\partial {\varrho }^2},$$

(23)

with the following properties

$${\ell }_f({\xi }_1{\varrho }^3+{\xi }_2{\varrho }^2+{\xi }_3\varrho +{\xi }_4)=0,$$

(24)

$${\ell }_{\theta }({\xi }_5\varrho +{\xi }_6)=0,$$

(25)

$${\ell }_H({\xi }_7\varrho +{\xi }_8)=0,$$

(26)

$${\ell }_G({\xi }_9\varrho +{\xi }_{10})=0,$$

(27)

where ${\xi }_1$, ${\xi }_2$, ${\xi }_3$, ${\xi }_4$, ${\xi }_5$, ${\xi }_6$, ${\xi }_7$, ${\xi }_8$, ${\xi }_9$ and ${\xi }_{10}$ are arbitrary constants.

The $zeroth$ order deformation problems can be obtained as:

$$(1;\varpi ){\ell }_f[\overline{f}(\varrho ;\varpi )-f_0\left(\varrho \right)]=q{\hslash }_f{\mathit{N}}_f[\overline{f}(\varrho ;\varpi ),\overline{\theta }(\varrho ;\varpi ),\overline{m}(\varrho ;\varpi ),\overline{n}(\varrho ;\varpi )],$$

(28)

$$(1;\varpi ){\ell }_{\theta }[\overline{\theta }(\varrho ;\varpi )-{\theta }_0\left(\varrho \right)]=q{\hslash }_{\theta }{\mathit{N}}_{\theta }[\overline{f}(\varrho ;\varpi ),\overline{\theta }(\varrho ;\varpi ),\overline{m}(\varrho ;\varpi ),\overline{n}(\varrho ;\varpi )],$$

(29)

$$(1;\varpi ){\ell }_H[\overline{H}(\varrho ;\varpi )-H_0\left(\varrho \right)]=q{\hslash }_H{\mathit{N}}_H[\overline{f}(\varrho ;\varpi ),\overline{H}(\varrho ;\varpi ),\overline{G}(\varrho ;\varpi )],$$

(30)

$$(1;\varpi ){\ell }_G[\overline{G}(\varrho ;\varpi )-G_0\left(\varrho \right)]=q{\hslash }_G{\mathit{N}}_G[\overline{f}(\varrho ;\varpi ),\overline{H}(\varrho ;q),\overline{G}(\varrho ;\varpi )].$$

(31)

The nonlinear operators of Equations (20)–(23) are defined as

$$\begin{array}{cc}\hfill {\mathit{N}}_f[\overline{f}(\varrho ;\varpi ),\overline{\theta }(\varrho ;\varpi )]=\frac{{\partial }^4\overline{f}(\varrho ;\varpi )}{\partial {\varrho }^4}& -S[\varrho \frac{{\partial }^3\overline{f}(\varrho ;\varpi )}{\partial {\varrho }^3}+2\frac{{\partial }^2\overline{f}(\varrho ;\varpi )}{\partial {\varrho }^2}\hfill \\ & +\frac{{\partial }^2\overline{f}(\varrho ;\varpi )}{\partial {\varrho }^2}\frac{\partial \overline{f}(\varrho ;\varpi )}{\partial \varrho }-f\frac{{\partial }^3\overline{f}(\varrho ;\varpi )}{\partial {\varrho }^3}]-H{a}^{2}L\frac{{\partial }^2\overline{f}(\varrho ;\varpi )}{\partial {\varrho }^2}\hfill \end{array}$$

(32)

$$\begin{array}{cc}\hfill {\mathit{N}}_{\theta }[\overline{f}(\varrho ;\varpi ),\overline{\theta }(\varrho ;\varpi )]=\frac{{\partial }^2\overline{\theta }(\varrho ;\varpi )}{\partial {\varrho }^2}(1+4/3R)& +PrS\left[f\frac{\partial \overline{\theta }(\varrho ;\varpi )}{\partial \varrho }-2\varrho \overline{\theta }(\varrho ;\varpi )+Q\overline{\theta }(\varrho ;\varpi )\right]\hfill \end{array}$$

(33)

$${\mathit{N}}_H[\overline{f}(\varrho ;\varpi ),\overline{H}(\varrho ;\varpi ),\overline{G}(\varrho ;\varpi )]=\frac{{\partial }^2\overline{H}(\varrho ;\varpi )}{\partial {\varrho }^2}$$

$$-M{S}_c{K}_1\overline{H}(\varrho ;\varpi )\frac{{\partial }^2\overline{G}(\varrho ;\varpi )}{\partial {\varrho }^2}-M{S}_{c}S\left[\varrho \frac{\partial \overline{H}(\varrho ;\varpi )}{\partial \varrho }-f\overline{H}(\varrho ;\varpi )\right],$$

(34)

$${\mathit{N}}_G[\overline{f}(\varrho ;\varpi ),\overline{H}(\varrho ;\varpi ),\overline{G}(\varrho ;\varpi )]=\frac{{\partial }^2\overline{G}(\varrho ;\varpi )}{\partial {\varrho }^2}\delta$$

$$-S_c{K}_1\overline{H}(\varrho ;\varpi )\frac{{\partial }^2\overline{G}(\varrho ;\varpi )}{\partial {\varrho }^2}-S_{c}S\left[\varrho \frac{\partial \overline{G}(\varrho ;\varpi )}{\partial \varrho }-f\frac{\partial \overline{G}(\varrho ;\varpi )}{\partial \varrho }\right],$$

(35)

where $\varpi$ is an embedding parameter, ${\hslash }_f$, ${\hslash }_{\theta }$, ${\hslash }_H$ and ${\hslash }_G$ are the nonzero auxiliary parameter and ${\mathit{N}}_f$, ${\mathit{N}}_{\theta }$, ${\mathit{N}}_H$ and ${\mathit{N}}_G$ are the nonlinear parameters.

For $\varpi =0$ and 1, we have

$$\begin{array}{cc}\hfill & \overline{f}(\varrho ,0)=f_o\left(\varrho \right),\overline{f}(\varrho ,1)=f\left(\varrho \right),\hfill \\ & \overline{\theta }(\varrho ,0)={\theta }_o\left(\varrho \right),\overline{\theta }(\varrho ,1)=\theta \left(\varrho \right),\hfill \\ & \overline{H}(\varrho ,0)=H_o\left(\varrho \right),\overline{H}(\varrho ,1)=H\left(\varrho \right),\hfill \\ & \overline{G}(\varrho ,0)=G_o\left(\varrho \right),\overline{G}(\varrho ,1)=G\left(\varrho \right),\hfill \end{array}$$

(36)

so we can say that as $\varpi$ varies from 0 to 1, $\overline{f}(\varrho ,0)$, $\overline{\theta }(\varrho ,0)$, $\overline{H}(\varrho ,0)$, $\overline{G}(\varrho ,0)$ varies from initial guesses $f_0\left(\varrho \right)$, ${\theta }_0\left(\varrho \right)$, $H_0\left(\varrho \right)$ and $G_0\left(\varrho \right)$ to exact solution $f\left(\varrho \right)$, $\theta \left(\varrho \right)$, $H\left(\varrho \right)$ and $G\left(\varrho \right)$, respectively.

Taylor’s series expansion of these functions yields:

$$f(\varrho ;\varpi )=f_0\left(\varrho \right)+\sum _{\mathsf{\Psi }=1}^{\infty }{\varpi }^{\mathsf{\Psi }}{f}_{\mathsf{\Psi }}\left(\varrho \right),$$

(37)

$$\theta (\varrho ;\varpi )={\theta }_0\left(\varrho \right)+\sum _{\mathsf{\Psi }=1}^{\infty }{\varpi }^{\mathsf{\Psi }}{\theta }_{\mathsf{\Psi }}\left(\varrho \right),$$

(38)

$$H(\varrho ;\varpi )=H_0\left(\varrho \right)+\sum _{\mathsf{\Psi }=1}^{\infty }{\varpi }^{\mathsf{\Psi }}{H}_{\mathsf{\Psi }}\left(\varrho \right),$$

(39)

$$G(\varrho ;\varpi )=G_0\left(\varrho \right)+\sum _{\mathsf{\Psi }=1}^{\infty }{\varpi }^{\mathsf{\Psi }}{G}_{\mathsf{\Psi }}\left(\varrho \right),$$

(40)

$$\begin{array}{cc}\hfill f_{\mathsf{\Psi }}\left(\varrho \right)=\frac{1}{\mathsf{\Psi }!}\frac{{\partial }^{\mathsf{\Psi }}f(\varrho ;\varpi )}{\partial {\varrho }^{\mathsf{\Psi }}}{|}_{\varpi =0},{\theta }_{\mathsf{\Psi }}\left(\varrho \right)=\frac{1}{\mathsf{\Psi }!}\frac{{\partial }^{\mathsf{\Psi }}\theta (\varrho ;\varpi )}{\partial {\varrho }^{\mathsf{\Psi }}}{|}_{\varpi =0},& H_{\mathsf{\Psi }}\left(\varrho \right)=\frac{1}{\mathsf{\Psi }!}\frac{{\partial }^{\mathsf{\Psi }}H(\varrho ;\varpi )}{\partial {\varrho }^{\mathsf{\Psi }}}{|}_{\varpi =0},\hfill \\ & G_{\mathsf{\Psi }}\left(\varrho \right)=\frac{1}{\mathsf{\Psi }!}\frac{{\partial }^{\mathsf{\Psi }}G(\varrho ;\varpi )}{\partial {\varrho }^{\mathsf{\Psi }}}{|}_{\varpi =0}\hfill \end{array}$$

(41)

it should be noticed that the convergence of above series strongly depends upon ${\hslash }_f$, ${\hslash }_{\theta }$, ${\hslash }_H$ and ${\hslash }_G$.

Assuming that these nonzero auxiliary parameters are chosen so that Equations (38)–(41) converges at $\varpi =1$. Therefore one can obtain

$$f\left(\varrho \right)=f_0\left(\varrho \right)+\sum _{\mathsf{\Psi }=1}^{\infty }{f}_{\mathsf{\Psi }}\left(\varrho \right),$$

(42)

$$\theta \left(\varrho \right)={\theta }_0\left(\varrho \right)+\sum _{\mathsf{\Psi }=1}^{\infty }{\theta }_{\mathsf{\Psi }}\left(\varrho \right),$$

(43)

$$H\left(\varrho \right)=H_0\left(\varrho \right)+\sum _{\mathsf{\Psi }=1}^{\infty }{H}_{\mathsf{\Psi }}\left(\varrho \right),$$

(44)

$$G\left(\varrho \right)=G_0\left(\varrho \right)+\sum _{\mathsf{\Psi }=1}^{\infty }{G}_{\mathsf{\Psi }}\left(\varrho \right),$$

(45)

Differentiating the deformation Equations (28)–(31) $\mathsf{\Psi }-times$ with respect to $\varpi$ and putting $\varpi =0$, we have

$${\ell }_f[f_{\mathsf{\Psi }}\left(\varrho \right)-{\chi }_{\mathsf{\Psi }}{f}_{\mathsf{\Psi }-1}\left(\varrho \right)]={\hslash }_f{R}_{f,\mathsf{\Psi }}\left(\varrho \right),$$

(46)

$${\ell }_{\theta }[{\theta }_{\mathsf{\Psi }}\left(\varrho \right)-{\chi }_{\mathsf{\Psi }}{\theta }_{\mathsf{\Psi }-1}\left(\varrho \right)]={\hslash }_{\theta }{R}_{g,\mathsf{\Psi }}\left(\varrho \right),$$

(47)

$${\ell }_H[H_{\mathsf{\Psi }}\left(\varrho \right)-{\chi }_{\mathsf{\Psi }}{H}_{\mathsf{\Psi }-1}\left(\varrho \right)]={\hslash }_H{R}_{H,\mathsf{\Psi }}\left(\varrho \right),$$

(48)

$${\ell }_G[G_{\mathsf{\Psi }}\left(\varrho \right)-{\chi }_{\mathsf{\Psi }}{G}_{\mathsf{\Psi }-1}\left(\varrho \right)]={\hslash }_{\varphi }{R}_{G,\mathsf{\Psi }}\left(\varrho \right),$$

(49)

subject to the boundary conditions

$$\begin{array}{cc}\hfill & f_{\mathsf{\Psi }}\left(0\right)=0,f_{\mathsf{\Psi }}^{\prime }\left(0\right)=0,{\theta }_{\mathsf{\Psi }}\left(0\right)=1,H_{\mathsf{\Psi }}^{\prime }\left(0\right)=k_{2}H\left(0\right),\delta {G^{\prime }}_{\mathsf{\Psi }}\left(0\right)=-k_{2}H\left(0\right),\hfill \\ & f_{\mathsf{\Psi }}\left(1\right)=0.5,f_{\mathsf{\Psi }}^{\prime }\left(1\right)=0,{\theta }_{\mathsf{\Psi }}\left(1\right)=0,H_{\mathsf{\Psi }}\left(1\right)=1,G_{\mathsf{\Psi }}^{\prime }\left(1\right)=0,\hfill \end{array}$$

(50)

where

$$R_{f,\mathsf{\Psi }}\left(\varrho \right)=f_{\mathsf{\Psi }-1}^{⁗}\left(\varrho \right)-S[\varrho f_{\mathsf{\Psi }-1}^{‴}\left(\varrho \right)+2f_{\mathsf{\Psi }-1}^{″}\left(\varrho \right)+2\sum _{j=0}^{\mathsf{\Psi }-1}{f}_j\left(\varrho \right)f_{\mathsf{\Psi }-j-1}^{″}\left(\varrho \right)$$

$$-\sum _{j=0}^{\mathsf{\Psi }-1}{f}_j\left(\varrho \right)f_{\mathsf{\Psi }-j-1}^{‴}\left(\varrho \right)]-H{a}^{2}L{f}_{\mathsf{\Psi }-1}^{″}\left(\varrho \right)$$

(51)

$$R_{\theta ,\mathsf{\Psi }}\left(\varrho \right)={\theta }_{\mathsf{\Psi }-1}^{″}\left(\varrho \right)(1+4/3R)+PrS\left[\sum _{j=0}^{\mathsf{\Psi }-1}{f}_j\left(\varrho \right){\theta }_{\mathsf{\Psi }-j-1}^{\prime }\left(\varrho \right)-2\varrho {\theta }_{\mathsf{\Psi }-1}\left(\varrho \right)+Q{\theta }_{\mathsf{\Psi }-1}\left(\varrho \right)\right]$$

(52)

$$R_{H,\mathsf{\Psi }}\left(\varrho \right)=H_{\mathsf{\Psi }-1}^{″}\left(\varrho \right)+M{S}_c{K}_1\sum _{j=0}^{\mathsf{\Psi }-1}{H}_j\left(\varrho \right)G_{\mathsf{\Psi }-j-1}^{″}\left(\varrho \right)$$

$$-M{S}_{c}S\left[\varrho H_{\mathsf{\Psi }-1}^{\prime }\left(\varrho \right)-\sum _{j=0}^{\mathsf{\Psi }-1}{f}_j\left(\varrho \right)H_{\mathsf{\Psi }-j-1}\left(\varrho \right)\right]$$

(53)

$$R_{G,\mathsf{\Psi }}\left(\varrho \right)=G_{\mathsf{\Psi }-1}^{″}\left(\varrho \right)\delta -S_c{K}_1\sum _{j=0}^{\mathsf{\Psi }-1}{H}_j\left(\varrho \right)G_{\mathsf{\Psi }-j-1}^2\left(\varrho \right)-S_{c}S\left[\varrho G_{\mathsf{\Psi }-1}^{\prime }\left(\varrho \right)-\sum _{j=0}^{\mathsf{\Psi }-1}{f}_j\left(\varrho \right)G_{\mathsf{\Psi }-j-1}^{\prime }\left(\varrho \right)\right]$$

(54)

and ${\chi }_{\mathsf{\Psi }}=\{1,if\mathsf{\Psi }>1,and0,\mathsf{\Psi }=1$.

Finally, the general solution of Equations (46)–(49) can be written as

$$f_{\mathsf{\Psi }}\left(\varrho \right)={\int }_0^{\varrho }{\int }_0^{\varrho }{\int }_0^{\varrho }{\int }_0^{\varrho }{\hslash }_f{R}_{f,\mathsf{\Psi }}\left(z\right)dzdzdzdzdz+{\chi }_{\mathsf{\Psi }}{f}_{\mathsf{\Psi }-1}+{\xi }_1{\varrho }^3+{\xi }_2{\varrho }^2+{\xi }_3\varrho +{\xi }_4,$$

(55)

$${\theta }_{\mathsf{\Psi }}\left(\varrho \right)={\int }_0^{\varrho }{\int }_0^{\varrho }{\hslash }_{\theta }{R}_{\theta ,\mathsf{\Psi }}\left(z\right)dzdz+{\chi }_{\mathsf{\Psi }}{\theta }_{\mathsf{\Psi }-1}+{\xi }_5\varrho +{\xi }_6,$$

(56)

$$H_{\mathsf{\Psi }}\left(\varrho \right)={\int }_0^{\varrho }{\int }_0^{\varrho }{\hslash }_H{R}_{H,\mathsf{\Psi }}\left(z\right)dzdz+{\chi }_{\mathsf{\Psi }}{H}_{\mathsf{\Psi }-1}+{\xi }_{11}\varrho +{\xi }_{12},$$

(57)

$$G_{\mathsf{\Psi }}\left(\varrho \right)={\int }_0^{\varrho }{\int }_0^{\varrho }{\hslash }_G{R}_{G,\mathsf{\Psi }}\left(z\right)dzdz+{\chi }_{\mathsf{\Psi }}{G}_{\mathsf{\Psi }-1}+{\xi }_{13}\varrho +{\xi }_{14},$$

(58)

and so the exact solution $f\left(\varrho \right)$, $\theta \left(\varrho \right)$, $H\left(\varrho \right)$ and $H\left(\varrho \right)$ becomes

$$\begin{array}{cc}\hfill & f\left(\varrho \right)\approx \sum _{n=0}^{\mathsf{\Psi }}{f}_n\left(\varrho \right),\theta \left(\varrho \right)\approx \sum _{n=0}^{\mathsf{\Psi }}{\theta }_n\left(\varrho \right),H\left(\varrho \right)\approx \sum _{n=0}^{\mathsf{\Psi }}{H}_n\left(\varrho \right),G\left(\varrho \right)\approx \sum _{n=0}^{\mathsf{\Psi }}{G}_n\left(\varrho \right).\hfill \end{array}$$

(59)

Optimal Convergence Control Parameters

It is important to note that the series solutions (56)–(59) include ${\hslash }_f$, ${\hslash }_{\theta }$, ${\hslash }_H$ and ${\hslash }_G$ non-zero auxiliary parameters that define the convergence area as well as rate of the homotopy series solutions. The so-called average residual error identified by the “Liao (2010)” to get optimal values of ${\hslash }_f$, ${\hslash }_{\theta }$, ${\hslash }_H$ and ${\hslash }_G$ was used as.

$${\epsilon }_{\mathsf{\Psi }}^f=\frac{1}{\xi +1}\sum _{j=0}^{\xi }[{\mathit{N}}_f{\left(\sum _{i=0}^{\mathsf{\Psi }}\overline{f}\left(\varrho \right),\sum _{i=0}^{\mathsf{\Psi }}\overline{\theta }\left(\varrho \right)\right]}^{2}d\varrho ,$$

(60)

$${\epsilon }_{\mathsf{\Psi }}^{\theta }=\frac{1}{\xi +1}\sum _{j=0}^{\xi }{\left[{\mathit{N}}_{\theta }{\left(\sum _{i=0}^{\mathsf{\Psi }}\overline{f}\left(\varrho \right),\sum _{i=0}^{\mathsf{\Psi }}\overline{\theta }\left(\varrho \right)\right)}_{n=j{D}_{u}n}\right]}^{2}d\varrho ,$$

(61)

$${\epsilon }_{\mathsf{\Psi }}^H=\frac{1}{\xi +1}\sum _{j=0}^{\xi }{\left[{\mathit{N}}_H{\left(\sum _{i=0}^{\mathsf{\Psi }}\overline{f}\left(\varrho \right),\sum _{i=0}^{\mathsf{\Psi }}\overline{H}\left(\varrho \right),\sum _{i=0}^{\mathsf{\Psi }}\overline{G}\left(\varrho \right)\right)}_{n=j{D}_{u}n}\right]}^{2}d\varrho ,$$

(62)

$${\epsilon }_{\mathsf{\Psi }}^G=\frac{1}{\xi +1}\sum _{j=0}^{\xi }{\left[{\mathit{N}}_G{\left(\sum _{i=0}^{\mathsf{\Psi }}\overline{f}\left(\varrho \right),\varpi \sum _{i=0}^{\mathsf{\Psi }}\overline{H}\left(\varrho \right),\sum _{i=0}^{\mathsf{\Psi }}\overline{G}\left(\varrho \right)\right)}_{n=j{D}_{u}n}\right]}^{2}d\varrho ,$$

(63)

Due to “Liao (2010)”

$${\epsilon }_{\mathsf{\Psi }}^t={\epsilon }_{\mathsf{\Psi }}^f+{\epsilon }_{\mathsf{\Psi }}^{\theta }+{\epsilon }_{\mathsf{\Psi }}^H+{\epsilon }_{\mathsf{\Psi }}^G,$$

(64)

where ${\epsilon }_{\mathsf{\Psi }}^t$ is the total squared residual error. Total average squared residual error is minimized by employing Mathematica package (BVPh 2.0, 2014.)

5. Error Analysis

An error analysis is carried out to ensure that the analysis is credible for minimum residual errors. The problem is solved both analytically and numerically by HAM and BVP4c, respectively. Analyses are carried out using 40th order of approximation. This analysis is also performed to ahead the reliability of HAM techniques using mathematica package BVPh 2.0 for minimum residual error ${10}^{-40}$. For authenticity and accuracy of the HAM solution, the results are compared with numerical solution of BVP4c using Matlab. HAM is an analytical method that offers solutions in series form. Although the BVP4c package is a numerical solver that uses non-uniform meshes for establishing error in each step size. Therefore, the BVP4c package is called collocation codes; they solve BVPs by computing a cube-spline on each subinterval of an interval mesh, so the method may be seen as a collocation or finite difference method with continuous extension. Error analyses shown in Figure 2 and

Table 1. Total residual error of $f\left(\eta \right)$, $\theta \left(\eta \right)$, $G\left(\eta \right)$ and $H\left(\eta \right)$ with fixed values of $Pr=0.5$, $S=-0.5$, $L=2$, $Q=4$, $R=0.5$, $\delta =2$, $Sc=K1=K2=Ha=0.5$.

m	${\mathit{ϵ}}^{\mathit{f}}\mathit{m}$	${\mathit{ϵ}}^{\mathit{\theta }}\mathit{m}$	${\mathit{ϵ}}^{\mathit{H}}\mathit{m}$	${\mathit{ϵ}}^{\mathit{G}}\mathit{m}$
1	$0.0112969$	$2.02804\times {10}^{-7}$	$3.20315\times {10}^{-7}$	$9.51647\times {10}^{-9}$
5	$5.95295\times {10}^{-15}$	$3.41791\times {10}^{-21}$	$1.21043\times {10}^{-21}$	$7.69318\times {10}^{-22}$
10	$1.98406\times {10}^{-29}$	$4.95817\times {10}^{-34}$	$1.03519\times {10}^{-34}$	$4.0926\times {10}^{-36}$
15	$1.88587\times {10}^{-32}$	$5.99659\times {10}^{-34}$	$1.05926\times {10}^{-34}$	$2.37732\times {10}^{-36}$
20	$1.78726\times {10}^{-32}$	$5.99659\times {10}^{-34}$	$1.05926\times {10}^{-34}$	$2.40741\times {10}^{-36}$
25	$1.73796\times {10}^{-32}$	$5.99659\times {10}^{-34}$	$1.05926\times {10}^{-34}$	$2.40741\times {10}^{-36}$
30	$1.73796\times {10}^{-32}$	$5.99659\times {10}^{-}34$	$1.05926\times {10}^{-34}$	$1.92593\times {10}^{-36}$
33	$1.73796\times {10}^{-32}$	$5.99659\times {10}^{-34}$	$1.05926\times {10}^{-34}$	$2.40741\times {10}^{-36}$
37	$1.73796\times {10}^{-32}$	$5.99659\times {10}^{-34}$	$1.05926\times {10}^{-34}$	$2.40741\times {10}^{-36}$
40	$1.73796\times {10}^{-32}$	$5.99659\times {10}^{-34}$	$1.05926\times {10}^{-34}$	$2.40741\times {10}^{-36}$

,

Table 2. Convergence of the homotopy solution for differing orders of calculation for $f^{″}\left(\eta \right)$, $-{\theta }^{\prime }\left(\eta \right)$, $-H^{\prime }\left(\eta \right)$ and $-G^{\prime }\left(\eta \right)$ when $S=-0.5$, $Pr=R=M=Sc=K1=K2=Ha=0.5$, $L=\delta =2$, $Q=4$.

$\mathit{\eta }$	${\mathit{f}}^{″}\left(\mathit{\eta }\right)$	$-{\mathit{\theta }}^{\prime }\left(\mathit{\eta }\right)$	$-{\mathit{G}}^{\prime }\left(\mathit{\eta }\right)$	$-{\mathit{H}}^{\prime }\left(\mathit{\eta }\right)$
0	$3.0185$	$-1.1595$	$0.7556$	$-0.3778$
$0.1001$	$2.4245$	$-1.1043$	$0.7556$	$-0.3778$
$0.2002$	$1.8116$	$-1.0582$	$0.7548$	$-0.3775$
$0.3003$	$1.2021$	$-1.0216$	$0.7533$	$-0.3769$
$0.4004$	$0.5795$	$-0.9927$	$0.7513$	$-0.3761$
$0.5005$	$-0.0320$	$-0.9717$	$0.7490$	$-0.3752$
$0.6006$	$-0.6448$	$-0.9572$	$0.7465$	$-0.3741$
$0.7007$	$-1.2437$	$-0.9488$	$0.7440$	$-0.3728$
$0.8008$	$-1.8148$	$-0.9460$	$0.7414$	$-0.3713$
$0.9009$	$-2.3714$	$-0.9478$	$0.7388$	$-0.3696$
1	$-2.8938$	$-0.9533$	$0.7360$	$-0.3675$

,

Table 3. Computations for $f\left(\eta \right)$, $\theta \left(\eta \right)$, $H\left(\eta \right)$ and $G\left(\eta \right)$ when $S=-2$, $Pr=R=Sc=K1=K2=Ha=0.5$, $L=\delta =2$, $Q=4$ and various values of $\eta$. Furthermore, $A=6.57027\times {10}^{-17}$ and $B=-4.0827\times {10}^{-19}.$

$\mathit{\eta }$	HAM Result				Numerical Result
	$\mathit{f}\left(\mathit{\eta }\right)$	$\mathit{\theta }\left(\mathit{\eta }\right)$	$\mathit{G}\left(\mathit{\eta }\right)$	$\mathit{H}\left(\mathit{\eta }\right)$	$\mathit{f}\left(\mathit{\eta }\right)$	$\mathit{\theta }\left(\mathit{\eta }\right)$	$\mathit{G}\left(\mathit{\eta }\right)$	$\mathit{H}\left(\mathit{\eta }\right)$
0	0	1	$0.6668$	$0.1650$	0	1	$0.6656$	$0.1648$
$0.1001$	$0.0140$	$0.8882$	$0.7001$	$0.1483$	$0.0141$	$0.8841$	$0.6989$	$0.1481$
$0.2002$	$0.0523$	$0.7808$	$0.7335$	$0.1317$	$0.0526$	$0.7743$	$0.7322$	$0.1315$
$0.3003$	$0.1086$	$0.6772$	$0.7670$	$0.1151$	$0.1092$	$0.6695$	$0.7655$	$0.1149$
$0.4004$	$0.1770$	$0.5764$	$0.8004$	$0.0985$	$0.1779$	$0.5686$	$0.7988$	$0.0983$
$0.5005$	$0.2513$	$0.4780$	$0.8339$	$0.0819$	$0.2524$	$0.4707$	$0.8321$	$0.0817$
$0.6006$	$0.3254$	$0.3811$	$0.8673$	$0.0654$	$0.3265$	$0.3750$	$0.8654$	$0.0652$
$0.7007$	$0.3933$	$0.2853$	$0.9007$	$0.0489$	$0.3942$	$0.2806$	$0.8989$	$0.0488$
$0.8008$	$0.4489$	$0.1899$	$0.9340$	$0.0325$	$0.4495$	$0.1868$	$0.9325$	$0.0324$
$0.9009$	$0.4864$	$0.0946$	$0.9672$	$0.0161$	$0.4866$	$0.0932$	$0.9663$	$0.0161$
1	$0.5$	A	1	B	$0.5000$	$-0.0000$	$1.0000$	0

,

Table 4. Computational for $f^{″}\left(0\right)$, $-{\theta }^{\prime }\left(0\right)$, $-H^{\prime }\left(0\right)$ and $-G^{\prime }\left(0\right)$ with $P_r=R=Sc=K1=K2=Ha=0.5$, $L=\delta =2$, $Q=4$ and differing values of S.

S	HAM Result				BVP4c Result
	${\mathit{f}}^{″}\left(\mathbf{0}\right)$	$-{\mathit{\theta }}^{\prime }\left(\mathbf{0}\right)$	$-{\mathit{H}}^{\prime }\left(\mathbf{0}\right)$	$-{\mathit{G}}^{\prime }\left(\mathbf{0}\right)$	${\mathit{f}}^{″}\left(\mathbf{0}\right)$	$-{\mathit{\theta }}^{\prime }\left(\mathbf{0}\right)$	$-{\mathit{H}}^{\prime }\left(\mathbf{0}\right)$	$-{\mathit{G}}^{\prime }\left(\mathbf{0}\right)$
$-0.5$	$3.0186$	$-1.1595$	$0.3328$	$-0.1664$	$3.0186$	$-1.1595$	$0.3328$	$-0.1664$
$-1$	$3.0143$	$-1.3111$	$0.3324$	$-0.1662$	$3.0143$	$-1.3111$	$0.3324$	$-0.1662$
$-1.5$	$3.0122$	$-1.4555$	$0.3320$	$-0.1660$	$3.0122$	$-1.4555$	$0.3320$	$-0.1660$
$-2$	$3.0157$	$-1.7257$	$0.3311$	$-0.1656$	$3.0157$	$-1.7257$	$0.3311$	$-0.1656$

,

Table 5. Computational for $f^{″}\left(0\right)$, $-{\theta }^{\prime }\left(0\right)$, $-H^{\prime }\left(0\right)$ and $-G^{\prime }\left(0\right)$ with $S=-0.5$, $P_r=R=Sc=K1=K2=0.5$, $L=\delta =2$, $Q=4$ and differing values of $Ha$.

$\mathit{Ha}$	HAM Result				Numerical Result
	${\mathit{f}}^{″}\left(\mathbf{0}\right)$	$-{\mathit{\theta }}^{\prime }\left(\mathbf{0}\right)$	$-{\mathit{H}}^{\prime }\left(\mathbf{0}\right)$	$-{\mathit{G}}^{\prime }\left(\mathbf{0}\right)$	${\mathit{f}}^{″}\left(\mathbf{0}\right)$	$-{\mathit{\theta }}^{\prime }\left(\mathbf{0}\right)$	$-{\mathit{H}}^{\prime }\left(\mathbf{0}\right)$	$-{\mathit{G}}^{\prime }\left(\mathbf{0}\right)$
$0.5$	$3.0186$	$-1.1595$	$0.3328$	$-0.1664$	$3.0186$	$-1.1595$	$0.3328$	$-0.1664$
1	$3.0916$	$-1.3111$	$0.3324$	$-0.1662$	$3.0916$	$-1.3111$	$0.3324$	$-0.1662$
$1.5$	$3.2102$	$-1.4555$	$0.3320$	$-0.1660$	$3.2102$	$-1.4555$	$0.3320$	$-0.1660$
2	$3.3702$	$-1.7257$	$0.3311$	$-0.1656$	$3.3702$	$-1.7257$	$0.3311$	$-0.1656$

,

Table 6. Computational for $f^{″}\left(0\right)$, $-{\theta }^{\prime }\left(0\right)$, $-H^{\prime }\left(0\right)$ and $-G^{\prime }\left(0\right)$ with $S=-0.5$, $P_r=Ha=R=Sc=K1=K2=0.5$, $L=2$, $Q=4$ and differing values of $\delta$.

D	HAM Result				Numerical Result
	${\mathit{f}}^{″}\left(\mathbf{0}\right)$	$-{\mathit{\theta }}^{\prime }\left(\mathbf{0}\right)$	$-{\mathit{H}}^{\prime }\left(\mathbf{0}\right)$	$-{\mathit{G}}^{\prime }\left(\mathbf{0}\right)$	${\mathit{f}}^{″}\left(\mathbf{0}\right)$	$-{\mathit{\theta }}^{\prime }\left(\mathbf{0}\right)$	$-{\mathit{H}}^{\prime }\left(\mathbf{0}\right)$	$-{\mathit{G}}^{\prime }\left(\mathbf{0}\right)$
1	$3.0186$	$-1.1595$	$0.3322$	$-0.3322$	$3.0186$	$-1.1595$	$0.3322$	$-0.3322$
2	$3.0186$	$-1.1595$	$0.3328$	$-0.1664$	$3.0186$	$-1.1595$	$0.3328$	$-0.1664$
3	$3.0186$	$-1.1595$	$0.3329$	$-0.1110$	$3.0186$	$-1.1595$	$0.3329$	$-0.1110$
4	$3.0186$	$-1.1595$	$0.3329$	$-0.0832$	$3.0186$	$-1.1595$	$0.3329$	$-0.0832$

,

Table 7. Computational for $f^{″}\left(0\right)$, $-{\theta }^{\prime }\left(0\right)$, $-H^{\prime }\left(0\right)$ and $-G^{\prime }\left(0\right)$ with $S=-0.5$, $Ha=R=Sc=K1=K2=0.5$, $\delta =L=2$, $Q=4$ and differing values of $P_r$.

$\mathit{Pr}$	HAM Result				Numerical Result
	${\mathit{f}}^{″}\left(\mathbf{0}\right)$	$-{\mathit{\theta }}^{\prime }\left(\mathbf{0}\right)$	$-{\mathit{H}}^{\prime }\left(\mathbf{0}\right)$	$-{\mathit{G}}^{\prime }\left(\mathbf{0}\right)$	${\mathit{f}}^{″}\left(\mathbf{0}\right)$	$-{\mathit{\theta }}^{\prime }\left(\mathbf{0}\right)$	$-{\mathit{H}}^{\prime }\left(\mathbf{0}\right)$	$-{\mathit{G}}^{\prime }\left(\mathbf{0}\right)$
1	$3.0186$	$-1.1595$	$0.3322$	$-0.3322$	$3.0186$	$-1.1595$	$0.3322$	$-0.3322$
2	$3.0186$	$-1.1595$	$0.3328$	$-0.1664$	$3.0186$	$-1.1595$	$0.3328$	$-0.1664$
3	$3.0186$	$-1.1595$	$0.3329$	$-0.1110$	$3.0186$	$-1.1595$	$0.3329$	$-0.1110$
4	$3.0186$	$-1.1595$	$0.3329$	$-0.0832$	$3.0186$	$-1.1595$	$0.3329$	$-0.0832$

,

Table 8. Computational for $f^{″}\left(0\right)$, $-{\theta }^{\prime }\left(0\right)$, $-H^{\prime }\left(0\right)$ and $-G^{\prime }\left(0\right)$ with $S=-0.5$, $P_r=Ha=R=Sc=K1=K2=0.5$, $L=\delta =2$ and differing values of Q.

Q	HAM Result				Numerical Result
	${\mathit{f}}^{″}\left(\mathbf{0}\right)$	$-{\mathit{\theta }}^{\prime }\left(\mathbf{0}\right)$	$-{\mathit{H}}^{\prime }\left(\mathbf{0}\right)$	$-{\mathit{G}}^{\prime }\left(\mathbf{0}\right)$	${\mathit{f}}^{″}\left(\mathbf{0}\right)$	$-{\mathit{\theta }}^{\prime }\left(\mathbf{0}\right)$	$-{\mathit{H}}^{\prime }\left(\mathbf{0}\right)$	$-{\mathit{G}}^{\prime }\left(\mathbf{0}\right)$
1	$3.0186$	$-1.0137$	$0.3328$	$-0.1664$	$3.0186$	$-1.0137$	$0.3328$	$-0.1664$
2	$3.0186$	$-1.0632$	$0.3328$	$-0.1664$	$3.0186$	$-1.0632$	$0.3328$	$-0.1664$
3	$3.0186$	$-1.1118$	$0.3328$	$-0.1664$	$3.0186$	$-1.1118$	$0.3328$	$-0.1664$
4	$3.0186$	$-1.1595$	$0.3328$	$-0.1664$	$3.0186$	$-1.1595$	$0.3328$	$-0.1664$

,

Table 9. Computational for $f^{″}\left(0\right)$, $-{\theta }^{\prime }\left(0\right)$, $-H^{\prime }\left(0\right)$ and $-G^{\prime }\left(0\right)$ with $S=-0.5$, $P_r=Ha=Sc=K1=K2=0.5$, $L=\delta =2$, $Q=4$ and differing values of R.

R	HAM Result				Numerical Result
	${\mathit{f}}^{″}\left(\mathbf{0}\right)$	$-{\mathit{\theta }}^{\prime }\left(\mathbf{0}\right)$	$-{\mathit{H}}^{\prime }\left(\mathbf{0}\right)$	$-{\mathit{G}}^{\prime }\left(\mathbf{0}\right)$	${\mathit{f}}^{″}\left(\mathbf{0}\right)$	$-{\mathit{\theta }}^{\prime }\left(\mathbf{0}\right)$	$-{\mathit{H}}^{\prime }\left(\mathbf{0}\right)$	$-{\mathit{G}}^{\prime }\left(\mathbf{0}\right)$
$0.1$	$3.0186$	$-1.2318$	$0.3328$	$-0.1664$	$3.0186$	$-1.2318$	$0.3328$	$-0.1664$
$0.2$	$3.0186$	$-1.2082$	$0.3328$	$-0.1664$	$3.0186$	$-1.2082$	$0.3328$	$-0.1664$
$0.3$	$3.0186$	$-1.1890$	$0.3328$	$-0.1664$	$3.0186$	$-1.1890$	$0.3328$	$-0.1664$
$0.4$	$3.0186$	$-1.1730$	$0.3328$	$-0.1664$	$3.0186$	$-1.1730$	$0.3328$	$-0.1664$

,

Table 10. Computational for $f^{″}\left(0\right)$, $-{\theta }^{\prime }\left(0\right)$, $-H^{\prime }\left(0\right)$ and $-G^{\prime }\left(0\right)$ with $S=-0.5$, $P_r=Ha=R=K1=K2=0.5$, $L=\delta =2$, $Q=4$ and differing values of $Sc$.

$\mathit{Sc}$	HAM Result				Numerical Result
	${\mathit{f}}^{″}\left(\mathbf{0}\right)$	$-{\mathit{\theta }}^{\prime }\left(\mathbf{0}\right)$	$-{\mathit{H}}^{\prime }\left(\mathbf{0}\right)$	$-{\mathit{G}}^{\prime }\left(\mathbf{0}\right)$	${\mathit{f}}^{″}\left(\mathbf{0}\right)$	$-{\mathit{\theta }}^{\prime }\left(\mathbf{0}\right)$	$-{\mathit{H}}^{\prime }\left(\mathbf{0}\right)$	$-{\mathit{G}}^{\prime }\left(\mathbf{0}\right)$
$0.1$	$3.0185$	$-1.1595$	$0.3332$	$-0.1666$	$3.0185$	$-1.1595$	$0.3332$	$-0.1666$
$0.2$	$3.0185$	$-1.1595$	$0.3331$	$-0.1666$	$3.0185$	$-1.1595$	$0.3331$	$-0.1666$
$0.3$	$3.0185$	$-1.1595$	$0.3330$	$-0.1665$	$3.0185$	$-1.1595$	$0.3330$	$-0.1665$
$0.4$	$3.0185$	$-1.1595$	$0.3329$	$-0.1664$	$3.0185$	$-1.1595$	$0.3329$	$-0.1664$

,

Table 11. Computational for $f^{″}\left(0\right)$, $-{\theta }^{\prime }\left(0\right)$, $-H^{\prime }\left(0\right)$ and $-G^{\prime }\left(0\right)$ with $S=-0.5$, $P_r=Ha=Sc=R=K2=0.5$, $L=\delta =2$, $Q=4$ and differing values of $K1$.

$\mathit{K}1$	HAM Result				Numerical Result
	${\mathit{f}}^{″}\left(\mathbf{0}\right)$	$-{\mathit{\theta }}^{\prime }\left(\mathbf{0}\right)$	$-{\mathit{H}}^{\prime }\left(\mathbf{0}\right)$	$-{\mathit{G}}^{\prime }\left(\mathbf{0}\right)$	${\mathit{f}}^{″}\left(\mathbf{0}\right)$	$-{\mathit{\theta }}^{\prime }\left(\mathbf{0}\right)$	$-{\mathit{H}}^{\prime }\left(\mathbf{0}\right)$	$-{\mathit{G}}^{\prime }\left(\mathbf{0}\right)$
$0.5$	$3.0186$	$-1.1595$	$0.3328$	$-0.1664$	$3.0186$	$-1.1595$	$0.3328$	$-0.1664$
1	$3.0186$	$-1.1595$	$0.3326$	$-0.1663$	$3.0186$	$-1.1595$	$0.3326$	$-0.1663$
2	$3.0186$	$-1.1595$	$0.3323$	$-0.1661$	$3.0186$	$-1.1595$	$0.3323$	$-0.1661$
3	$3.0185$	$-1.1595$	$0.3317$	$-0.1659$	$3.0185$	$-1.1595$	$0.3317$	$-0.1659$

and

Table 12. Computational for $f^{″}\left(0\right)$, $-{\theta }^{\prime }\left(0\right)$, $-H^{\prime }\left(0\right)$ and $-G^{\prime }\left(0\right)$ with $S=-0.5$, $P_r=Ha=Sc=R=K1=0.5$, $L=\delta =2$, $Q=4$ and differing values of $K2$.

$\mathit{K}2$	HAM Result				Numerical Result
	${\mathit{f}}^{″}\left(\mathbf{0}\right)$	$-{\mathit{\theta }}^{\prime }\left(\mathbf{0}\right)$	$-{\mathit{H}}^{\prime }\left(\mathbf{0}\right)$	$-{\mathit{G}}^{\prime }\left(\mathbf{0}\right)$	${\mathit{f}}^{″}\left(\mathbf{0}\right)$	$-{\mathit{\theta }}^{\prime }\left(\mathbf{0}\right)$	$-{\mathit{H}}^{\prime }\left(\mathbf{0}\right)$	$-{\mathit{G}}^{\prime }\left(\mathbf{0}\right)$
$0.5$	$3.0186$	$-1.1595$	$0.3328$	$-0.1664$	$3.0186$	$-1.1595$	$0.3328$	$-0.1664$
1	$3.0186$	$0.5010$	$0.5010$	$-0.2505$	$3.0186$	$0.5010$	$0.5010$	$-0.2505$
2	$3.0185$	$-1.1595$	$0.6704$	$-0.3352$	$3.0185$	$-1.1595$	$0.6704$	$-0.3352$
3	$3.0185$	$-1.1595$	$0.7556$	$-0.3778$	$3.0185$	$-1.1595$	$0.7556$	$-0.3778$

are provided to investigated the accuracy of both the methods for different involved physical parameters.

Figure 2 depicts the residual error at various orders of approximation for $f\left(\eta \right)$, $\theta \left(\eta \right)$, $H\left(\eta \right)$ and $G\left(\eta \right)$. It is obvious from the sub-figures that the error is persistently decreased up to 15th order of approximation. The cumulative residual error for the various approximation orders of fixed $Pr=0.5$, $S=-0.5$, $L=2$, $Q=4$, $R=0.5$, $\delta =2$, $Sc=K1=K2=Ha=0.5$ values is shown in Table 1. Table 2 has been shown at various orders of approximation with distinct average squared residual error. The comparison of the analytical and numerical values by HAM and BVP4c shown in Table 3 for various values of $\eta$ and fixed values of other parameters $Pr=-0.5$, $S=-0.5$, $L=0.01$, $Q=5$, $R=0.7$, $\delta =1$, $Sc=2.5$, $K1=0.1$, $K2=1$ and $Ha=0.1$.

6. Results and Discussions

Here, Equations (10)–(13) subject to the boundary conditions given in Equation (14) are used to model the flow between squeezing plates using the governing model for unsteady Newtonian fluid. These equations are solved and compared for numerical investigations through HAM and BVP$4c$. To investigate and analyze the problem graphically, Figure 3, Figure 4, Figure 5, Figure 6, Figure 7, Figure 8, Figure 9 and Figure 10 are drawn. The squeeze Reynolds number, $S=\frac{\alpha l^2}{2v_f}$, is the ratio between the normal velocity of the upper plate and the kinematic viscosity of the fluid. It is important to remember that the slow or fast vertical velocity of an upper plate to the lower plate is mean the small and large values of S. Furthermore, positive values of S implies that the upper plate moves out from the lower plate or increases the gap between the parallel plates, while negative values of S means that the top plate moves in the direction of the lower plate or decreases the distance between the plates. Figure 3 shows the influence of S on $f^{\prime }\left(\eta \right)$, $f\left(\eta \right)$, $\theta \left(\eta \right)$ and $H\left(\eta \right)$. The increase in the distance between the plates, i.e., $S=-0.5$, $-2.5$, $-4.5$, $-6.5$ with fixed values of $P_r=0.1$, $Ha=1$, $R=0.05$, $L=0.01$, $K1=K2=Sc=0.5$, $Q=10$, $\delta =5$ will allow the fluid to move in the x-axis direction, which reduces the radial velocity $f^{\prime }\left(\eta \right)$. As the fluid moves through the central region, it begins to increase as seen in Figure 3a. This occurrence also reduces the velocity when the velocity leads in the squeezing effect. It is also obvious from Figure 3b, that the fluid friction makes squeezing of plates, which generate heat and this heat raises the fluid temperature.

Impacts of relevant parameters on velocity, temperature distributions are noticed graphically in this part. The squeezing velocity distribution $f\left(\eta \right)$, is noticed for the parameter $H_a$ in Figure 4a. The magnetic field makes the resistive force act between the particles causing about the smaller velocity. Thus, diminishing Ha expands the velocity. Ha oppositely affects the temperature $\theta \left(\eta \right)$ (see Figure 4b). This resistive force produces heat that outcomes in the reducing in temperature. Figure 6 shows the impact of Prandtl number on $f\left(\eta \right)$ and $\theta \left(\eta \right)$. It is noticed that an expansion in $P_r$ causes a decline in temperature. Most extreme reduction is found in the fluid field. Furthermore, a opposite conduct is noted for $\theta \left(\eta \right)$. Figure 7 portrays the impact of heat generation Q on non-dimensional temperature and velocity models on $f\left(\eta \right)$ and $\theta \left(\eta \right)$, individually. We saw that the temperature increment with an increase in Q. The result is based on the choices that any heat generated would increase in the liquid’s nuclear power. We also saw that the fluid temperature is higher to a smaller than the fluid temperature. Figure 8 depicts the impact of Radiation parameter R on the velocity and temperature profiles of the flow. It is obvious to say that the radiation parameter have inclination to diminish the velocity and temperature profiles of the flow. This may occur because of the way that lessens in the radiation produces heat energy to the flow, which leads to decline in the velocity and temperature of the liquid. The influence on concentration profiles $f\left(\eta \right)$, $H\left(\eta \right)$ and $G\left(\eta \right)$ of the equations to describe of a homogenous parameter $k1$ and the equations to describe of a heterogeneous parameter $k2$ can be observed in Figure 9, respectively. The increase in the concentration profiles $f\left(\eta \right)$, $H\left(\eta \right)$ and $G\left(\eta \right)$, is observed due to increase in $k1$. According to this, an increase in the homogeneous chemical reaction parameter which decreases viscosity. However, the heterogeneous $k2$ parameter indicates the reverse of the $k1$ above can seen in Figure 9. Because of the diffusion decreases with an increase in $k2$ and a smaller amount of diffused particles will reduce in the concentration of $f\left(\eta \right)$ and $H\left(\eta \right)$ but it will increase in the concentration of $G\left(\eta \right)$. The Schmidt number $Sc$’s effects can be seen in Figure 10, in which the value of $Sc$ increases. As the fluid concentration profiles $f\left(\eta \right)$ and $\theta \left(\eta \right)$ increases because of the ratio viscus diffusion Sc to the molecular diffusion rate. Table 4, Table 5, Table 6, Table 7, Table 8, Table 9, Table 10, Table 11 and Table 12 are made to numerically examine the impact of various definite parameters. It can be seen from given tables that all consequences have excellent agreement with the results found through BVP4c and HAM. It is seen that the impacts of skin friction coefficient, velocity, temperature, both homogeneous and heterogeneous parameters have the tendency to increase the mass transfer rate. However, it does not depict a huge change in friction factor and heat transfer rate. An enhancement of the internal heat generation depreciates the coefficient of skin friction and the rate of heat. We have found similar results with an increase in the Squeeze parameter S. The friction factor decreases because of an increase in the squeeze parameter S and increasing in the local Nusselt and Sherwood numbers. A decrease in the skin friction coefficient was found, indicating that the surface drags on the fluid. Table 4, Table 5, Table 6, Table 7, Table 8, Table 9, Table 10, Table 11 and Table 12 are made to show the effects of $f^{″}\left(0\right)$, $f\left(0\right)$, $\theta \left(0\right)$, $G\left(0\right)$, $H^{\prime }\left(0\right)$, $-{\theta }^{\prime }\left(0\right)$, $-G^{\prime }\left(0\right)$ and $-H^{\prime }\left(0\right)$, respectively.

7. Conclusions

This paper offers the numerical solution for analyzing the consequence of the flow between two squeezing plates with the homogeneous and heterogeneous reaction in the presence of internal heat generation/absorption and thermal radiation. The effect of non-dimensional leading parameters on velocity, temperature and concentration profiles with the friction factor, local Nusselt, and Sherwood numbers are explained and displayed using graphs and tables. The conclusions of the present analysis are as follows:

• It is found that increasing the distance between discs, decrease the radial velocity and increase tangential velocity, while decreasing the distance between discs have shown opposite behavior.
• Homogeneous and heterogeneous parameters assist in the monitoring of the stream’s utility profiles.
• It is also found that increase in prandtl number increasing skin friction, mass flux while it decrease heat flux.
• Table 2 depicts a drop in the skin friction and a decreases in the local Nusselt and Sherwood numbers.
• Table 4 depict that as the friction factor, local Nusselt and Sherwood numbers all are reduces as the squeezing parameter.
• As heat generation $Ha$ increases the friction factor decreases, it can be seen from Table 5 that the local Nusselt and Sherwood numbers decrease in the fluid viscosity.
• It is also observed from Table 8 and Table 9 that Heat generation Q and radiation R parameter tend to reduce the friction factor.