**2. Methods**

A rotating flow of hydromagnetic, time independent and incompressible hybrid nanofluid between two parallel disks in three dimensions is analyzed. Homogeneous-heterogeneous chemical reactions are also considered. The lower disk is supposed to locate at *z* = 0 while the upper disk is at a constant distance *H* apart. The velocities and stretching on these disks are (Ω1, Ω2) and (*a*1, *a*2), respectively while the temperatures on these disks are *T*<sup>1</sup> and *T*2, respectively. A magnetic field of strength *B*<sup>0</sup> is applied in the direction of *z*-axis (please see Figure 1). Ethylene glycol is chosen for the base fluid in which zinc oxide and gold nanoparticles are added.

For cubic auto-catalysis, the homogeneous reaction is

$$\mathfrak{C}\mathfrak{C} + \mathcal{B} \to \mathfrak{BC}, \quad rate = c^2 b k\_c. \tag{1}$$

The first order isothermal reaction on the surface of catalyst is

$$B \to \mathbb{C}, \quad rate = b\mathbf{k}\_{\rm s} \tag{2}$$

where *B* and *C* denote the chemical species with concentrations *b* and *c*, respectively. *k<sup>c</sup>* and *k<sup>s</sup>* are the rate constants.

**Figure 1.** Geometry of the problem.

Cylindrical coordinates (*r*, *ϑ*, *z*), are applied to provide the thermodynamics of hybrid nanofluid as [57–59]

$$
\frac{\partial w}{\partial z} + \frac{\partial u}{\partial r} + \frac{u}{r} = 0,\tag{3}
$$

$$\rho\_{\rm Inf} \left( -\frac{v^2}{r} + \frac{\partial u}{\partial r} u + \frac{\partial u}{\partial z} w \right) = \mu\_{\rm Int} \left( \frac{\partial^2 u}{\partial z^2} + \frac{\partial^2 u}{\partial r^2} - \frac{u}{r^2} + \frac{\partial u}{\partial r} \frac{1}{r} \right) - \sigma\_{\rm Int} B\_0^2 u - \frac{\mu\_{\rm Int}}{S} u^2 - S\_1 u^2 - \frac{\partial P}{\partial r'} \tag{4}$$

$$\rho\_{\rm Inf} \left( \frac{w}{r} + w \frac{\partial v}{\partial z} + u \frac{\partial v}{\partial r} \right) = \mu\_{\rm Int} \left( \frac{\partial^2 v}{\partial z^2} + \frac{\partial^2 v}{\partial r^2} - \frac{v}{r^2} + \frac{1}{r} \frac{\partial v}{\partial r} \right) - \sigma\_{\rm Int} B\_0^2 v - \frac{\mu\_{\rm Inf}}{S} v^2 - S\_1 v^2 \tag{5}$$

$$
\rho\_{\ln f} \left( w \frac{\partial w}{\partial z} + u \frac{\partial w}{\partial r} \right) = -\frac{\partial P}{\partial z} + \mu\_{\ln f} \left( \frac{\partial^2 w}{\partial z^2} + \frac{\partial^2 w}{\partial r^2} + \frac{1}{r} \frac{\partial w}{\partial r} \right) - \frac{\mu\_{\ln f}}{S} w^2 - S\_1 w^2 \tag{6}
$$

$$(\rho c\_p)\_{\text{hnf}} \left( w \frac{\partial T}{\partial \overline{z}} + u \frac{\partial T}{\partial \overline{r}} \right) = \left( k\_{\text{hnf}} + \frac{16T\_1^3 \sigma\_1}{3k\_0} \right) \left( \frac{\partial^2 T}{\partial z^2} + \frac{\partial^2 T}{\partial r^2} + \frac{1}{r} \frac{\partial T}{\partial \overline{r}} \right) + \sigma\_{\text{hnf}} B\_0^2 (v^2 + u^2), \tag{7}$$

$$4w\frac{\partial b}{\partial z} + \mu \frac{\partial b}{\partial r} = -c^2 b k\_c + D\_B \left( \frac{\partial^2 b}{\partial z^2} + \frac{\partial^2 b}{\partial r^2} + \frac{1}{r} \frac{\partial b}{\partial r} \right) \, , \tag{8}$$

$$w\frac{\partial c}{\partial z} + u\frac{\partial c}{\partial r} = c^2 b k\_c + D\_\mathcal{C} \left( \frac{\partial^2 c}{\partial z^2} + \frac{\partial^2 c}{\partial r^2} + \frac{1}{r} \frac{\partial c}{\partial r} \right). \tag{9}$$

The boundary conditions are

$$\text{at } \ z = 0, \quad D\_{\mathbb{C}} \frac{\partial c}{\partial z} = k\_{\mathfrak{s}} c, \quad D\_{\mathbb{B}} \frac{\partial b}{\partial z} = k\_{\mathfrak{s}} b, \quad T = T\_1, \quad w = 0, \qquad v = r \Omega\_1, \quad u = n\_1. \tag{10}$$

$$\begin{array}{ccccccccc}\text{at} & z = H, & c \to 0, & b \to b\_{0\prime} & T = T\_{\Sigma}, & P \to \infty, & w = 0, & \upsilon = r\Omega\_{\Sigma}, & u = m\_{\Sigma} \end{array} \tag{11}$$

where *u*(*r, ϑ, z*), *v*(*r, ϑ, z*) and *w*(*r, ϑ, z*) are the components of velocity, *P* is the pressure. *S* is the permeability of porous medium, *S*<sup>1</sup> = *Cb rS* 1 2 is the non-uniform inertia coefficient of porous medium with *C<sup>b</sup>* as the drag coefficient. Temperature of hybrid nanofluid is *T* and *B* = (0, 0, *B*0) is the magnetic field. *σ*<sup>1</sup> is the Stefan Boltzmann constant and *k*<sup>0</sup> is the absorption coefficient. For the hybrid nanofluid, the important quantities are *ρhn f* (density), *µhn f* (dynamic viscosity), *σhn f* (electrical conductivity), (*cp*)*hn f* (heat capacity) and *khn f* (thermal conductivity). The subscript "hnf" shows the hybrid nanofluid. For the thermal conductivity, the mathematical formulation is obtained via Hamilton–Crosser model [9] as

$$\frac{k\_{nf}}{k\_f} = \frac{k\_1 + (n\_1 - 1)k\_f - (n\_1 - 1)(k\_f - k\_1)\phi\_1}{k\_1 + (n\_1 - 1)k\_f + (k\_f - k\_1)\phi\_1},\tag{12}$$

where *n* is the empirical shape factor for the nanoparticle whose value is given in Table 1.


**Table 1.** Values of shape factor of different shapes of nanoparticles.

The subscript "*f* " denotes the base fluid namely ethylene glycol and the subscript "*nf* " is used for nanofluid. *ρ<sup>s</sup>* and *(cP)<sup>s</sup>* are the density and heat capacity at specified pressure of nanoparticles, respectively. *φ*<sup>1</sup> is the first nanoparticle volume fraction while *φ*<sup>2</sup> is the second nanoparticle volume fraction which can be formulated as [57].

$$\rho\_s = \frac{(\rho\_1 \times m\_1) + (\rho\_2 \times m\_2)}{m\_1 + m\_2},\tag{13}$$

$$(c\_P)\_s = \frac{((c\_P)\_1 \times m\_1) + ((c\_P)\_2 \times m\_2)}{m\_1 + m\_2} \,\,\,\tag{14}$$

$$\phi\_1 = \frac{\frac{m\_1}{\rho\_1}}{\frac{m\_1}{\rho\_1} + \frac{m\_2}{\rho\_2} + \frac{m\_f}{\rho\_f}} \,\,\,\tag{15}$$

$$\phi\_2 = \frac{\frac{m\_2}{\rho\_2}}{\frac{m\_1}{\rho\_1} + \frac{m\_2}{\rho\_2} + \frac{m\_f}{\rho\_f}} \tag{16}$$

$$
\phi = \phi\_1 + \phi\_{2'} \tag{17}
$$

where *m*1, *m*<sup>2</sup> and *m<sup>f</sup>* are, respectively the mass of first nanoparticle, mass of the second nanoparticle and mass of the base fluid. *φ* is the total volume fraction of zinc oxide and gold nanoparticles.

The thermophysical properties of C2H6O<sup>2</sup> as well as nanoparticles are given in Table 2.


**Table 2.** Thermophysical properties of ethylene glycol and nanoparticles.

The mathematical formulations for *ρhn f* (density), *µhn f* (dynamic viscosity), *σhn f* (electrical conductivity), (*cp*)*hn f* (heat capacity) are given in Table 3 where *φ<sup>s</sup>* shows the particle concentration.


**Table 3.** Mathematical forms of thermophysical properties.

Following transformations are used

$$f'(\zeta)\Omega\_1 r = u, \quad v = r\Omega\_1 g(\zeta), \quad -2f(\zeta)H\Omega\_1 = w, \quad \frac{-T\_2 + T}{-T\_2 + T\_1} = \theta(\zeta), \quad \Omega\_1 \rho\_f v\_f \left(\frac{\varepsilon r^2}{2H^2} + P(\zeta)\right) = P,$$

$$\varrho b\_0 = b, \quad c = b\_0 \varrho\_1, \quad \frac{z}{H} = \zeta,\tag{18}$$

where *ν<sup>f</sup>* = *µf ρf* is the kinematic viscosity and *e* is the pressure parameter.

Using the values from Equation (18) in Equations (4)–(11), the following eight Equations (19)–(26) are obtained

$$\mathrm{Re}f'''' + \mathrm{Re}\left(2f'' - f^2 + g^2 - \mathrm{MB}\_2f\right) - \varepsilon - k\_2 \mathrm{Re}\mathrm{B}\_1 f - k\_3 \mathrm{Re}\frac{1}{\rho\_{\mathrm{Imf}}}(f')^2 = 0,\tag{19}$$

$$\left(B\_1 \mathbf{g}^{\prime\prime} + \text{Re}\left(2\mathfrak{f}\mathbf{g}^{\prime} - MB\_2 \mathbf{g}^{\prime}\right) - k\_2 \text{Re}B\_1 \mathbf{g} - k\_3 \text{Re}\frac{1}{\rho\_{\text{mf}f}}(\mathbf{g})^2 = 0,\tag{20}$$

$$P' = \frac{2}{k\_2} - 4\text{Ref}f' - f',\tag{21}$$

$$B\_3 \frac{k\_{\rm Inf}}{k\_f} \theta'' + \frac{1}{\rm Rd} Pr \text{Re} \left[ 2f \theta' + M \text{EcB}\_4 \left( g^2 + (f')^2 \right) \right] = 0,\tag{22}$$

$$\text{ScRe}\left(2\varphi'f - k\_4\varphi\varphi\_1^2\right) + \varphi'' = 0,\tag{23}$$

$$
\varphi\_1^{\prime\prime} + \mathrm{ScRe}\left(2\varphi\_1^{\prime}f + k\_4\varphi\varphi\_1^2\right)\frac{1}{k\_5} = 0,\tag{24}
$$

$$f = 0, \quad f' = k\_{\theta}, \quad g = 1, \quad \theta = 1, \quad \phi' = k\_7 \phi, \quad k\_4 \phi'\_1 = -k\_7 \phi, \quad P = 0 \quad \text{at} \quad \zeta = 0,\tag{25}$$

$$f = 0, \quad f' = k\_{\theta}, \quad g = \Omega, \quad \theta = 0, \quad \varphi = 1, \quad \varphi\_1 = 0 \quad \text{at} \quad \zeta = 1,\tag{26}$$

where (0 ) represents the derivative with respect to *ζ*. *B*<sup>1</sup> = " 1 − *m*1 *ρ*1 *m*1 *ρ*1 + *m*2 *ρ*2 + *mf ρ f* #−2.5 × " 1 − *m*1 *ρ*1 *m*1 *ρ*1 + *m*2 *ρ*2 + *mf ρ f* + *m*1 *ρ*1 *m*1 *ρ*1 + *m*2 *ρ*2 + *mf ρ f ρs ρ f* #−<sup>1</sup> , *B*<sup>2</sup> = 1 + 3 *σ*1*φ*1+*σ*2*φ*<sup>2</sup> *σf* − (*φ*<sup>1</sup> + *φ*2) 2+ *σ*1*φ*1+*σ*2*φ*<sup>2</sup> (*φ*1+*φ*2)*σ<sup>f</sup>* − *σ*1*φ*1+*σ*2*φ*<sup>2</sup> *σf* − (*φ*<sup>1</sup> + *φ*2) , *B*<sup>3</sup> = (*ρcP*)*<sup>f</sup>* " 1 − *m*1 *ρ*1 *m*1 *ρ*1 + *m*2 *ρ*2 + *mf ρ f* # *ρ<sup>f</sup>* + " 1 − *m*1 *ρ*1 *m*1 *ρ*1 + *m*2 *ρ*2 + *mf ρ f* # *ρs* × " 1 − *m*1 *ρ*1 *m*1 *ρ*1 + *m*2 *ρ*2 + *mf ρ f* # (*cP*)*<sup>f</sup>* + " 1 − *m*1 *ρ*1 *m*1 *ρ*1 + *m*2 *ρ*2 + *mf ρ f* # (*cP*)*<sup>s</sup>* , *B*<sup>4</sup> = *σhn f ρhn f* . *k*<sup>2</sup> = *νf S*Ω<sup>1</sup> is the porosity parameter, *k*<sup>3</sup> = *Cb S* 1 2 is the inertial parameter due to Darcy Forchheimer effect. The other non-dimensional parameters are Ω = Ω<sup>2</sup> Ω<sup>1</sup> , *Re* = Ω1*H*<sup>2</sup> *νf* , *M* = *σ<sup>f</sup> B* 2 0 *ρ <sup>f</sup>* Ω<sup>1</sup> , *Rd* = 16*σ*1*T* 3 1 3*k <sup>f</sup> k*0 , *Pr* = (*ρcP* )*hn f ν<sup>f</sup> k f* , *Ec* = *r* <sup>2</sup>Ω<sup>2</sup> 1 *cP*(*T*1−*T*2) , *Sc* = *νf D<sup>B</sup>* , *k*<sup>4</sup> = *kcb* 2 0 Ω<sup>1</sup> , *k*<sup>5</sup> = *D<sup>C</sup> D<sup>B</sup>* , *k*<sup>6</sup> = *a*1 Ω<sup>1</sup> , *k*<sup>7</sup> = *ksH D<sup>B</sup>* and *k*<sup>8</sup> = *a*2 Ω<sup>1</sup> which are known as rotation parameter, Reynolds number, magnetic field parameter, thermal radiation parameter, Prandtl number, Eckert number, Schmidt number, homogeneous chemical reaction parameter, diffusion

stretching parameter at upper disk, respectively. Regarding the homogeneous-heterogeneous chemical reaction, the quantities *B* and *C* may be considered in a special case, i.e., if *D<sup>B</sup>* is equal to *DC*, then in such a case *k*<sup>5</sup> equals unity, which leads to

coefficient ratio, stretching parameter for lower disks, heterogeneous chemical reaction parameter and

$$1 = \varrho\_1(\zeta) + \varrho(\zeta). \tag{27}$$

Using Equation (27), Equations (23) and (24) generate

$$0 = \text{ScRe}\left[\left(1 - \varphi\right)^2 k\_4 \varphi + 2\varphi' f\right] + \varphi'',\tag{28}$$

whose corresponding boundary conditions become

$$
\varphi' = k \varphi \quad \text{for} \quad \zeta = 0 \quad \text{while} \quad \varphi = 1 \quad \text{for} \quad \zeta = 1. \tag{29}
$$

By taking derivative of Equation (19) with respect to *ζ*, it becomes

$$\mathrm{Re}f^{\prime\prime\prime} + \mathrm{Re}\left(2f^{\prime\prime\prime} + 2\mathrm{gg}^{\prime} - \mathrm{MB}\_{2}f^{\prime}\right) - k\_{2}\mathrm{Re}\mathrm{B}\_{1}f^{\prime} - 2k\_{3}\mathrm{Re}\frac{1}{\rho\_{\mathrm{lnf}}}f^{\prime} = 0. \tag{30}$$

Considering Equation (21), Equations (25) and (26), the quantity *e* is computed as

$$\epsilon = f'''(0) - \text{Re}\left[ -\left(g(0)\right)^2 + \left(f'(0)\right)^2 + MB\_2f'(0) + \frac{1}{B\_1k\_2f'(0)} \right].\tag{31}$$

Integrating Equation (21) with respect to *ζ* by using the limit 0 to *ζ* for evaluating *P* as

$$P = -2\left[\text{Re}\left( (f)^2 + \frac{1}{k\_2} \int\_0^\zeta f\right) \left( f - f'(0) \right) \right]. \tag{32}$$

*Skin Frictions and Nusselt Numbers*

The important physical quantities are defined as

 $\mathbb{C}\_{f\_1}(\text{Local skin friction at lower disk}) = \frac{\tau|\_{z=0}}{\rho\_{\text{Inf}}(r\Omega\_1)^2},$ 
$$\mathbb{C}\_{f\_2}(\text{Local skin friction at upper disk}) = \frac{\tau|\_{z=H}}{\rho\_{\text{Inf}}(r\Omega\_1)^2}, \quad \text{(33)}$$

*Crystals* **2020**, *10*, 1086

where

$$
\pi = \sqrt{(\tau\_{zr})^2 + (\tau\_{z\theta})^2} \,\,\,\tag{34}
$$

denotes the sum of shear stress of tangential forces *τzr* and *τz<sup>θ</sup>* along radial and tangential directions which are defined as

$$\begin{aligned} \text{tr}\_{\text{2r}} \left( \text{Shear stress } friction \text{ at lower disk} \right) &= \mu\_{\text{lnf}} \frac{\partial u}{\partial z}|\_{z=0} = \frac{\mu\_{\text{lnf}} r \Omega\_1 f''(0)}{H} \text{ and} \\ \text{tr}\_{\text{2\theta}} &= \mu\_{\text{lnf}} \frac{\partial v}{\partial z}|\_{z=0} = \frac{\mu\_{\text{lnf}} r \Omega\_1 g'(0)}{H} . \end{aligned} \text{ and }$$

Using the information of Equations (34) and (35), Equation (33) proceeds to

$$\mathcal{C}\_{f\_1} = \frac{1}{\mathcal{Re}\_r} \left[ 1 - \frac{\frac{m\_1}{\rho\_1}}{\frac{m\_1}{\rho\_1} + \frac{m\_2}{\rho\_2} + \frac{m\_f}{\rho\_f}} \right]^{-2.5} \left[ \left( f''(0) \right)^2 + \left( g'(0) \right)^2 \right]^{\frac{1}{2}} \tag{36}$$

$$\mathcal{C}\_{f\_2} = \frac{1}{\mathcal{Re}\_r} \left[ 1 - \frac{\frac{m\_1}{\rho\_1}}{\frac{m\_1}{\rho\_1} + \frac{m\_2}{\rho\_2} + \frac{m\_f}{\rho\_f}} \right]^{-2.5} \left[ \left( f''(1) \right)^2 + \left( g'(1) \right)^2 \right]^{\frac{1}{2}} \tag{37}$$

where *Re<sup>r</sup>* = *r*Ω1*H νhn f* is the Reynolds number.

Another important physical quantity is

$$\text{Nu}\_{\text{T}\_1} \left( \text{Local Nusselt number at lower disk} \right) \\ = \frac{Hq\_w}{k\_f (T\_1 - T\_2)} |\_{z=0}$$

$$\text{Nu}\_{\text{T}\_2} \left( \text{Local Nusselt number at upper disk} \right) \\ = \frac{Hq\_w}{k\_f (T\_1 - T\_2)} |\_{z=H\_V} \quad \text{(38)}$$

where *q<sup>w</sup>* is the surface temperature defined as

$$q\_{\rm pw} \left( At \, lower \, disk \right) \, = -k\_{\rm hnf} \frac{\partial T}{\partial z}|\_{z=0} = -k\_{\rm hnf} \frac{T\_1 - T\_2}{H} \theta'(0). \tag{39}$$

Taking information from Equation (39), Equation (38) becomes

$$Nu\_{r\_1} = -\frac{k\_{\ln f}}{k\_f} \theta'(0), \quad Nu\_{r\_2} = -\frac{k\_{\ln f}}{k\_f} \theta'(1). \tag{40}$$
