Non-Linear Dynamic Movements of CNT/Graphene/Aluminum Oxide and Copper/Silver/Cobalt Ferrite Solid Particles in a Magnetized and Suction-Based Internally Heated Surface: Sensitivity and Response Surface Optimization

Abstract

Hybrid nanofluids combine two or more nano properties with a base fluid such as water ethylene. Usually, this helps enhance the heat transfer rate; in this article, using new similarity transformations created by Lie group analysis, the governing nonlinear partial differential equations are transformed into a system of connected nonlinear ordinary differential equations. The resulting design is numerically solved using a BVP4C solver with the shooting method (MATLAB). The magneto hydrodynamic flow of an incompressible fluid and the rate of heat and mass transfer were investigated for two cases, with various nanoparticle shapes including cylindrical, spherical, and platelet. Case 1 was CNT (1%), graphene (1%), and aluminum oxide (1%), and Case 2 was copper (1%), silver (1%), and cobalt ferrite (1%). When the Hartmann number rises, velocity and temperature exhibit inverse behavior: the velocity profile increases, and the temperature profile decreases. When the suction rises, the velocity and temperature profiles both increase. Optimization techniques were used from response surface methodology (RSM) to set factorial variables so that the response met the desired maximum or minimum value. Factorial methods like ANOVA were used to model the response, but they were expanded to simulate the effects in terms of extrapolation.

1. Introduction

A potent, comprehensive, and systematic approach to obtaining group-invariant solutions, also known as self-similarity transformations, is provided by Lie group methods and associated invariants. The independent variable numbers of a group of partial differential equations were reduced using self-similarity transformations, which resulted in the change of the non-linear controlling partial differential equations into ordinary differential equations. The fundamental instrument for this investigation was the application of the matching ‘infinitesimal’ Lie algebra representations. Continual series of transformations applied to a more extensive set of variables, including an equation’s parameters and independent and dependent variables, is what we mean when we refer to a group of developed Lie transformations of partial differential equations. Recent studies on the Lie group transformations theory and its uses in other disciplines are established in [1,2,3,4]. When the temperature and concentration gradients are significant, the impact of thermal-diffusion spread in chemical reactions becomes more effective. The convective transport in porous material examined by Li et al. [5] also significantly benefits from the dispersion, inertial, and injection/suction effects. For a hydromagnetic flow, i.e., an electrically conducting fluid over an expanding vertical sheet with a chemical reaction at varying stream conditions, Kandaswamy et al. [6] employed the scaling group of transformation. To evaluate the increase of the steady border layer flow, nanoparticle volume fraction, and heat transmitted over a porous expanding surface in a nanofluid for different parameters, Abdul et al. [7] used a scaling group of changes. Using the Lie symmetry group of transformations, Siva et al. [8] examined the continuous border layer flow and heat transmission evolution for different parameters using a porous wedge sheet in a nanofluid. Kandasamy et al. [9] explored the border layer fluid flow of seawater and water-based nano-sized particles through a porous wedge under the uniform transverse thermal energy radiation and magnetic field. Alessa et al. [10] explored the mixed effects of radiative heat transmission and a non-Newtonian fluid with viscous and elastic properties through a channel with a porous media for optically thin liquid.

Hybrid nanofluids that are useful in cooling, microelectronics, temperature enhancement and reduction, vehicle thermal management, and pharmaceutical processes have been identified. Through a rigorous investigation of the liquid above, a base fluid with two different nanosized particles has been found (i.e., the mixture of nanofluid). The dynamic viscosity of aluminium-based nanofluids varies according to the morphologies of nanoparticles at various temperatures, as shown in Timofeeva et al. [11]. Depending on the surface charge, different types of nanoparticles (cylinders, blades, bricks, and platelets) collect and interact with the base fluid differently. This is supported by Sahu and Sarkar’s [12] finding that “Nanoparticle shapes affect both the energy and energetic performance”. The dynamics of nanofluids produced by thermal-capillary consecutiveness caused by five different nanoparticle morphologies were characterized by Jiang et al. [13] (brick, platelet, sphere, cylinder, and blade). The proper operation of numerous industrial components depends on heat transmission. The natural convection of a nanofluid composed of magnetized carbon nanotubes in a curved cage has not been thoroughly studied, according to researchers. Research conducted by Izadi et al. [14] took a variety of flow and geometric factors into account when analyzing the natural convection of a nanofluid interior space. Nayak [15] concluded that thermal radiation and viscous dissipation are to blame for the decrease in heat transmission after considering the impact of thermal radiative on the depth of molecules. Rasool et al. [16] explored the numerical investigation of EMHD nanofluid flows over a convectively heated riga pattern positioned horizontally in a Darcy-Forchheimer porous medium: Application of passive control strategy and generalized transfer laws. Hassan et al. [17] tested the Buongiorno model by increasing the thermal conduction of the Falkner–Skan magnetic nanofluid while microorganisms were present. Numerical simulation of a thermally enhanced EMHD flow of a heterogeneous micropolar mixture comprising (60%)-ethylene glycol (EG), (40%)-water (W), and copper oxide nanomaterials (CuO) was studied by Shah et al. [18]. Using slip boundary conditions and convection, Khilap et al. [19] investigated the rate-of-momentum profile and the thermal heat energy transfer in the flow of a micropolar fluid on a diminishing porous surface. Sowmya et al. [20] studied the impact of perpendicular shapes on the flow of hybrid kerosene/Fe₃O₄–Ag and water/Fe₃O₄–Ag nanofluids under slip conditions and nonlinear thermal radiation circumstances. Lou et al. [21] studied micropolar dusty fluid: Coriolis force effects on dynamics of MHD rotating fluid when Lorentz force is significant.

In the present article, we employed the Lie group scaling transformations to convert from dimensional to non-dimensional equations and analyzed the heat and velocity transfer rate under the effect of magneto-hydrodynamics, suction, and nanoparticles with different shapes such as spherical, cylindrical, and platelet.

2. Mathematical Formulation

An incompressible tangent hyperbolic fluid flows in a continuous, two-dimensional flow. We suppose that the liquid is perpendicular to the plane at $y=0$, limiting the flow to the region at $y>0$. We also suppose that linear stretching is responsible for fluid production. The liquid’s governing equation is given by the following equation, as described in [22]:

$$\overline{\tau }=\left[{\mu }_{\infty }+\left({\mu }_0+{\mu }_{\infty }\right)\tanh {\left(\Gamma \overline{\dot{\gamma }}\right)}^n\right]\overline{\dot{\gamma }}.$$

(1)

The excess stress tensor, “infinite shear rate viscosity” and “zero shear rate viscosity”, time-dependent material fixed value, and flow behavior index are denoted by $\overline{\tau }$, ${\mu }_0$, ${\mu }_{\infty }$, $\Gamma$, and $n$. $\overline{\dot{\gamma }}$ is defined by

$$\overline{\dot{\gamma }}=\sqrt{\frac{1}{2}{\displaystyle \sum {}_i{\displaystyle \sum {}_j{\overline{\dot{\gamma }}}_{ij}{\overline{\dot{\gamma }}}_{ji}}}}=\sqrt{\frac{1}{2}\Pi }$$

(2)

in collaboration with $\Pi =\frac{1}{2}tr{\left(gradV+(grad{V}^T)\right)}^2$. We will limit ourselves to the case ${\mu }_{\infty }=0$ in the following, because we are dealing with infinite shear rate viscosity difficulties. We must assume that $\Gamma \overline{\dot{\gamma }}<1$, because the fluid under investigation is described as shear thinning.

Then, Equation (1) becomes

$$\begin{gathered} \overline{\tau }={\mu }_0\left[\left(\Gamma \overline{\dot{\gamma }}\right)\right]\overline{\dot{\gamma }}={\mu }_0\left[{\left(1+\Gamma \overline{\dot{\gamma }}-1\right)}^n\right]\overline{\dot{\gamma }}, \\ ={\mu }_0\left[1+n\left(\Gamma \overline{\dot{\gamma }}-1\right)\right]\overline{\dot{\gamma }}. \end{gathered}$$

(3)

The suggested model’s governing equations for continuity, energy, and momentum are as described in Ullah et al. [22].

In the equations below, subscripts denote PDE.

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

(4)

$${\rho }_{hnf}\left(\overline{\overline{u}}\frac{\partial \overline{\overline{u}}}{\partial \overline{\overline{x}}}+\overline{\overline{v}}\frac{\partial \overline{\overline{v}}}{\partial \overline{\overline{y}}}\right)={\mu }_{hnf}\frac{{\partial }^2\overline{\overline{u}}}{\partial {\overline{\overline{y}}}^2}-\sigma B^2\overline{\overline{u}}$$

(5)

$${\left(\rho c_p\right)}_{hnf}\left(\overline{\overline{u}}\frac{\partial T}{\partial \overline{\overline{x}}}+\overline{\overline{v}}\frac{\partial T}{\partial \overline{\overline{y}}}\right)=K_{hnf}{T}_{\overline{\overline{y}}\overline{\overline{y}}}\frac{{\partial }^{2}T}{\partial {\overline{\overline{y}}}^2}+{\theta }_0\left(T-T_{\infty }\right)$$

(6)

The velocity ingredient $u,v$ is in the direction $x,y$, while the symbols $\upsilon$ and $\rho$ stand for fluid density and kinematic fluid viscosity, respectively. The term $\sigma$ denotes the fluid’s electrical conductivity, i.e., the applied uniform magnetic field. T and $T_{\infty }$ stand for free stream temperature and temperature, respectively. $c_p$, $k$, and $Q_0$ stand for thermal conductivity, the specific heat of the fluid, and volumetric appraise of generating heat, respectively.

The temperature and velocity components associative boundary conditions are given by

$$T=T_w,\overline{\overline{v}}=f_0,\overline{\overline{u}}=a\overline{\overline{x}},\mathrm{for} \overline{\overline{y}}=0.$$

(7)

$$T\to T_{\infty },\overline{\overline{u}}\to 0, \mathrm{for} \overline{\overline{y}}\to \infty .$$

(8)

Here, the stretching rate is $a$.

The first step is to prepare a non-dimensional transform of the system that we have been provided. We introduce the dimensionless quantities below for this purpose.

$$\theta =\frac{T-T_{\infty }}{T_w-T_{\infty }},\overline{\overline{u}}=\overline{u}\sqrt{a\nu },\overline{\overline{v}}=\overline{v}\sqrt{a\nu },\overline{\overline{y}}=\sqrt{\frac{\nu }{a}}\overline{y},\overline{\overline{x}}=\sqrt{\frac{\nu }{a}}\overline{x}.$$

(9)

When the momentum, continuity, and energy equations are applied to the system given by Equations (4)–(6) and the bars are dropped, the continuity, speed, and energy equations become

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

(10)

$${\rho }_{hnf}\left(\overline{u}\frac{\partial \overline{u}}{\partial \overline{x}}+\overline{v}\frac{\partial \overline{v}}{\partial \overline{y}}\right)={\mu }_{hnf}\frac{{\partial }^2\overline{u}}{\partial {\overline{y}}^2}-\frac{\sigma B^2}{a}\overline{u}$$

(11)

$${\left(\rho c_p\right)}_{hnf}\left(\overline{u}\frac{\partial T}{\partial \overline{x}}+\overline{v}\frac{\partial T}{\partial \overline{y}}\right)=K_{hnf}\frac{{\partial }^{2}T}{\partial {\overline{y}}^2}+{\theta }_0\left(T-T_{\infty }\right)$$

(12)

The boundary conditions (7) and (8) in the scaling scenario described in (9) are as follows:

$$\overline{u}=\overline{x},\overline{v}=\frac{f_0}{\sqrt{\nu a}}, \theta =1 \mathrm{for}\overline{y}=0.$$

(13)

$$\theta \to 0,\overline{u}\to \infty \mathrm{for}\overline{y}\to \infty .$$

(14)

The stream function $\overline{v}=-\frac{\partial \psi }{\partial \overline{x}}$, $\overline{u}=\frac{\partial \psi }{\partial \overline{y}}$ is then used to decrease the dependent variables and equations. The continuity Equation (10) is naturally satisfied, as are Equations (12) and (11), which are expressed in terms of the stream function $\psi$, and have the form

$${\rho }_{hnf}\left(\frac{{\partial }^2\psi }{\partial \overline{x}\partial \overline{y}}\frac{\partial \psi }{\partial \overline{y}}-\frac{{\partial }^2\psi }{\partial {\overline{y}}^2}\frac{\partial \psi }{\partial \overline{x}}\right)={\mu }_{hnf}\frac{{\partial }^3\psi }{\partial {\overline{y}}^3}-\frac{\sigma B^2}{a}\frac{\partial \psi }{\partial \overline{y}}$$

(15)

$$\left(\frac{\partial \psi }{\partial \overline{y}}\frac{\partial \theta }{\partial \overline{x}}-\frac{\partial \psi }{\partial \overline{x}}\frac{\partial \theta }{\partial \overline{y}}\right)=\frac{K_{hnf}}{{\left(\mu c_p\right)}_{hnf}}\frac{{\partial }^2\theta }{\partial {\overline{y}}^2}+\frac{Q_0\theta }{{\left(\rho c_p\right)}_{hnf}a}$$

(16)

The stream function’s induction translates the boundary conditions (13) and (14) to

$$\theta =1, \frac{\partial \psi }{\partial \overline{x}}=-\frac{f_0}{\sqrt{\nu a}}, \frac{\partial \psi }{\partial \overline{y}}=\overline{x} \mathrm{for} \overline{y}=0.$$

(17)

$$\theta \to 0, \frac{\partial \psi }{\partial \overline{y}}\to 0 \mathrm{for} \overline{y}\to \infty .$$

(18)

3. Analysis

We use Lie group analysis in this section to develop novel similarity conventions for Equations (15) and (16). The nonlinear PDE will be reduced to nonlinear ODE. We evaluate the following transformation scaling group for this purpose.

$$\Gamma :{\Gamma }^{*}=\Gamma e^{\epsilon {\gamma }_5},{\theta }^{*}=\theta e^{\epsilon {\gamma }_4},{\psi }^{*}=\psi e^{\epsilon {\gamma }_3},{\overline{x}}^{*}=\overline{x}{e}^{\epsilon {\gamma }_1},{\overline{y}}^{*}=\overline{y}{e}^{\epsilon {\gamma }_2}.$$

(19)

Then, Group $\Gamma$ and ${\gamma }_i$, $\left(i=1,2,3,4,5\right)$ in terms of the parameter, and $\epsilon$ are the actual numbers to be determined. The coordinates are transformed by point transformation $\left(\overline{x},\overline{y},\psi ,\theta ,\Gamma \right)$ to $\left({\overline{x}}^{*},{\overline{y}}^{*},{\psi }^{*},{\theta }^{*},{\Gamma }^{*}\right)$ by (19).

By plugging (19) into Equations (15) and (16), we obtain

$${\rho }_{hnf}{e}^{\epsilon \left({\gamma }_1+2{\gamma }_2-2{\gamma }_3\right)}\left(\frac{\partial {\psi }^{\ast }}{\partial {\overline{y}}^{\ast }}\frac{{\partial }^2{\psi }^{*}}{\partial {\overline{x}}^{*}\partial {\overline{y}}^{*}}-\frac{\partial {\psi }^{*}}{\partial {\overline{x}}^{*}}\frac{{\partial }^2{\psi }^{*}}{\partial {\overline{y}}^{{*}^2}}\right)={\mu }_{hnf}{e}^{\epsilon \left(3{\gamma }_2-{\gamma }_3\right)}\frac{{\partial }^3\psi }{\partial {\overline{y}}^{{*}^3}}-e^{\epsilon \left({\gamma }_2-{\gamma }_3\right)}\frac{\sigma B^2}{a}\frac{\partial {\psi }^{*}}{\partial {\overline{y}}^{*}}$$

(20)

$$e^{\epsilon \left({\gamma }_1+{\gamma }_2-{\gamma }_3-{\gamma }_4\right)}\left(\frac{\partial {\psi }^{*}}{\partial {\overline{y}}^{*}}\frac{\partial {\theta }^{*}}{\partial {\overline{x}}^{*}}-\frac{\partial {\psi }^{*}}{\partial {\overline{x}}^{*}}\frac{\partial {\theta }^{*}}{\partial {\overline{y}}^{*}}\right)=\frac{K_{hnf}}{{\left(\mu c_p\right)}_{hnf}}{e}^{\epsilon \left(2{\gamma }_2-{\gamma }_4\right)}\frac{{\partial }^2\theta }{\partial {\overline{y}}^{{*}^2}}+e^{-\epsilon {\gamma }_4}\frac{Q_0{\theta }^{*}}{{\left(\rho c_p\right)}_{hnf}a}$$

(21)

If the coefficients of the above equations are equal, the changed systems (20) and (21) will remain invariant under the group of transformations.

$${\gamma }_1+2{\gamma }_2-2{\gamma }_3=3{\gamma }_2-{\gamma }_3={\gamma }_2-{\gamma }_3$$

(22)

$${\gamma }_1+{\gamma }_2-{\gamma }_3-{\gamma }_4=2{\gamma }_2-{\gamma }_4=-{\gamma }_4$$

(23)

We have, from the border condition, the following:

$${\gamma }_4=0$$

(24)

To obtain, we solve Equations (22) and (23) together.

$${\gamma }_1={\gamma }_3,{\gamma }_2=0,{\gamma }_3={\gamma }_1,{\gamma }_4=0.$$

(25)

The transformations are transformed into the following one-parameter group of changes by embedding (25) into the scaling (19).

$$\Gamma :{\overline{x}}^{*}=\overline{x}{e}^{\epsilon {\gamma }_1},{\overline{y}}^{*}=\overline{y},{\psi }^{*}=\psi e^{\epsilon {\gamma }_1},\theta =0.$$

(26)

The following simplified version is obtained by expanding the one-parameter group of (26) with Taylor’s series and maintaining terms up to first-order $\epsilon$.

$${\theta }^{*}-\theta =0,{\psi }^{*}-\psi =\overline{x}\epsilon {\gamma }_1,{\overline{y}}^{*}-\overline{y}=0,{\overline{x}}^{*}-\overline{x}=\overline{x}\epsilon {\gamma }_1.$$

(27)

The collection of transformations can readily be represented as characteristic equations using Equation (27).

$$\frac{d\overline{x}}{\overline{x}{\gamma }_1}=\frac{d\overline{y}}{0}=\frac{d\psi }{\overline{x}{\gamma }_1}=\frac{d\theta }{0}$$

(28)

The similarity transformations can be found in Equation (28). We obtain $\frac{d\overline{x}}{\overline{x}{\gamma }_1}=\frac{d\overline{y}}{0}$ from the first two terms in Equation (28), which we can integrate to obtain

$$\overline{y}=\mathrm{constant}=\zeta \left(\mathrm{say}\right).$$

(29)

From Equation (28), using 1,3 terms yields $\frac{d\overline{x}}{\overline{x}{\gamma }_1}=\frac{d\psi }{\overline{x}{\gamma }_1}$.

$$\frac{\psi }{\overline{x}}=\mathrm{constant}=f\left(\zeta \right)\left(\mathrm{say}\right), \mathrm{then} \psi =\overline{x}f\left(\zeta \right).$$

(30)

By integrating both sides of Equation (28) and equating the first and fourth terms, we obtain

$$\theta =\theta \left(\zeta \right).$$

(31)

As a result, the new similarity transformations can be written as

$$\zeta =\overline{y},\psi =\overline{x}f\left(\zeta \right),\theta =\theta \left(\zeta \right).$$

(32)

We obtain the following nonlinear ODE by substituting Equation (32) into Equations (15) and (16):

$$\frac{{\mu }_{hnf}}{{\rho }_{hnf}}{f}^{‴}\left(\zeta \right)-M^{2}f\left(\zeta \right)+f\left(\zeta \right)f^{″}\left(\zeta \right)-{\left(f^{\prime }\left(\zeta \right)\right)}^2=0$$

(33)

$$\left(\frac{k_{hnf}}{{\left(\mu c_p\right)}_{hnf}}\right){\theta }^{″}\left(\zeta \right)+f\left(\zeta \right){\theta }^{\prime }\left(\zeta \right)+Q\theta \left(\zeta \right)=0$$

(34)

The Hartmann number is $M^2=\frac{\sigma B^2}{\rho a}$, the heat source/sink parameter is $Q=\frac{Q_0}{\rho c_{p}a}$, and the Prandtl number is $\mathrm{Pr}=\frac{\mu c_p}{k}$.

Primes are used to represent differentials using $\zeta$. Subject to the boundary conditions, we solve the system by Equations (33) and (34).

$$\theta \left(0\right)=1,f\left(0\right)=f_w,f^{\prime }\left(0\right)=1. \mathrm{for} \zeta \to 0.$$

(35)

$$\theta (\infty )\to 0,f^{\prime }(\infty )\to 0 \mathrm{at} \zeta \to \infty .$$

(36)

where $f\left(\zeta \right)=f_w,f_w=\frac{f_0}{\sqrt{\nu a}}$, $f_w\succ 0$ corresponds to suction-permeable, and $f_w\prec 0$ corresponds to injection-permeable.

The skin friction coefficient $C_f$ and the local Nusselt number $N{u}_x$ are defined as follows:

$$C_f=\frac{{\tau }_w}{\rho {\left(a\overline{x}\right)}^2}, N{u}_{\overline{x}}=\frac{\overline{x}{q}_w}{k\left(T_w-T_{\infty }\right)}$$

(37)

where ${\tau }_w$ represents skin friction, and $q_w$ represents heat transfer from the plate:

$${\tau }_w={\mu }_{hnf}{\left(\frac{\partial \overline{u}}{\partial \overline{y}}\right)}_{\zeta =0}, q_w=-k_{hnf}{\left(\frac{\partial T}{\partial \overline{y}}\right)}_{\zeta =0}.$$

(38)

When we plug Equations (9) and (32) into (37), the skin friction and local Nusselt become dimensionless:

$${\mathrm{Re}}^{\frac{1}{2}}Cf=\frac{{\mu }_{hnf}}{{\mu }_f}{f}^{″}\left(0\right), {\mathrm{Re}}^{-\frac{1}{2}}N{u}_{\overline{x}}=-\frac{k_{hnf}}{k_f}{\theta }^{\prime }\left(0\right).$$

where $\mathrm{Re}$ is the Reynolds number.

Next, we apply the hybrid nanofluid thermophysical properties.

${\mathsf{\phi }}_1,{\mathsf{\phi }}_2$, and ${\mathsf{\phi }}_3$ represent the volume fraction of solid particles for Case 1 (Platelet $\left(A{l}_2{O}_3\right)$, Cylindrical (GNT), and Spherical (CNT)) and solid particles for Case 2 (Spherical (copper), Cylindrical (silver), Platelet (cobalt ferrite)) nanoparticles. The ternary hybrid nanoparticles with different shapes have platelet, cylindrical and spherical shapes.

Thermal conductivity and viscosity are

$${\mu }_{hnf}={\phi }^{-1}({\mu }_{nf3}{\mathsf{\phi }}_3+{\mu }_{nf2}{\mathsf{\phi }}_2+{\mu }_{nf1}{\mathsf{\phi }}_1)$$

(39)

$$k_{hnf}={\phi }^{-1}(k_{nf3}{\mathsf{\phi }}_3+k_{nf2}{\mathsf{\phi }}_2+k_{nf1}{\mathsf{\phi }}_1)$$

(40)

The ${\rho }_{hnf}$ density of ternary hybrid nanoparticles (platelet, cylindrical, and spherical form) is given by

$${\rho }_{hnf}=(1-{\mathsf{\phi }}_1-{\mathsf{\phi }}_2-{\mathsf{\phi }}_3){\rho }_{bf}+{\mathsf{\phi }}_1{\rho }_{sp1}+{\mathsf{\phi }}_2{\rho }_{sp2}+{\mathsf{\phi }}_3{\rho }_{sp3}$$

(41)

The ${(\rho c_p)}_{hnf}$ heat capacity of ternary hybrid nanoparticles (platelet, cylindrical, spherical, and shape) is calculated as follows:

$${(\rho c_p)}_{hnf}=(1-{\mathsf{\phi }}_1-{\mathsf{\phi }}_2-{\mathsf{\phi }}_3){(\rho c_p)}_{bf}+{\mathsf{\phi }}_1{(\rho c_p)}_{sp1}+{\mathsf{\phi }}_2{(\rho c_p)}_{sp2}+{\mathsf{\phi }}_3{(\rho c_p)}_{sp3}$$

(42)

Thermal conductivity and viscosity for spherical nanoparticles are expressed by

$${\mu }_{nf1}=\left({\mu }_{bf}\right)\left(1+2.5\mathsf{\phi }+6.2{\mathsf{\phi }}^2\right)\hspace{0.17em}$$

(43)

$$\frac{k_{nf1}}{k_{bf}}=\left[\frac{k_{sp1}+2k_{bf}-2\mathsf{\phi }(k_{bf}-k_{sp1})}{k_{sp1}+2k_{bf}+\mathsf{\phi }(k_{bf}-k_{sp1})}\right]$$

(44)

For cylindrical nanoparticles, the thermal conductivity and viscosity are expressed by

$$\frac{{\mu }_{nf2}}{{\mu }_{bf}}=\left(1+13.5\mathsf{\phi }+904.4{\mathsf{\phi }}^2\right)\hspace{0.17em}$$

(45)

$$\frac{k_{nf2}}{k_{bf}}=\left[\frac{k_{sp2}+3.9k_{bf}-3.9\mathsf{\phi }(k_{bf}-k_{sp2})}{k_{sp2}+3.9k_{bf}+\mathsf{\phi }(k_{bf}-k_{sp2})}\right]$$

(46)

Platelet nanoparticles have the following thermal conductivity and thermal viscosity:

$$\frac{{\mu }_{nf3}}{{\mu }_{bf}}=\left(1+37.1\mathsf{\phi }+612.6{\mathsf{\phi }}^2\right)\hspace{0.17em}$$

(47)

$$\frac{k_{nf3}}{k_{bf}}=\left[\frac{k_{sp3}+4.7k_{bf}-4.7\mathsf{\phi }(k_{bf}-k_{sp3})}{k_{sp3}+4.7k_{bf}+\mathsf{\phi }(k_{bf}-k_{sp3})}\right]$$

(48)

The above properties for Equations (33) and (34) are as follows:

$$\left(\frac{A_1}{A_2\mathsf{\phi }}\right)f^{‴}\left(\zeta \right)-{\left(f^{\prime }\left(\zeta \right)\right)}^2-M^{2}f\left(\zeta \right)+f\left(\zeta \right)f^{″}\left(\zeta \right)=0$$

(49)

$$\left(\frac{A_4}{A_3}\right)\left(\mathrm{Pr}\mathsf{\phi }\right){\theta }^{″}\left(\zeta \right)+f\left(\zeta \right){\theta }^{\prime }\left(\zeta \right)+Q\theta \left(\zeta \right)=0$$

(50)

In the above equations, the transformations are as follows:

$$A_1=B_3{\mathsf{\phi }}_3+B_2{\mathsf{\phi }}_2+B_1{\mathsf{\phi }}_1$$

$$A_2=1-{\mathsf{\phi }}_1-{\mathsf{\phi }}_2-{\mathsf{\phi }}_3+{\mathsf{\phi }}_1\frac{{\rho }_{sp1}}{{\rho }_{bf}}+{\mathsf{\phi }}_2\frac{{\rho }_{sp2}}{{\rho }_{bf}}+{\mathsf{\phi }}_3\frac{{\rho }_{sp3}}{{\rho }_{bf}}$$

$$A_3=B_4{\mathsf{\phi }}_1+B_5{\mathsf{\phi }}_2+B_6{\mathsf{\phi }}_3$$

$$A_4=1-{\mathsf{\phi }}_1-{\mathsf{\phi }}_2-{\mathsf{\phi }}_3+{\mathsf{\phi }}_1\frac{{(\rho c_p)}_{sp1}}{{(\rho c_p)}_{bf}}+{\mathsf{\phi }}_2\frac{{(\rho c_p)}_{sp2}}{{(\rho c_p)}_{bf}}+{\mathsf{\phi }}_3\frac{{(\rho c_p)}_{sp3}}{{(\rho c_p)}_{bf}}$$

$$A_5=\left(1-{\mathsf{\phi }}_1-{\mathsf{\phi }}_2-{\mathsf{\phi }}_3\right)+{\mathsf{\phi }}_1\frac{{\left(\rho {\beta }_0\right)}_{sp1}}{{\left(\rho {\beta }_0\right)}_{bf}}+{\mathsf{\phi }}_2\frac{{\left(\rho {\beta }_0\right)}_{sp2}}{{\left(\rho {\beta }_0\right)}_{bf}}+{\mathsf{\phi }}_3\frac{{\left(\rho {\beta }_0\right)}_{sp3}}{{\left(\rho {\beta }_0\right)}_{bf}}$$

$$A_6=\left(1-{\mathsf{\phi }}_1-{\mathsf{\phi }}_2-{\mathsf{\phi }}_3\right)+{\mathsf{\phi }}_1\frac{{\left(\rho {\beta }_1\right)}_{sp1}}{{\left(\rho {\beta }_1\right)}_{bf}}+{\mathsf{\phi }}_2\frac{{\left(\rho {\beta }_1\right)}_{sp2}}{{\left(\rho {\beta }_1\right)}_{bf}}+{\mathsf{\phi }}_3\frac{{\left(\rho {\beta }_1\right)}_{sp3}}{{\left(\rho {\beta }_1\right)}_{bf}}$$

$$B_1=1+2.5\hspace{0.17em}\mathsf{\phi }+6.2{\mathsf{\phi }}^2$$

$$B_2=1+13.5\hspace{0.17em}\mathsf{\phi }+904.4{\mathsf{\phi }}^2$$

$$B_3=1+37.1\hspace{0.17em}\mathsf{\phi }+612.6{\mathsf{\phi }}^2$$

$$B_4=\frac{k_{sp1}+2k_{bf}-2\mathsf{\phi }(k_{bf}-k_{sp1})}{k_{sp1}+2k_{bf}+\mathsf{\phi }(k_{bf}-k_{sp1})}$$

$$B_5=\frac{k_{sp2}+3.9k_{bf}-3.9\mathsf{\phi }(k_{bf}-k_{sp2})}{k_{sp2}+3.9k_{bf}+\mathsf{\phi }(k_{bf}-k_{sp2})}$$

$$B_6=\frac{k_{sp3}+4.7k_{bf}-4.7\mathsf{\phi }(k_{bf}-k_{sp3})}{k_{sp3}+4.7k_{bf}+\mathsf{\phi }(k_{bf}-k_{sp3})}$$

4. Numerical Procedure

The changed Equations (49) and (50) with the given circumstances can be solved using bvp4c (35–36), a built-in function in MATLAB. In this problem, we used the bvp4c solver because it is a built-in function, informed by Mamatha and Raju [23].

We utilized the bvp4c solver via the Runge–Kutta approach.

We used the following assumptions as a preprocedural before coding:

$$f=J_{1,}f\prime =J_{2,}f″=J_{3,}\theta =J_4,\theta \prime =J_{5.}$$

Then, using the Equations (49) and (50) and conditions (35–36), we may construct a first-order system of ODEs:

$$\begin{array}{l}{J}_1^{\prime }=J_2\\ J^{\prime }{}_2=J_3\\ J^{\prime }{}_3=\mathsf{\phi }\left(\frac{A_2}{A_1}\right)\left[{\left(J^{\prime }\left(\zeta \right)\right)}^2-J\left(\zeta \right)J^{″}\left(\zeta \right)+M^2{J}^{\prime }\left(\zeta \right)\right]\\ J_4^{\prime }=J_5\\ J_5^{\prime }=\mathrm{Pr}\mathsf{\phi }\left(\frac{A_4}{A_3}\right)\left[-J\left(\zeta \right)\theta \prime \left(\zeta \right)-Q\theta \left(\zeta \right)\right]\end{array}\}$$

with the circumstances

$$\begin{array}{l}{f}_a\left(1\right)=f_w\\ f_a\left(2\right)=1\\ f_a\left(4\right)=1\\ f_b\left(2\right)=0\\ f_b\left(4\right)=0\end{array}\}$$

5. Results and Discussion

The obtained findings show how non-dimensional parameters like Prandtl number ($\mathrm{Pr}$), Hartmann number ($M$), heat sink/source ($Q$), and suction ($f_w$) affect non-dimensional velocity and temperature profiles (see Figure 1, Figure 2, Figure 3, Figure 4, Figure 5, Figure 6, Figure 7, Figure 8, Figure 9 and Figure 10). $f_w=0.2,\mathrm{Pr}=6.2,M=0.5,Q=0.5$ are the same parameter values used in the simulations; these numbers were consistent throughout the inquiry, except for the various discounts given in the tables and figures. In

Table 1. Thermophysical properties of a hybrid nanofluid.

	Nomenclature of Solid Particles and Base Fluid	ρ (kg/m3)	Cp (J/kgK)	K (W/mK)	Nanoparticle Shapes
Base fluid	Water H2O	997.1	4.179	0.623
Ternary hybrid nanofluid 1	Graphene (1%)	2200	5000	790	Platelet
	Carbon nanotubes (1%)	5100	410	3007	Cylindrical
	Aluminum oxide (Al2O3) (1%)	3970	765	40	Spherical
Ternary hybrid nanofluid 2	Copper (1%)	10,500	235	429	Spherical
	Silver (1%)	8933	385	400	Cylindrical
	Cobalt ferrite (1%)	4907	700	3.7	Platelet

mentioned the thermophysical properties of nanoparticles for each case (Case 1: aluminum oxide/CNT/graphene; Case 2: copper/silver/cobalt ferrite), using water as the base fluid with different nanoparticle shapes such as spherical, cylindrical, and platelet.

Table 2. The values of the skin friction coefficient and the local Nusselt number as a physical parameter.

					Skin Friction		Nusselt Number
Q	M	$\mathit{\phi }$	Pr	fw	Case-1	Case-2	Case-1	Case-2
0.2					2.402399	2.540362	3.314880	3.292131
0.4					2.402398	2.540362	2.613102	2.575383
0.6					2.402399	2.540362	1.709249	1.639159
0.8					2.402399	2.540362	0.366501	0.198966
	0.6				2.489423	2.633932	2.160871	2.106608
	1.2				3.239851	3.436586	1.827122	1.732227
	1.8				4.220439	4.478802	1.252630	1.070368
	2.4				5.297751	5.621093	0.303402	0.084177
		0.002			6.215759	6.569464	3.003631	2.818563
		0.004			4.367599	4.631485	2.512344	2.370365
		0.006			3.578756	3.798511	2.284544	2.174230
		0.008			3.112536	3.305193	2.142240	2.052131
			0.7		2.797875	2.971118	0.184786	0.159446
			1.4		2.797877	2.971122	0.025169	0.117266
			2.1		2.797877	2.971122	0.094425	0.065792
			2.8		2.797877	2.971122	0.456514	0.289438
				1	3.295625	3.466278	2.173269	2.172145
				2	3.959843	4.225607	4.089725	4.095232
				3	4.702199	5.079428	6.032918	6.044240
				4	5.513988	6.015086	7.987920	8.004670

shows skin friction and the Nusselt number for the local area.

Figure 1, Figure 2, Figure 3, Figure 4, Figure 5, Figure 6, Figure 7, Figure 8, Figure 9 and Figure 10 depict the differences in the profiles for various $f_w=0.2,\mathrm{Pr}=6.2,M=0.5,Q=0.3.$ values. From Figure 1 and Figure 2, we can observe the effect of the suction on temperature and velocity profiles: when the suction rises, the temperature and velocity of both profiles increase. From Figure 3 and Figure 4, we can observe the effect of the Hartmann number on temperature and velocity profiles: when the Hartmann number rises, the velocity increases, and the temperature profile decreases. From Figure 5 and Figure 6, we can observe the effect of the volume fraction on temperature and velocity profiles: when the volume fraction rises, the temperature and velocity of both profiles increase. From Figure 7 and Figure 8, we can observe the effect of the Prandtl number on temperature and velocity profiles: when the Prandtl number rises, the temperature and velocity of both profiles increase. From Figure 9 and Figure 10, we can observe the effect of the heat source/sink on temperature and velocity profiles: when the heat source/sink rises, the temperature and velocity of both profiles decrease.

From Table 2, when the heat source/sink rises, skin friction is steady in both Cases 1 and 2, and we can observe that Case 2 has a greater skin friction transfer rate than Case 1. When the heat source/sink rises, the Nusselt number transfer rate decreases in both Cases 1 and 2, and we can observe that Case 1 has a greater Nusselt number transfer rate than Case 2.

When the Hartmann number rises, the skin friction transfer rate decreases in both Cases 1 and 2, and we can see that Case 2 has a greater skin friction transfer rate than Case 1. When the Hartmann number increases, the Nusselt number transfer rate decreases in both Cases 1 and 2, and we can observe that Case 1 has a greater Nusselt number transfer rate than Case 2.

When the volume fraction rises, the skin friction transfer rate decreases in both Cases 1 and 2, and we can see that Case 2 has a greater skin friction transfer rate than Case-1. When the volume fraction increases, the Nusselt number transfer rate decreases in both Cases 1 and 2, and we can observe that Case 1 has a greater Nusselt number transfer rate than Case 2.

When the Prandtl number rises, the skin friction transfer rate is steady in both Cases 1 and 2, and we can see that Case 2 has a greater skin friction transfer rate than Case 1. When the Prandtl number increases, the Nusselt number transfer rate has mixed behavior in both Cases 1 and 2, even though we can observe that Case 1 has a greater Nusselt number transfer rate compared to Case 2.

When suction values rise in both Cases 1 and 2, the skin friction transfer rate also increases, and we can observe that Case 2 has a greater skin friction transfer rate than Case 1. When suction values rise, the Nusselt number transfer rate also grows in both Cases 1 and 2, and we can see that Case 2 has a greater Nusselt number transfer rate than Case 1.

It can be observed that, in the coefficients of the models, Case 1 $M$ and $\mathrm{Pr}$ have a positive impact on $Nu{s}_1$, but the reverse was true with $f_w$; in Case 2, $Nu{s}_2$ showed the same trend. Surfaces were used to illustrate the interactive influence of the important factors on the response variable, as seen in Figure 11, Figure 12, Figure 13, Figure 14, Figure 15 and Figure 16. In Case 1, low-level values of $M$ and high-level values of $\mathrm{Pr}$ enhance the heat transfer rate, medium values of $\mathrm{Pr}$ and high-level values of $f_w$ also enhance the heat transfer rate, but low-level values of $M$ and $f_w$ decrement the heat transfer rate.

6. RSM: Response Surface Methodology

To review, Case 1 consisted of aluminum oxide, CNT, and graphene particles, while Case 2 consisted of copper, silver, and cobalt ferrite particles. In both cases, particles exhibited different shapes such as sphere, cylindrical, and platelet.

The RSM was utilized to investigate the continual impacts of important parameters (i.e., independent parameters) on the dependent-variable response variables and to identify the key parameter’s ideal levels. Additionally, the approach is consistent with an observed relationship for the (independent-variable) response variables, which is useful for forecasting and optimization. Due to its computing efficiency, this approach was first introduced by Box and Wilson [24] and has since been applied in other disciplines. For accurate modeling, a suitable experimental design was selected. The experimental design made use of the (FC-CCD) face-cantered central composite design described in [25]. A fractional factorial or factorial design with center points, enhanced with a collection of axial points that permit curvature estimation, is known as a central composite design.

Central composite design (CCD) can be used to estimate 1st- and 2nd-order terms and to model a (dependent-variable) response variable with curvature by including center and axial points to a prior factorial design. Having more degrees of freedom and less error than other designs, the FC-CCD design is preferable. Also taken into account in this design are the high levels of the factors. Using $F+2N+2^N$, where $F$ is the number of faces in the design and $N$ is the number of independent variables, the number of runs for the experimental design using CCD can be obtained. The following equation gives the complete quadratic model for three variables that include the interactional, linear, and square terms of the important factors:

$$Z=b_0+b_1{Y}_1+b_2{Y}_2+b_3{Y}_3+b_4{Y}_1{Y}_2+b_5{Y}_1{Y}_3+b_6{Y}_2{Y}_3+b_7{Y}_1^2+b_8{Y}_2^2+b_9{Y}_3^2$$

(51)

where the regression coefficients $b_i$ are indicated. The RSM analysis, which is focused on the goal of maximizing attractiveness, is used to optimize multiple replies. The value of desirability ranges from 0 to 1. Although greater desirability values are preferable, simultaneous optimization of various factors may result in lower desirability values. The following relationship is used to determine the desirability, as described in [26]:

$$\begin{array}{ll}d=0,& \mathrm{Z}\le {\mathrm{Z}}_L\\ d=\left[\frac{Z-Z_L}{Z_H-Z_L}\right],& {\mathrm{Z}}_L\le Z\le Z_H\\ d=1,& \mathrm{Z}\ge {\mathrm{Z}}_H\end{array}$$

(52)

Here, $Z_{H,}$ and $Z_L$ stand for, respectively, the highest and lowest response value.

6.1. RSM and Optimization Outcomes

The study’s three main components include the Hartmann number ($M$), Prandtl number ($\mathrm{Pr}$), and a suction parameter ($f_w$) to investigate the ideal amounts of these elements to improve the rate of heat transmission. They are chosen based on three levels (medium (0), high (+1), low (−1)); the key variables, their levels, and their codes values are listed in

Table 3. Key RSM parameters, their levels, and their symbols.

Key Factors	Symbols	Levels
		−1 (Low)	0 (Medium)	1 (High)
$M$	$X_1$	1.2	1.8	2.4
$\mathrm{Pr}$	$X_2$	1.4	2.1	2.8
$f_w$	$X_3$	2	3	4

. In accordance with

Table 4. Design and results of heat transfer rate experiments (Nur).

Runs	Coded Values			Real Values			Response
	${\mathit{X}}_1$	${\mathit{X}}_2$	${\mathit{X}}_3$	$\mathit{M}$	$\mathbf{Pr}$	${\mathit{f}}_{\mathit{w}}$	Nus-1	Nus-2
1	−1	−1	−1	1.2	1.4	2	2.516313	2.716616
2	1	−1	−1	2.4	1.4	2	2.701862	2.534297
3	−1	1	−1	1.2	2.8	2	5.711829	5.744885
4	1	1	−1	2.4	2.8	2	5.567382	5.602613
5	−1	−1	1	1.2	1.4	4	5.550822	5.586327
6	1	−1	1	2.4	1.4	4	5.495897	5.532817
7	−1	1	1	1.2	2.8	4	11.197125	11.270286
8	1	1	1	2.4	2.8	4	11.15152	11.225642
9	−1	0	0	1.2	2.1	3	6.321479	6.360699
10	1	0	0	2.4	2.1	3	6.238506	6.279378
11	0	−1	0	1.8	1.4	3	4.103701	4.130272
12	0	1	0	1.8	2.8	3	8.420175	8.474129
13	0	0	−1	1.8	2.1	2	5.639609	5.673726
14	0	0	1	1.8	2.1	4	11.175048	11.248695
15	0	0	0	1.8	2.1	3	8.420175	8.474129
16	0	0	0	1.8	2.1	3	8.420175	8.474129
17	0	0	0	1.8	2.1	3	8.420175	8.474129
18	0	0	0	1.8	2.1	3	8.420175	8.474129
19	0	0	0	1.8	2.1	3	8.420175	8.474129
20	0	0	0	1.8	2.1	3	8.420175	8.474129

, this design requires 20 runs, using the formula $F+2N+2^N$.

Analysis of variance, shown in

Table 5. For Nus, the variance was examined.

Source	Degrees of Freedom	Adjusted Sum of Squares		Adjusted Mean Square		F-Value		p-Value
		Case-1	Case-2	Case-1	Case-2	Case-1	Case-2	Case-1	Case-2
Model	9	121.385	123.026	13.4872	13.6695	50.34	50.64	0	0
Linear	3	12.507	11.961	4.169	3.9869	15.56	14.77	0	0
$M$	1	4.524	4.248	4.5241	4.2476	16.89	15.73	0.002	0.003
$\mathrm{Pr}$	1	4.823	4.726	4.8234	4.7261	18	17.51	0.002	0.002
$f_w$	1	1.161	1.254	1.161	1.2542	4.33	4.65	0.064	0.057
Square	3	20.608	20.871	6.8695	6.9571	25.64	25.77	0	0
$M\ast M$	1	5.135	5.202	5.135	5.202	19.17	19.27	0.001	0.001
$\mathrm{Pr}\ast \mathrm{Pr}$	1	5.272	5.338	5.2716	5.3378	19.67	19.77	0.001	0.001
$f_w\ast f_w$	1	1.592	1.613	1.592	1.6128	5.94	5.97	0.035	0.035
2-Way Interaction	3	3.449	3.492	1.1496	1.1639	4.29	4.31	0.034	0.034
$M\ast \mathrm{Pr}$	1	0.013	0	0.0129	0.0003	0.05	0	0.831	0.974
$\mathrm{Pr}\ast f_w$	1	0.003	0.006	0.0025	0.0064	0.01	0.002	0.925	0.881
$f_w\ast M$	1	3.433	3.485	3.4334	3.4851	12.81	12.91	0.005	0.005
Error	10	2.679	2.7	0.2679	0.27
Lack of Fit	5	2.679	2.7	0.5359	0.5399
Pure Error	5	0	0	0	0
Total	19	124.064	125.725
$R^2=$	97.85%

, is used to determine the precision of the estimated RSM model. The coefficient of determination ($R^2$ = 97.85%) quantifies the amount of variation in the (dependent-variable) response variable that is evaluated by the important variables. It also serves as a gauge for the model’s degree of accuracy. $R^2$ = 97.85% is considered high and indicates the precision of the model. The model’s correctness is further confirmed by the residual plots in Figure 17 and Figure 18. The data points in the normal probability plot along the straight line indicate the residuals’ normality, and the residuals’ histogram has a bell-shaped distribution. All of these results indicate the model’s correctness. According to the coded coefficients (X₁, X₂, and X₃), the quadratic model for $Nus$ for two cases is as follows, in accordance with the Taguchi model [27]:

$$\begin{array}{l}Nu{s}_1=-17.62+13.93M+12.33\mathrm{Pr}-4.23f_w-3.796M^2-2.826{\mathrm{Pr}}^2+0.761f_w^2\\ -0.095M\mathrm{Pr}-0.030M{f}_w+0.936\mathrm{Pr}{f}_w.\\ Nu{s}_2=-16.78+13.50M+12.20\mathrm{Pr}-4.40f_w-3.820M^2-2.843{\mathrm{Pr}}^2\\ +0.766f_w^2+0.015M\mathrm{Pr}+0.047M{f}_w+0.943\mathrm{Pr}{f}_w.\end{array}\}$$

(53)

In Case 2, low-level values of $M$ and high-level values of $\mathrm{Pr}$ enhance the heat transfer rate, medium-level values of $\mathrm{Pr}$ and high-level values $f_w$ also enhance the heat transfer rate, but low-level values of $M$ and $f_w$ decrement the heat transfer rate.

The curvature in the surface plots also amply demonstrates the non-linear link between the important elements. To increase process efficiency, it is crucial to identify the critical variables’ ideal values that maximize the rate of heat transfer. RSM estimates the maximal heat transfer at this ideal state to be $Nu{s}_1$ = 11.197125 and $Nu{s}_2$ = 11.270286.

6.2. Sensitivity Analysis

The sensitivity of the heat transmission rate is evaluated in this subsection. The sensitivity analysis gives details about how the primary factors’ increments affect the response variable. Designers and engineers can learn about the relative influence of the important elements by looking at the proportions of sensitivity values at particular value levels for the parameters. A positive correlation follows from a positive sensitivity, and the reverse from a negative sensitivity. The following calculations represent the heat transmission rate’s sensitivity functions:

$$\begin{array}{l}\frac{\partial Nus}{\partial M}=13.93-7.592\ast M-0.095\ast \mathrm{Pr}-0.030\ast fw\\ \frac{\partial Nus}{\partial \mathrm{Pr}}=12.33-5.652\ast \mathrm{Pr}-0.095\ast M+0.936\ast fw\\ \frac{\partial Nus}{\partial fw}=-4.23+1.522\ast fw-0.030\ast M+0.936\ast \mathrm{Pr}\end{array}\}$$

(54)

$$\begin{array}{l}\frac{\partial Nus}{\partial M}=13.50-7.64\ast M+0.015\ast \mathrm{Pr}+0.047\ast fw\\ \frac{\partial Nus}{\partial \mathrm{Pr}}=12.20-5.686\ast \mathrm{Pr}+0.015\ast M+0.943\ast fw\\ \frac{\partial Nus}{\partial fw}=-4.40+1.532\ast fw+0.047\ast M+0.943\ast \mathrm{Pr}\end{array}\}$$

(55)

The sensitivity values for each level of the significant components are listed in

Table 6. Sensitivity values at all levels of the key factors.

			$\frac{\partial \mathit{N}\mathit{u}\mathit{s}}{\partial \mathit{M}}$		$\frac{\partial \mathit{N}\mathit{u}\mathit{s}}{\partial \mathbf{Pr}}$		$\frac{\partial \mathit{N}\mathit{u}\mathit{s}}{\partial \mathit{f}\mathit{w}}$
${\mathit{X}}_1$	${\mathit{X}}_2$	${\mathit{X}}_3$	Case-1	Case-2	Case-1	Case-2	Case-1	Case-2
−1	−1	−1	4.6266	4.447	6.1752	6.1436	0.0884	0.0406
−1	−1	0	4.5966	4.494	7.1112	7.0866	1.6104	1.5726
−1	−1	1	4.5666	4.541	8.0742	8.0296	3.1324	3.1046
−1	0	−1	4.5601	4.4575	2.2188	2.1634	0.7436	0.7707
−1	0	0	4.5301	4.5045	3.1548	3.1064	2.2656	2.2327
−1	0	1	4.5001	4.5515	4.0908	4.0494	3.7876	3.7647
−1	1	−1	4.4936	4.468	−1.7376	−1.8168	1.3988	1.0608
−1	1	0	4.4636	4.515	−0.8016	−0.8738	2.9208	2.8928
−1	1	1	4.4336	4.562	0.1344	6.822	4.4428	4.4248
0	−1	−1	0.0714	−0.137	6.1182	6.1526	0.0704	0.0688
0	−1	0	0.0414	−0.09	7.0542	7.0956	1.5924	1.6008
0	−1	1	0.0114	−0.043	7.9902	8.0386	3.1144	3.1328
0	0	−1	0.0049	−0.1265	2.1618	2.1724	0.7256	0.7289
0	0	0	−0.0251	−0.0795	3.0978	3.1154	2.2476	2.2609
0	0	1	−0.0551	−0.0325	4.0338	4.0584	3.7696	3.7929
0	1	−1	−0.0616	−0.116	−1.7946	−1.8078	1.3808	1.389
0	1	0	−0.0916	−0.069	−0.8586	−0.8648	0.9028	2.921
0	1	1	−0.1216	−0.022	0.0774	0.0782	4.4248	4.453
1	−1	−1	−4.4838	−4.721	6.0612	6.1616	0.0524	0.097
1	−1	0	−4.5138	−4.674	6.9972	7.1046	1.5744	1.629
1	−1	1	−4.5438	−4.627	7.9332	8.0476	3.0964	3.161
1	0	−1	−4.5503	−4.7105	2.1048	2.1814	0.7076	0.7571
1	0	0	−4.5803	−4.6635	3.0408	3.1244	2.2296	2.2891
1	0	1	−4.6103	−4.6165	3.9768	4.0674	3.7516	3.8211
1	1	−1	−4.6168	−4.7	−1.8516	−1.7988	1.3628	1.4172
1	1	0	−4.6468	−4.653	−0.9516	−0.8558	2.8848	2.9492

. The Nusselt number’s sensitivity is shown in Table 6, where $\frac{\partial Nus}{\partial M}$ showed a mixed nature in Case 1 and Case 2, $\frac{\partial Nus}{\partial \mathrm{Pr}}$ showed a positive nature in Case 1 and Case 2, and $\frac{\partial Nus}{\partial fw}$ showed a positive nature in Case 1 and Case 2.

The following outcomes are produced by an increase in the key factor levels:

• Raising $M$ levels (from −1 to 1) results in decreased sensitivity for the Case 1 $\frac{\partial Nus}{\partial M}$ and the same sensitivity nature for Case 2.
• Raising $\mathrm{Pr}$ levels (from −1 to 1) results in mixed sensitivity in both Case 1 and Case 2.
• Raising $f_w$ levels (from −1 to 1) results in mixed sensitivity in Case 1 and increased sensitivity in Case 2.

7. Conclusions

Hybrid nanofluids are generally employed as coolants in heat transfer equipment such as electronic cooling, heat exchange systems (such as flat-panel systems), and radiators due to their improved thermal characteristics.

Here, the effects of MHD and suction were studied for two hybrid nanofluids (Case 1: CNT/aluminum oxide/graphene; Case 2: copper/silver/cobalt ferrite), using water as a base fluid and nanoparticles with different shapes such as spherical, cylindrical, and platelet. The effects of four parameters ($\mathrm{Pr},M,Q,f_w$) on heat and velocity transfer rate were also analyzed. The scaling group of transformations was used to convert from dimensional to non-dimensional equations to conduct the analysis.

This research yielded the following conclusions:

• When suction rises, the velocity and temperature profiles both increase.
• When the Hartmann number rises, velocity and temperature exhibit inverse behaviors: the velocity profile increases, while the temperature profile decreases.
• When the volume fraction rises, the temperature and velocity profiles both increase.
• When the Prandtl number rises, the temperature and velocity profiles both increase.
• When the heat source/sink rises, the temperature and velocity profiles both decrease.