Multiple-Relaxation-Time Lattice Boltzmann Simulation of Soret and Dufour Effects on the Thermosolutal Natural Convection of a Nanofluid in a U-Shaped Porous Enclosure

Abstract

This article reports an investigation of the Soret and Dufour effects on the double-diffusive natural convection of $A{l}_2{O}_3$-$H_{2}O$ nanofluids in a U-shaped porous enclosure. Numerical problems were resolved using the multiple-relaxation-time (MRT) lattice Boltzmann method (LBM). The indented part of the U-shape was cold, and the right and left walls were heated, while the bottom and upper walls were adiabatic. The experimental data-based temperature and nanoparticle size-dependent correlations for the $A{l}_2{O}_3$-water nanofluids are used here. The benchmark results thoroughly validate the graphics process unit (GPU) based in-house compute unified device architecture (CUDA) C/C++ code. Numeral simulations were performed for a variety of dimensionless variables, including the Rayleigh number, ($Ra$ = ${10}^4,{10}^5,{10}^6$), the Darcy number, ($Da$ = ${10}^{-2},{10}^{-3},{10}^{-4}$), the Soret number, ($Sr$ = $0.0,0.1,0.2$), the Dufour number, ($D_f$ = $0.0,0.1,0.2$), the buoyancy ratio, ($-2\le Br\le 2$), the Lewis number, ($Le$ = $1,3,5$), the volume fraction, ($0\le \varphi \le 0.04$), and the porosity, $ϵ$ = ($0.2-0.8$), and the Prandtl number, $Pr$ = $6.2$ (water) is fixed to represent the base fluid. The numerical results are presented in terms of streamlines, isotherms, isoconcentrations, temperature, velocity, mean Nusselt number, mean Sherwood number, entropy generation, and statistical analysis using a response surface methodology (RSM). The investigation found that fluid mobility was enhanced as the $Ra$ number and buoyancy force increased. The isoconcentrations and isotherm density close to the heated wall increased when the buoyancy force shifted from a negative magnitude to a positive one. The local $Nu$ increased as the Rayleigh number increased but reduced as the volume fraction augmented. Furthermore, the mean $Nu$ ($\overline{Nu}$) decreased by $3.12\%$ and $6.81\%$ and the $\overline{Sh}$ increased by $83.17\%$ and $117.91\%$ with rising Lewis number for ($Ra={10}^5$ and $Da={10}^{-3}$) and ($Ra={10}^6$ and $Da={10}^{-4}$), respectively. Finally, the $Br$ and $Sr$ demonstrated positive sensitivity, and the $Ra$ and $\varphi$ showed negative sensitivity only for higher values of $\varphi$ based on the RSM.

1. Introduction

A colloidal suspension of nanoparticles is referred to as a nanofluid. Usually, the nanoparticles are metals, oxides, carbon nanotubes, base fluids, or carbides, typically, oil, ethylene glycol, or water. Nanofluids have become an essential component of industrial applications in recent years. The base fluid’s thermal conductivity is increased when nanoparticles are dissolved in it [1,2]. This particular characteristic of nanofluids has increased their use in various technical services, such as machining and manufacturing [3,4], refrigeration [5,6], and heat and mass transfer [7,8,9]. The type of numerical simulations and experimental analyses employed vary based on the applications and other relevant requirements.

Thermosolutal natural convection, which is also known as double-diffusive natural convection, is a phenomenon where buoyancy-driven fluid flow is induced by both the temperature gradients and the concentration gradients of elements concurrently [10,11]. The combined influence of the gradients produces two major effects within a system. A diffusion-thermo effect caused by the concentration gradient during the heat flux production is known as the Dufour effect. On the other hand, if the temperature gradient is responsible for producing the mass flux through a thermo-diffusion process, it is referred to as the Soret effect. The Dufour and Soret effects are not always considered in fluid dynamics due to their second-order genre and because they are of a comparatively smaller order of magnitude than the effects demonstrated by Fourier’s and Fick’s laws [12]. However, these two effects cannot be neglected in some exceptional cases, such as where there are steep temperature and concentration gradients. Thermosolutal natural convection is one of these exceptional cases in fluid dynamics where a coupled heat and mass transfer process is observed [10]. This governed coupling process often creates intriguing and complicated fluid flow patterns requiring meticulous analyses. Therefore, a coherent numerical technique is required to meet the abovementioned requirements. Ghachem et al. [13] studied the Soret and Dufour aspects of a viscoelastic fluid due to the movement of cylinders with viscous dissipation and convective boundary conditions, and they found increasing velocity profile results due to the curvature parameter and the viscoelastic parameter.

The lattice Boltzmann method (LBM) is a popular and powerful simulation technique to model complex flow patterns [2,14] numerically. The LBM considers the fluid density on a special lattice formation and develops the analysis through simulation with streaming and collision [15,16]. Among different relaxation variants, multiple-relaxation-time (MRT) is considered to be one of the most precise and positive options. The MRT-LBM option has been considered the most efficient for instability reduction. However, the selection of MRT-LBM combinations varies based on the geometric configuration, the types of transport phenomena, and the computational resources. While industries responsible for the thermal efficiency enhancement of thermo-mechanical devices will undoubtedly benefit from the inclusion of nanoparticles, certain issues require consideration. One is the need to impose control on the augmented heat transfer rate. Thermosolutal natural convection, heat, and mass transfer must be observed. Magnetohydrodynamics (MHD) is the study of fluid flow and electrically conducting fluid in a heat transfer process and is a strong candidate to serve as the controller [17]. A strong magnetic field works against the heat transfer process and could be used as a restriction if a heat transfer process needs to be stable to a certain extent. However, further analyses are required to investigate the effect of using LBM-MRT-MHD techniques on thermosolutal natural convection where both Dufour and Soret effects exist. Nevertheless, the proper geometry that can represent a certain heat and mass transfer unit needs to be selected.

Researchers have opted for various shaped cavities with certain boundary conditions to understand fluid flow phenomena. Some of them are C-shaped, U-shaped (by rotating the C-shaped geometry ${90}^{\circ }$ anticlockwise once), L-shaped, and T-shaped, to name a few. These cavities can help develop electronic heating and cooling systems where the optimal use of space is a necessity. A U-shaped cavity is an excellent option due to its symmetric nature and simpler design portfolio. There has been evidence of power enhancement by utilizing a U-shaped cavity microring resonator [18] and fracture evolution [19]. Considering the use of porous cavities can be beneficial to the flow control mechanism as the porosity can be varied for further rectification.

A number of studies have recently been published on U-shaped cavities. Ali et al. [20] investigated the natural convection of a non-Newtonian nanofluid in a U-shaped cavity in relation to MHD using the Galerkin finite element method. The power-law index was reported to impact the Rayleigh ($Ra$) number by ${10}^6$. The authors also emphasized the need to keep the magnetic parameter Hartmann ($Ha$) number below 30 to enhance the heat transfer mechanism. As described above, their well-written work focused on simple natural convection—the process would have been much more complicated under double-diffusive conditions. In addition, the analyses did not consider the correlations among the pertinent parameters. A similar report to that of Ali et al. [20] was also produced by Ma et al. [21], who focused more on the aspect ratio (AR) of the cavity in the presence of nanoparticles. However, the possible threshold of the concurrent effect of the buoyancy and heat transfer mechanisms needed further clarification. Unsteady natural convection was the area of interest in the study by Nabwey et al. [22], where the usefulness of the magnetic field in heat transfer control was demonstrated. Despite the numerical precision and accuracy of the different numerical analyses, no evidence was found for the thermosolutal natural convection of a nanofluid inside a U-shaped porous cavity. However, a number of well-cited published works that consider different cavities using various numerical methods are worth discussing. An overview is provided in tabular format in Table 1 for easier comparison.

According to our literature review, there have been several investigations on ferrofluids or nanofluids and natural convection in different enclosures using magnetic fields. Most of these studies were performed using non-Newtonian fluids and magnetic fields within a U-shaped cavity. To the best of the authors’ knowledge, there has never been an investigation that used an MRT-LB simulation of the Soret ($Sr$) and Dufour ($D_f$) effects on the thermosolutal natural convection of a Newtonian nanofluid in a U-shape porous enclosure, where $A{l}_2{O}_3$ or alumina is used as the nanoparticle and $H_{2}O$ is provided as the base fluid. Alumina has a relatively higher thermal conductivity and is not difficult to obtain. Therefore, the combination of alumina and water is relatively inexpensive and could be utilized inside any heat-exchange device. Furthermore, this work investigated the correlations among the relevant numerical parameters applying rigorous statistical analysis to provide a better insight into thermosolutal natural convection in the presence of the $Sr$ and $D_f$ effects. All the numerical simulations were carried out using GPU-accelerated computing.

2. Theoretical Formulation

2.1. Problem Statement

Figure 1a depicts a schematic and Figure 1b the mesh composition for the explanation domain chamber. This simulation investigates the natural heat transfer of an $A{l}_2{O}_3$-$H_{2}O$ nanofluid in a U-shaped enclosure with a steady laminar flow. The density is also computed using the Boussinesq approximation. The flow was steady-state and assumed to be composed of water and nanoparticles; the properties are listed in Table 2. The solution has a two-dimensional domain. A solid matrix having interconnected cavities is referred to as a porous medium. The matrix is considered inflexible or only slightly deformable under normal circumstances [28,29]. The fluid can flow due to the connections and interconnections of the cavities. At a stable temperature, $T_h$ heats the channel from below (the curved wall), while the U-shaped fence exposes it to low-temperature $T_c$. The temperatures of the upper and lower sides were likewise assumed to be adiabatic. Furthermore, the porous media were considered to be homogeneous. The grid size for the entire simulation is $256\times 256$. The grid density in the mesh composition figure would be too high with a $256\times 256$ grid size, which would not be clearly visible. Because of this, a grid of $128\times 128$ dimensions has been used to represent the mesh composition in Figure 1b.

2.2. Properties of Nanofluids

Some of the mathematical parameters need to be defined. For example, ${\rho }_{nf}$ is the effective density, the heat capacity is ${\left(\rho C_p\right)}_{nf}$, and the expansion of the thermal coefficient of the nanofluid is ${\left(\rho \beta \right)}_{nf}$, which is provided in the following ways, as described in [17,33]:

(1)${\rho }_{nf}=\varphi {\rho }_s+{\rho }_f(1-\varphi )$

(2)${\left(\rho c_p\right)}_{nf}={\left(\rho c_p\right)}_s\varphi +{\left(\rho c_p\right)}_f(1-\varphi )$

and

(3)${\left(\rho \beta \right)}_{nf}=\varphi {\left(\rho \beta \right)}_s+(1-\varphi ){\left(\rho \beta \right)}_f$

where $\varphi$ is the volume fraction of the nanoparticles, the subscripts f and s are the base fluid and the solid particle, respectively. The viscosity of the nanofluid $\left({\overline{\mu }}_{nf}\right)$ using the power-law model for a Newtonian fluid is defined by following the Corcione et al. model [34,35]:

(4)$\begin{array}{c}\hfill {\overline{\mu }}_{nf}=\frac{{\mu }_f}{\left(1-34.87{\left(\frac{d_p}{d_f}\right)}^{-0.3}{\varphi }^{1.03}\right)}\end{array}$

The consistency factor for the power-law model is ${\mu }_f$. Here, the Newtonian pure fluid is represented by $\varphi =0$, and the Newtonian fluid viscosity is ${\overline{\mu }}_{nf}={\mu }_f$. The diameter of the base fluid molecule is given by $d_f$, the nanoparticle size, $d_p=25$ nm, is uniform, and the shape is spherical.

(5)$\begin{array}{c}\hfill d_f=0.1{\left(\frac{6M}{N_A\pi {\rho }_f}\right)}^{\frac{1}{3}}\end{array}$

where the molecular weight of the base fluid is $M=18.01528\times {10}^{-3}$ kg mol${}^{-1}$ and the Avogadro number is $N_A$ = $6.022\times {10}^{23}$ mol${}^{-1}$. At temperature $T_0$ = 293 K, the mass density of the base fluid ${\rho }_f$ is calculated. The nanofluid’s temperature-dependent effective thermal conductivity, $k_{nf}$, is described in Equation (6), according to [34,35]:

(6)$\begin{array}{c}\hfill \frac{k_{nf}}{k_f}=1+4.4R{e}_p^{0.4}{\left(\frac{k_s}{k_f}\right)}^{0.03}{\varphi }^{0.66}{\left(\frac{T_f}{T_{fr}}\right)}^{10}P{r}^{0.66}\end{array}$

(7)$\begin{array}{c}\hfill R{e}_p=\frac{2{\rho }_f{k}_b{T}_f}{\pi {\mu }_f^2{d}_P}\end{array}$

where the Reynolds number of the nanofluid is $R{e}_p$ and the Prandtl number for the base fluid is $Pr=6.2$. The freezing temperature of the base fluid is $T_{fr}=273$ K, and $k_f$ is the effective thermal conductivity of the fluid. $k_s$ is the effective thermal conductivity of the nanoparticles. For the present simulation, $T_f$ is calculated according to the relation $T_f=T_0+T(t,\overline{x},\overline{y})\Delta T$ and $\Delta T=50\mathrm{K}$ is considered for this study.

2.3. Macroscale Dimensional Governing Equations

Continuity equation:

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

u-momentum equation:

(9)$\begin{array}{cc}\hfill {\rho }_{nf}\left(\frac{\partial \overline{u}}{\partial \overline{t}}+\frac{\overline{u}}{ϵ}\frac{\partial \overline{u}}{\partial \overline{x}}+\frac{\overline{v}}{ϵ}\frac{\partial \overline{u}}{\partial \overline{y}}\right)=& -ϵ\frac{\partial \overline{P}}{\partial \overline{x}}+{\overline{\mu }}_{nf}\frac{{\partial }^2\overline{u}}{\partial {\overline{x}}^2}+{\overline{\mu }}_{nf}\frac{{\partial }^2\overline{u}}{\partial {\overline{y}}^2}+F_x\hfill \end{array}$

v-momentum equation:

(10)$\begin{array}{cc}\hfill {\rho }_{nf}\left(\frac{\partial \overline{v}}{\partial \overline{t}}+\frac{\overline{u}}{ϵ}\frac{\partial \overline{v}}{\partial \overline{x}}+\frac{\overline{v}}{ϵ}\frac{\partial \overline{v}}{\partial \overline{y}}\right)=& -ϵ\frac{\partial \overline{P}}{\partial \overline{y}}+{\overline{\mu }}_{nf}\frac{{\partial }^2\overline{v}}{\partial {\overline{x}}^2}+{\overline{\mu }}_{nf}\frac{{\partial }^2\overline{v}}{\partial {\overline{y}}^2}+gϵ{\left(\rho {\beta }_T\right)}_{nf}\left(T-T_c\right)+\hfill \\ \hfill & gϵ{\left(\rho {\beta }_C\right)}_{nf}\left(C-C_c\right)+F_y\hfill \end{array}$

Energy equation:

(11)$\begin{array}{c}\hfill \frac{\partial T}{\partial \overline{t}}+\overline{u}\frac{\partial T}{\partial \overline{x}}+\overline{v}\frac{\partial T}{\partial \overline{y}}=\frac{1}{{\left(\rho C_p\right)}_{nf}}\left[\frac{\partial }{\partial \overline{x}}\left(k_{nf}\frac{\partial T}{\partial \overline{x}}\right)+\frac{\partial }{\partial \overline{y}}\left(k_{nf}\frac{\partial T}{\partial \overline{y}}\right)\right]+D_{TC}\left(\frac{{\partial }^{2}C}{\partial {\overline{x}}^2}+\frac{{\partial }^{2}C}{\partial {\overline{y}}^2}\right)\end{array}$

Concentration equation:

(12)$\begin{array}{c}\hfill \frac{\partial C}{\partial \overline{t}}+\overline{u}\frac{\partial C}{\partial \overline{x}}+\overline{v}\frac{\partial C}{\partial \overline{y}}=D\left(\frac{{\partial }^{2}C}{\partial {\overline{x}}^2}+\frac{{\partial }^{2}C}{\partial {\overline{y}}^2}\right)+D_{CT}\left(\frac{{\partial }^{2}T}{\partial {\overline{x}}^2}+\frac{{\partial }^{2}T}{\partial {\overline{y}}^2}\right)\end{array}$

Here, T denotes the nanofluid’s temperature, and $\overline{u}$ and $\overline{v}$ are the velocity components’ horizontal and vertical dimensions, respectively. The gravitational acceleration is g and C is the concentration of the fluid. $F_x$ and $F_y$ of the Brinkman–Forchheimer extended Darcy model (BFEDM) for the Newtonian nanofluids are, respectively, as follows [36]:

(13)$F_x=-\frac{ϵ{\mu }_{nf}}{K}\overline{u}-\frac{1.75{\rho }_{nf}}{\sqrt{150ϵK}}\left|\mathbf{u}\right|\overline{u}$

(14)$F_y=-\frac{ϵ{\mu }_{nf}}{K}\overline{v}-\frac{1.75{\rho }_{nf}}{\sqrt{150ϵK}}\left|\mathbf{u}\right|\overline{v}$

where ${\beta }_0$ is the applied uniform magnetic field along the $\overline{x}$ direction. $\rho$ is the density, $C_p$ is the heat capacity, and the thermal expansion coefficient of the fluid is $\beta$. The appropriate boundary conditions to solve Equations (8)–(12) are:

Left wall:

(15)$\begin{array}{c}\hfill \overline{u}=\overline{v}=0,T=T_h,C=C_h,x=0,0\le y\le H;\end{array}$

Right wall:

(16)$\begin{array}{c}\hfill \overline{u}=\overline{v}=0,T=T_h,C=C_h,x=L,0\le y\le H;\end{array}$

Lower wall:

(17)$\begin{array}{c}\hfill \overline{u}=\overline{v}=0,y=0,\frac{\partial T}{\partial y}=\frac{\partial C}{\partial y}=0,x=L,0\le x\le L;\end{array}$

Upper wall:

(18)$\begin{array}{cc}\hfill & Aty=H,0\le x\le \frac{3H}{10},\frac{\partial T}{\partial y}=\frac{\partial C}{\partial y}=0,u=v=0;\hfill \\ \hfill & Aty=H,\frac{7L}{10}\le x\le L,\frac{\partial T}{\partial y}=\frac{\partial C}{\partial y}=0,u=v=0;\hfill \\ \hfill & At\frac{3H}{10},\frac{3L}{10}\le x\le \frac{7L}{10},T=T_c,u=v=0;\hfill \\ \hfill & Atx=\frac{3L}{10},\frac{3H}{5}\le y\le H,T=T_c,u=v=0;\hfill \\ \hfill & Atx=\frac{7H}{10},\frac{3H}{5}\le y\le H,T=T_c,u=v=0;\hfill \end{array}$

Here, the height is H and the length is L. H is equal to L, $(H=L)$.

2.4. Non-Dimensional Equation

To non-dimensionalize the governing Equations (8)–(12), the following dimensionless variables [37,38] are used:

(19)$x=\frac{\overline{x}}{L},y=\frac{\overline{y}}{L},u=\frac{\overline{u}L}{{\alpha }_f\sqrt{Ra}},v=\frac{\overline{v}L}{{\alpha }_f\sqrt{Ra}},t=\frac{\alpha \overline{t}\sqrt{Ra}}{L^2}$

(20)$p=\frac{\overline{p}{L}^2}{{\rho }_{nf}{\alpha }_f^{2}Ra},\theta =\frac{T-T_c}{T_h-T_c},\mathsf{\Phi }=\frac{C-C_c}{C_h-C_c},{\nu }_f=\frac{{\mu }_f}{{\rho }_f},{\alpha }_f=\frac{k_f}{{\left(\rho C_p\right)}_f}$

(21)$D_f=\frac{D_{TC}\Delta C}{{\alpha }_f\Delta T},Sr=\frac{D_{CT}\Delta T}{D\Delta C},Pr=\frac{{\nu }_f}{{\alpha }_f},Le=\frac{{\alpha }_f}{D}$

(22)$Ra=\frac{g{\beta }_{T_f}\Delta T{L}^3}{{\nu }_f{\alpha }_f},Br=\frac{{\beta }_{C_{nf}}\Delta C}{{\beta }_{T_{nf}}\Delta T},Da=\frac{K}{H^2}$

where u and v denote the nanofluid’s non-dimensional velocity along the x and y’s non-dimensional directions. The dimensionless forms of the temperature, time, solutal concentration, and pressure are denoted by the symbols $\theta$, t, $\mathsf{\Phi }$, and p, respectively.${\nu }_f$ is the kinematic viscosity, D is the mass diffusivity, and ${\alpha }_f$ is the thermal diffusivity of the base fluid. $Br$ is the buoyancy ratio and the Rayleigh number is $Ra$. The Prandtl number is $Pr$, $Le$ is the Lewis number, and K is the permeability of the porous medium.

(23)$\frac{\partial u}{\partial x}+\frac{\partial v}{\partial y}=0$

(24)$\begin{array}{c}\hfill \frac{\partial u}{\partial t}+\frac{u}{ϵ}\frac{\partial u}{\partial x}+\frac{v}{ϵ}\frac{\partial u}{\partial y}=-ϵ\frac{\partial p}{\partial x}+\left(\frac{{\mu }_{nf}}{{\mu }_f}\right)\left(\frac{{\rho }_f}{{\rho }_{nf}}\right)\frac{Pr}{\sqrt{Ra}}\left[\frac{{\partial }^{2}u}{\partial x^2}+\frac{{\partial }^{2}u}{\partial y^2}-ϵ\frac{u}{Da}\right]-\frac{1.75}{\sqrt{150ϵDa}}\left|u\right|u\end{array}$

(25)$\begin{array}{cc}\hfill \frac{\partial v}{\partial t}+\frac{u}{ϵ}\frac{\partial v}{\partial y}+\frac{v}{ϵ}\frac{\partial v}{\partial y}=& -ϵ\frac{\partial p}{\partial y}+\left(\frac{{\mu }_{nf}}{{\mu }_f}\right)\left(\frac{{\rho }_f}{{\rho }_{nf}}\right)\frac{Pr}{\sqrt{Ra}}\left[\frac{{\partial }^{2}v}{\partial x^2}+\frac{{\partial }^{2}v}{\partial y^2}-ϵ\frac{v}{Da}\right]+ϵ\left(\frac{{\beta }_{T_{nf}}}{{\beta }_{T_f}}\right).\hfill \\ \hfill & Pr\left[\theta +Br\mathsf{\Phi }\right]-\frac{1.75}{\sqrt{150ϵDa}}\left|u\right|u\hfill \end{array}$

(26)$\begin{array}{c}\hfill \frac{\partial \theta }{\partial t}+u\frac{\partial \theta }{\partial x}+v\frac{\partial \theta }{\partial y}=\frac{{\left(\rho C_p\right)}_f}{{\left(\rho C_p\right)}_{nf}}\frac{1}{\sqrt{Ra}}\left[\frac{\partial }{\partial x}\left(\frac{k_{nf}}{k_f}\frac{\partial \theta }{\partial x}\right)+\frac{\partial }{\partial y}\left(\frac{k_{nf}}{k_f}\frac{\partial \theta }{\partial \overline{y}}\right)\right]+\frac{D_f}{\sqrt{Ra}}\left(\frac{{\partial }^2\mathsf{\Phi }}{\partial x^2}+\frac{{\partial }^2\mathsf{\Phi }}{\partial y^2}\right)\end{array}$

(27)$\begin{array}{c}\hfill \frac{\partial \mathsf{\Phi }}{\partial t}+u\frac{\partial \mathsf{\Phi }}{\partial x}+v\frac{\partial \mathsf{\Phi }}{\partial y}=\frac{1}{Le\sqrt{Ra}}\left[\left(\frac{{\partial }^2\mathsf{\Phi }}{\partial x^2}+\frac{{\partial }^2\mathsf{\Phi }}{\partial y^2}\right)+Sr\left(\frac{{\partial }^2\theta }{\partial x^2}+\frac{{\partial }^2\theta }{\partial y^2}\right)\right]\end{array}$

${\left(\rho C_p\right)}_{nf}$ represents the effective heat capacity and ${\rho }_{nf}$ the effective density. The thermal expansion coefficient of the nanofluid is ${\left(\rho {\beta }_T\right)}_{nf}$.

2.5. Numerical Methods

The physical quantities of the fluid, such as the velocity, temperature, and solute concentration, were computed for the current study using the D2Q9 and D2Q5 MRT-LBM lattice models [39,40,41]. The equations’ formulations are briefly described below:

(28)$\mathbf{f}(\mathbf{r}+{\mathbf{e}}_i\Delta t,t+\Delta t)-\mathbf{f}(\mathbf{r},t)=-{\mathbf{M}}^{-\mathbf{1}}\mathbf{S}\left[\mathbf{m}(\mathbf{r},t)-{\mathbf{m}}^{eq}(\mathbf{r},t)\right]+{\mathbf{M}}^{-\mathbf{1}}\left(\mathbf{I}-\frac{\mathbf{S}}{2}\right)\mathbf{F}(\mathbf{r},t)$

(29)$\mathbf{g}(\mathbf{r}+{\mathbf{e}}_i\Delta t,t+\Delta t)-\mathbf{g}(\mathbf{r},t)=-{\mathbf{N}}^{-\mathbf{1}}{\mathbf{S}}_{\mathbf{T}}\left[{\mathbf{m}}_{\mathbf{T}}(\mathbf{r},t)-{\mathbf{m}}_T^{eq}(\mathbf{r},t)\right]+{\mathbf{N}}^{-\mathbf{1}}\left(\mathbf{I}-\frac{{\mathbf{S}}_{\mathbf{T}}}{2}\right){\mathbf{Q}}_{\mathbf{T}}(\mathbf{r},t)$

(30)$\mathbf{h}(\mathbf{r}+{\mathbf{e}}_i\Delta t,t+\Delta t)-\mathbf{h}(\mathbf{r},t)=-{\mathbf{N}}^{-\mathbf{1}}{\mathbf{S}}_{\mathbf{C}}\left[{\mathbf{m}}_{\mathbf{C}}(\mathbf{r},t)-{\mathbf{m}}_C^{eq}(\mathbf{r},t)\right]+{\mathbf{N}}^{-\mathbf{1}}\left(\mathbf{I}-\frac{{\mathbf{S}}_{\mathbf{C}}}{2}\right){\mathbf{Q}}_{\mathbf{C}}(\mathbf{r},t)$

where $\mathbf{r}=(\overline{x},\overline{y})$, $\mathbf{f}={(f_0,f_1,f_2,\dots ,f_8)}^T$. Here, the vector moments are $\mathbf{m}=(m_0,m_1,$$m_2,\dots ,m_8{)}^T$ and ${\mathbf{m}}^{eq}={({m_0}^{eq},{m_1}^{eq},{m_2}^{eq},\dots ,{m_8}^{eq})}^T$ and the forcing components are $\mathbf{F}={(F_0,F_1,F_2,\dots ,F_8)}^T$. Similarly, for the D2Q5 model, $\mathbf{g}={(g_0,g_1,g_2,\dots ,g_4)}^T$ and $\mathbf{h}={(h_0,h_1,h_2,\dots ,h_4)}^T$ are the density distribution processes for the temperature and the solutal concentration, receptively. Here, ${\mathbf{Q}}_{\mathbf{T}}$ and ${\mathbf{Q}}_{\mathbf{C}}$ are source terms for the Soret and Dufour effects in the energy and solutal concentrations equations, respectively.

The velocity and the moment spaces’ mapping for the momentum equation are modified by a linear transformation, considering $\mathbf{m}=\mathbf{M}\mathbf{f}$, which implies $\mathbf{f}={\mathbf{M}}^{-1}\mathbf{m}$, where

(31)$\mathbf{M}=\left[\begin{array}{ccccccccc}1& 1& 1& 1& 1& 1& 1& 1& 1\\ -4& -1& -1& -1& -1& 2& 2& 2& 2\\ 4& -2& -2& -2& -2& 1& 1& 1& 1\\ 0& 1& 0& -1& 0& 1& -1& -1& 1\\ 0& -2& 0& 2& 0& 1& -1& -1& 1\\ 0& 0& 1& 0& -1& 1& 1& -1& -1\\ 0& 0& -2& 0& 2& 1& 1& -1& -1\\ 0& 1& -1& 1& -1& 0& 0& 0& 0\\ 0& 0& 0& 0& 0& 1& -1& 1& -1\end{array}\right]$

In Equation (28), There are nine eigenvalues in the diagonal collision matrix S, with values ranging from 0 to 2 and given in Equation (32):

(32)$\mathbf{S}=\mathbf{diag}[s_0,s_1,s_2,s_3,s_4,s_5,s_6,s_7,s_8]$

where $s_0=s_3=s_5=1.0,s_1=s_2=1.4,s_4=s_6=1.2$ and $s_7=s_8=\frac{1}{\tau }$. Here, $\tau =3\nu (\mathbf{r},t)+\frac{1}{2}$ is the relaxation time, where $\nu (\mathbf{r},t)$ is the kinematic viscosity computed by the Newtonian power-law viscosity model:

(33)$\nu (\mathbf{r},t)={\nu }_{nf}=\frac{{\rho }_{f}Pr}{{\rho }_{nf}\left(1-34.87{\left(\frac{d_p}{d_f}\right)}^{-0.3}{\varphi }^{1.03}\right)\sqrt{Ra}}$

The discrete velocity distribution function ${\mathbf{e}}_i$ for the D2Q9 model is defined as in [41] along the velocity direction.

In Equation (28), ${\mathbf{m}}^{\mathbf{eq}}$ is defined by

(34)${\mathbf{m}}^{eq}=\left[\rho ,2\rho +3\left|j\right|,\rho -3\left|j\right|,j_x,-j_x,j_y,-j_y,(j_x^2-j_y^2),j_x{j}_y\right]$

where $(j_x,j_y)=(\rho \overline{u},\rho \overline{v})$ and $\left|j\right|=j_x^2+j_y^2$. Here, $\rho$ and $\mathbf{u}$ are given by [42],

(35)$\rho =\sum _{i=0}^8{f}_i$

(36)$\mathbf{u}=\frac{1}{\rho }\sum _{i=0}^8{f}_i{\mathbf{e}}_i$

The collision matrix N for D2Q5 is shown below,

(37)$\mathbf{N}=\left[\begin{array}{ccccc}1& 1& 1& 1& 1\\ 0& 1& 0& -1& 0\\ 0& 0& 1& 0& -1\\ -4& 1& 1& 1& 1\\ 0& 1& -1& 1& -1\end{array}\right]$

The five eigenvalues of $S_T$ and $S_C$ for the D2Q5 model are each represented as a flow as follows:

(38)$S_T=diag[s_0,s_{\alpha },s_{\alpha },s_e,s_{\nu }]$

and

(39)$S_C=diag[s_0,s_D,s_D,s_e,s_{\nu }]$

As indicated in Equations (40) and (43), the formulation and explanation of the $s_i$ values were exhaustively provided in [40,42]:

(40)$s_0=1,\frac{1}{s_e}-\frac{1}{2}=\frac{1}{s_{\nu }}-\frac{1}{2}=\frac{1}{6}$

Here, the formula presented in Mezrhad et al. [43] is used to calculate $s_{\alpha }$:

(41)$s_{\alpha }=\frac{1}{\frac{1}{2}+5\alpha }$

The thermal diffusivity is computed as follows for the issue in this study:

(42)$\alpha =\frac{{\alpha }_{nf}}{{\alpha }_f}\frac{1}{\sqrt{Ra}}$

Similarly, $s_D$ is calculated as:

(43)$s_D=\frac{1}{\frac{1}{2}+5D}$

The mass diffusivity, D, is calculated for the current Newtonian nanofluid flow simulation as follows:

(44)$D=\frac{1}{Le\sqrt{Ra}}$

The equilibrium moments ${\mathbf{m}}_T^{eq}$ and ${\mathbf{m}}_C^{eq}$ for the distribution functions $g_i$ and $h_i$ are, respectively:

(45)${\mathbf{m}}_T^{eq}=\left[T,\overline{u}T,\overline{v}T,aT,0\right]$

(46)${\mathbf{m}}_C^{eq}=\left[C,\overline{u}C,\overline{v}C,aC,0\right]$

The discrete velocities for the D2Q5 model are defined in [41], where $a=-2$ is less than one to prevent numerical instability [43].

The following definitions are given for the fluid temperature, T, and C is the solute concentration:

(47)$T=\sum _{i=0}^4{g}_i$

and

(48)$C=\sum _{i=0}^4{h}_i$

2.6. Wall Boundary Conditions

Application of the bounce-back condition to the wall is performed as follows:

(49)$f_{\overline{i}}(r_f,t+\Delta t)=f_i^c(r_f,t)$

The following boundary conditions were used at a fixed temperature and concentration ($T_w$ and $C_w$) wall for the D2Q5 model MRT thermal and solutal LB simulations:

(50)$g_{\overline{i}}({\mathbf{r}}_f,t+\Delta t)=2T_w\frac{(1+\frac{a}{4})}{5}-g_i^c({\mathbf{r}}_f,t)$

and

(51)$h_{\overline{i}}({\mathbf{r}}_f,t+\Delta t)=2C_w\frac{(1+\frac{a}{4})}{5}-h_i^c({\mathbf{r}}_f,t)$

The adiabatic walls were operated under the following zero-gradient, temperature and concentration conditions:

(52)$g_i(r_f,t+\Delta t)=g_i^c(r_f,t)$

and

(53)$h_i(r_f,t+\Delta t)=h_i^c(r_f,t)$

2.7. Heat and Mass Transfer Rate

The mean heat transfer rate is calculated by,

(54)$\overline{Nu}=\frac{1}{H}{\int }_0^{H}Nu\left(y\right)dy{{|}_{leftwall}+\frac{1}{H}{\int }_0^{H}Nu\left(y\right)dy|}_{rightwall}$

and the mean mass transfer rate is then computed by

(55)$\overline{Sh}=\frac{1}{H}{\int }_0^{H}Sh\left(y\right)dy{{|}_{leftwall}+\frac{1}{H}{\int }_0^{H}Sh\left(y\right)dy|}_{rightwall}$

where the local Nusselt number is defined by

(56)$Nu\left(y\right)=-\frac{k_{nf}}{k_f}\frac{L}{\Delta T}\frac{\partial T}{\partial \overline{x}}{|}_{wall}-D_f\frac{L}{\Delta C}\frac{\partial C}{\partial \overline{x}}{|}_{wall}$

and the local Sherwood number is finally determined by

(57)$Sh\left(y\right)=-\frac{L}{\Delta C}\frac{\partial C}{\partial \overline{x}}{|}_{wall}-S_r\frac{k_{nf}}{k_f}\frac{L}{\Delta T}\frac{\partial T}{\partial \overline{x}}{|}_{wall}$

2.8. Formulation of Entropy Production

The current study examines the irreversibilities or entropy creation resulting from magnetic fields, heat transfer, and fluid friction. The overall entropy is the sum of the irreversibilities caused by the thermal gradients, viscous dissipation, and the concentration gradients [44,45]:

(58)${\overline{M}}_S={\overline{M}}_F+{\overline{M}}_T+{\overline{M}}_D$

where the irreversibilities due to fluid friction (${\overline{M}}_F$), heat transfer (${\overline{M}}_T$) and mass transfer (${\overline{M}}_D$) are given by:

(59)${\overline{M}}_F=\frac{{\mu }_f}{T_0(1-34.87{\left(\frac{d_p}{d_f}\right)}^{-0.3}{\varphi }^{1.03})}\left[2{\left(\frac{\partial \overline{u}}{\partial \overline{x}}\right)}^2+2{\left(\frac{\partial \overline{v}}{\partial \overline{y}}\right)}^2+{\left(\frac{\partial \overline{v}}{\partial \overline{x}}+\frac{\partial \overline{u}}{\partial \overline{y}}\right)}^2\right]$

(60)${\overline{M}}_T=\frac{k_{nf}}{T_0^2}\left[{\left(\frac{\partial T}{\partial \overline{x}}\right)}^2+{\left(\frac{\partial T}{\partial \overline{y}}\right)}^2\right]$

(61)${\overline{M}}_D=\frac{RD}{C_0}\left[{\left(\frac{\partial C}{\partial \overline{x}}\right)}^2+{\left(\frac{\partial C}{\partial \overline{y}}\right)}^2\right]+\frac{RD}{T_0}\left[\left(\frac{\partial C}{\partial \overline{x}}\right)\left(\frac{\partial T}{\partial \overline{x}}\right)+\left(\frac{\partial C}{\partial \overline{y}}\right)\left(\frac{\partial T}{\partial \overline{y}}\right)\right]$

Bejan numbers measure the proportion of entropy generations resulting from the irreversible mass and heat transfers to the overall entropy generation and can be expressed as follows:

(62)$\overline{Be}=\frac{{\overline{M}}_T+{\overline{M}}_D}{{\overline{M}}_S}$

The non-dimensionalized equations are denoted as follows:

(63)$M_S=M_F+M_T+M_D$

(64)$M_F={\overline{M}}_F\times \frac{T_0^2{L}^2}{k_f\Delta T^2}$

(65)$M_F=\frac{{\lambda }_1}{(1-34.87{\left(\frac{d_p}{d_f}\right)}^{-0.3}{\varphi }^{1.03})}\left[2{\left(\frac{\partial u}{\partial x}\right)}^2+2{\left(\frac{\partial v}{\partial y}\right)}^2+{\left(\frac{\partial v}{\partial x}+\frac{\partial u}{\partial y}\right)}^2\right]$

Here,

(66)${\lambda }_1=\frac{T_0{\nu }_{f}Ra{\alpha }_f^2}{L^2{k}_f\Delta T^2}$

(67)$M_T={\overline{M}}_T\times \frac{T_0^2{L}^2}{k_f\Delta T^2}$

(68)$M_T=\frac{k_{nf}}{k_f}\left[{\left(\frac{\partial \theta }{\partial x}\right)}^2+{\left(\frac{\partial \theta }{\partial y}\right)}^2\right]$

(69)$M_D={\overline{M}}_D\times \frac{T_0^2{L}^2}{k_f\Delta T^2}$

(70)$M_D={\lambda }_2\left[{\left(\frac{\partial \mathsf{\Phi }}{\partial x}\right)}^2+{\left(\frac{\partial \mathsf{\Phi }}{\partial y}\right)}^2\right]+{\lambda }_3\left[\left(\frac{\partial \mathsf{\Phi }}{\partial x}\right)\left(\frac{\partial \theta }{\partial x}\right)+\left(\frac{\partial \mathsf{\Phi }}{\partial y}\right)\left(\frac{\partial \theta }{\partial y}\right)\right]$

where,

(71)${\lambda }_2=\frac{RD}{C_0}\frac{T_0^2\Delta C^2}{k_f\Delta T^2}$

(72)${\lambda }_3=\frac{RD}{k_f}\frac{\Delta C}{\Delta T}{T}_0$

The definition of the dimensionless local Bejan number is,

(73)$Be=\frac{M_T+M_D}{M_S}$

The following are the formulas that provide the total quantities of the irreversibilities across the entire domain:

(74)$M{F}_t={\int }_0^1{\int }_0^1{M}_{F}dxdy$

(75)$M{T}_t={\int }_0^1{\int }_0^1{M}_{T}dxdy$

(76)$M{S}_t={\int }_0^1{\int }_0^1{M}_{S}dxdy$

(77)$M{D}_t={\int }_0^1{\int }_0^1{M}_{D}dxdy$

(78)$B{e}_{avg}=\frac{{\int }_0^1{\int }_0^{1}Bedxdy}{{\int }_0^1{\int }_0^{1}dxdy}$

For the present problem, the parameters ${\lambda }_1={10}^{-4}$, as in Illis et al. [44], and ${\lambda }_2={\lambda }_3=0.5$.

3. Grid Independent Test

For acceptable accuracy of the numerical results, considering the proper grid size is essential to any investigation. For the present problem, three diffident grid sizes were chosen and the corresponding results in terms of the mean Nusselt and Sherwood numbers are given in the Table 3 at $Ra={10}^6$, $Da={10}^{-4}$, $Le=5$, $\varphi =0.04$, $Br=1$, $D_f=0.2$ and $ϵ=0.4$. After considering the results, a $256\times 256$ grid size was selected for the whole simulation.

4. Code Validation

In this section, another essential test was performed for the numerical validity. Two validation tests were conducted in the present study, as discussed in the following subsections.

4.1. Code Validation for the Porous Media

The code was also tested for the porous media with different $Ra$, as reported in Table 4. The present MRT-LBM results were compared with the results of Nithiarasu et al. [46] and Guo and Zhao [47]. This comparison ensures the present code can simulate the fluid flow and heat transfer in porous media.

4.2. Code Validation for the Soret and Dufour Effects

Table 5 shows another validation for the cross-diffusion with the results of Xu et al. [11] and Ren et al. [37] for the different Rayleigh numbers while $Pr=1.0$, $Sr=0.1$, $Le=2.0$, $Dr=0.1$, $\varphi =0.0$, $ϵ=1$ and $Da={10}^9$. Here, the computational domain was a side-heated vertical cavity with an aspect ratio of $A=H/L=2.0$, where the height is H and L is the width of the cavity. This table also ensures that the present LBM code was validated for the cross-diffusion with the Soret and Dufour effects.

A qualitative comparison is shown in Figure 2 for the Soret and Dufour effects in terms of the isotherms and the isoconcentration lines obtained by Ren and Chan [37] using the single-relaxation-time (SRT) LBM while $Ra={10}^5$, $Pr=1.0$, $Le=2.0$, $Sr=0.1$, $D_f=0.1$, $Br=-10$, $\varphi =0.0$, $ϵ=1$ and $Da={10}^9$. By inspection, this qualitative comparison also confirmed the code validation.

5. Results and Discussion

The MRT-LB Simulation was applied to simulate the cross-diffusion natural convection of a nanofluid in a porous U-shape enclosure. The parameters specified here are the Rayleigh number, ${10}^4\le Ra\le {10}^6$; the volume fraction of the nanoparticle, $0.0\le \varphi \le 0.04$; the buoyancy ratio, $-2\le Br\le 2$; the Darcy number, ${10}^{-4}\le Da\le {10}^{-2}$; the Dufour number, $0\le D_f\le 0.2$; the Soret number, $0\le Sr\le 0.2$; and the Prandtl number, $Pr$ = $6.2$. In the velocity and temperature distribution, we vary the porosity, $ϵ$ = $0.2$ to $0.8$. The variant of the Lewis number, $Le$ = 1 to 5, is also shown in $\overline{Nu}$ and $\overline{Sh}$ as a bar chart.

5.1. Observation of Streamlines, Isotherms and Isoconcentrations

Figure 3 represents the isotherms, streamlines, and isoconcentrations for various $Ra$ numbers and $Da$ numbers at $Sr=0$, $Br=1$, $D_f=0.2$, $ϵ$ = $0.4$, and $\varphi$ = $0.0$ (solid line), $0.04$ (long-dashed line). The fluid responded by passing through the heated wall first, then the cold wall, moving upward and lower. For $Ra$ = ${10}^4$, ${10}^5$, and ${10}^6$, two vortexes were created on the left and heated right sides of the enclosure. The velocity changes across the vertical walls were greater than in the core zone of the left vortex, resulting in dense streamlines. Nevertheless, the concentrations were comparatively lower along the steep walls than in the middle area of the right vortex. The configurations of the streamlines and the isotherms were organized and distributed equally in all circumstances. The dimensions of the vortex grew rapidly as the $Ra$ numbers increased. The buoyancy effect caused the flow to circulate at the center of the enclosure and to descend along each side wall, creating mirror-image structures that revolved in the right and left half counter-clockwise. Both cells presented elliptical shapes. The symmetry ensured that the flows in the enclosure’s left and right halves were identical. As the fluid permeability across a cross-sectional area decreased, reduction in the $Da$ number resulted in less buoyancy-induced flow, progressively thickening the boundary layer. The heat transfer decreased slightly as the $Ra$ number increased and the Darcy number decreased. The transfer of mass changed as shown in Figure 3g–i. A higher $Ra$ number and a lower $Da$ number resulted in a synergistic rise in the mass exchange rate. The lower Darcy number complemented the higher flow velocity and convective transport generated by the higher $Ra$ number to produce a much more effective mass transfer process in the fluid.

In addition, the buoyancy ratio was calculated using the density difference between the two fluids relative to their mean densities. This is essential for predicting and analyzing fluid behavior when the density differences are significant. The streamlines, isotherms, and isoconcentrations for the different buoyancy ratios with some fixed parameters are shown in Figure 4, Figure 5 and Figure 6. Two vortexes formed at the cavity’s upward portion in the streamline while the buoyancy force was opposing. The $Ra$ number is comparatively lower in Figure 4a–c. As the fluid flow increased, the flow structure became irregular due to the $Ra$ number’s augmentation and the reduction in the $Da$ number, as shown in Figure 4d–i. The flow structure became erratic when the buoyancy force was applied from negative to positive due to the gravitational forces interacting with the density differences, leading to more pronounced and dominant buoyancy-driven convection. Positive buoyancy encourages vigorous fluid motion, making it a critical factor in natural convection, stratified flows, and other buoyancy-dominated fluid phenomena. Unsurprisingly, the most dominant fluid flow was observed when the $Ra$ number was at a maximum, and the buoyancy ratio was favorable, as shown in Figure 4i. As Figure 5 and Figure 6 show, the isotherms and isoconcentrations were consistent at opposing the buoyancy forces. The isoconcentrations and isotherm density close to the heated wall grew when the buoyancy force shifted from a negative to a positive value, indicating a rise in mass and heat transfer from the hot wall. Additionally, as the $Ra$ number and buoyancy ratio increased, the isotherms and isoconcentrations became more inconsistent.

5.2. Velocity and Temperature Effects

The changes in the velocity and temperature profiles under variable input parameters will be discussed here. In all cases, the u-velocity and the v-velocity will be represented along the x-axis and the y-axis, respectively, as shown in Figure 7a,b. The temperature profile along the y-axis for the different buoyancy effects (Br) is shown in Figure 7c. The other input parameters were as follows: $Ra$ = ${10}^4$, $Da={10}^{-2}$, $D_f=0.2$, $Sr=0.2$, $\varphi =0.04$ and $ϵ=0.4$. Both the $Sr$ and $D_f$ impacts are present in this case. When the buoyancy effect was increased, the u-velocity grew in the region $0<y<0.3$ but plummeted in the area $0.3<y<0.6$, as shown in Figure 7a, forming a symmetric shape. Unsurprisingly, the symmetrical shape was only visible at $Br=2$ as the higher buoyancy influenced the proper velocity development. At $y>0.6$, the flow met the adiabatic wall of the cavity; hence, the fluid was immobile regardless of the $Br$ values. Meanwhile, as per Figure 7b, the v-velocity went up rapidly at $Br=2$ in the regions $0<x<0.2$ and $0.8<x<1$. A similar pattern was observed at $Br=0$ as well, but not as sharp as that of $Br=2$. However, as $Br=-2$ implied, the opposite pattern was observed regarding the v-velocity due to the negative magnitude of the buoyancy effect. In all cases, the symmetrical patterns were observed due to the symmetric configuration of a U-shaped cavity. Finally, the temperature distribution exhibited consistent patterns across all $Br$ values, as shown in Figure 7c. The temperature distribution declined rapidly as the enhancing buoyancy ratio increased due to the buoyancy’s opposing force to the gravitational force. The strength of the buoyancy force grew relative to the gravitational force as the buoyancy ratio was augmented, resulting in a more consistent temperature distribution.

Figure 8 shows the temperature and the velocity distribution for the various porosities, $ϵ$ ($0.2,0.4,0.6,0.8$), with some fixed parameters. The u-velocity increased with increasing porosity in the lowest half of the indicated circle but reduced in the upper half, as shown in Figure 8a. Meanwhile, as the rate of porosity increased, the v-velocity increased at both the heated walls but declined in the center of the cavity (Figure 8b). As the fluid attempted to pass through the restricted pathways that the increased density had generated, the flow resistance caused by the higher porosity closer to the heated walls caused a rise in the velocity. The higher pressure produced by the higher porosity caused this rise in velocity. Finally, the temperature distribution also decreased as the rate of porosity increased; this scenario was similar to that shown in Figure 8c. The temperature distribution decreased as the porosity increased due to increased friction within the system. As the fluid’s porosity increased, the interactions between its molecules became more intense, increasing the friction and heat release. As a consequence, the temperature gradient was diminished, causing a reduction in the temperature distribution.

Figure 9 describes the temperature and velocity distributions for varying nanoparticle volume fractions with some constant parameters, which are mentioned in the caption. The u-velocity declined with increasing volume fraction within the $0<y<0.3$ region, but it increased in the $0.3<y<0.6$ area, as shown in Figure 9a. According to Figure 9b, the v-velocity increased as a function of $\varphi$ near the heated walls on the left- and right-sides. However, the velocity grew in the cavity’s center. Meanwhile, the temperature distribution shown in Figure 9c was diametrically opposed to that in Figure 7c and Figure 8c. In this case, the temperature distribution increased with rising $\varphi$. The increased number of molecules accessible to transfer energy improved the temperature distribution as the volume fraction increased. Since there were more molecules, more thermal energy raised the temperature as the volume fraction increased. As the number of molecules increases, a higher temperature gradient will result from more effective energy transfer.

The effect of various $Ra$ numbers on the velocity and temperature distributions is depicted in Figure 10. In the region $0<y<0.3$, the u-velocity increased for $Ra$ = ${10}^5$, but decreased for a higher $Ra$ number at ${10}^6$. On the other hand, in the region $0.3<y<0.6$, the exact opposite scenario occurred. Here, the velocity decreased for $Ra$ = ${10}^5$ but increased for $Ra$ = ${10}^6$. At $Ra={10}^5$, the buoyancy still did not fully dominate the flow. However, $Ra={10}^6$ took over control of the velocity profile. In Figure 10b, the v-velocity increased near both the heated walls as the $Ra$ number increased. As the $Ra$ number augmented, the buoyancy forces became greater than the viscous forces, resulting in more freedom for fluid mobility. The increased fluid mobility led to high velocities near the heated walls. The temperature distribution was reduced as the $Ra$ number was augmented. The decreasing rate was greater at $Ra$ = ${10}^5$, but the decreasing rate was reduced when $Ra$ = ${10}^6$ at the region $0.4\le x\le 0.6$. As per Figure 10c, the temperature increased for $Ra$ = ${10}^5$ and decreased for $Ra$ = ${10}^6$ near the heated walls. In the region between $0.2\le x\le 0.4$ and $0.6\le x\le 0.8$, the temperature distribution increased as the $Ra$ number increased.

The local $Nu$ and local $Sh$ are illustrated in Figure 11 for different $Ra$ numbers when $Da={10}^{-2}$, $D_f=0.1$, $Sr=0.2$, $Br=1$, $Le=5$, $\varphi =0.04$, and $ϵ=0.4$. Figure 11a,c for the left wall, and Figure 11b,d for the right wall. It was observed that the $Nu$ and $Sh$ for the heated walls were quite similar. This indicated that the heat and mass transfers near both side walls were identical. When a fluid is heated from below during natural convection, it becomes buoyant and elevates while the colder fluid sinks. This circulation produces a convective flow pattern known as a thermal boundary layer. The $Ra$ number affects the flow pattern and the heat transmission properties. Unlike viscous forces, buoyant forces gain strength as the $Ra$ number rises. Convection cells inside the fluid grow more prominent and active as a consequence. These cells’ convective heat transfer enhancement produces a more significant convective heat transfer coefficient. As a result, increased values of the Nusselt numbers were observed in the area. Growth in the $Ra$ number also implied an excellent temperature differential across the fluid layer. The Nusselt number increased when there was a more significant temperature differential as it encouraged improved and efficient heat transfer rates. Consequently, the $Sh$ also increased with rising $Ra$ number, as presented in Figure 11c,d. The $Ra$ number, as previously stated, indicates the ratio of buoyant to viscous strengths and is not directly connected to the Sherwood number. Diffusion, rather than convection, is the principal driver of mass transfer in natural convection. The convective mass transfer coefficient is often substantially lower than the diffusive mass transfer coefficient. As a result, the Sherwood number is controlled more by the fluid’s diffusive features, such as the diffusion coefficient, than by the convective flow patterns associated with the $Ra$ number. Therefore, the local Sherwood number in natural convection is not directly affected by the rise in the $Ra$ number. However, it is essential to remember that when convective mass transfer becomes considerable, in some circumstances, modifications in the heat transfer and flow patterns features associated with variations in the $Ra$ number may indirectly impact the mass transfer rate and, subsequently, the Sherwood number.

Figure 12 shows the local $Nu$ and local $Sh$ for various volume fractions. The $Nu$ reduces with rising volume fractions near the walls. Particle migration can occur when the volume proportion of particles close to the hot wall increases. This is because of the greater particle concentration on the heated surface, which might form particle-rich layers near the wall. The particles, which are usually less thermally conductive than the fluid, might behave as insulators, lowering the effective thermal conductivity of the combination. A lower $Nu$ was produced due to the decline in thermal conductivity, which hinders heat transfer from the heated border to the fluid. The scenario is opposite in the middle of the cavity. The local $Nu$ increased with an increasing volume fraction in that place. The local $Sh$ falls with enhancing the $\varphi$. A thicker concentration boundary barrier, decreased diffusivity, hindered fluid mixing, and concentration polarization can all result from an increased particle or solute volume fraction. These elements work together to lower the local Sherwood number, which denotes a slower mass transfer rate as the volume fraction rises. Considering these drawbacks, while optimizing the flow conditions, mixing techniques, and concentration gradients, is critical to increasing the mass transfer efficiency in such situations.

5.3. Mean Nusselt Number $\left(\overline{Nu}\right)$ and Mean Sherwood Number $\left(\overline{Sh}\right)$ for Different Parameters

The $\overline{Nu}$ and the $\overline{Sh}$ for different $Ra$ numbers, $Da$, and $Sr$ are shown in Table 6. It was found that the $\overline{Nu}$ decreased as the flow became more laminar, and the fluid velocity dropped as $Da$ decreased. The equation of $Da$ represents the flow regime in a porous medium. Consequently, the $\overline{Nu}$ also dropped as the convective heat exchange factor dropped. The $\overline{Nu}$ increased with increase in the $Ra$ number. As $Ra$ increased, the fluid flow experienced more mobility due to improved strength in the buoyancy that enhanced the heat transfer. The ratio of the convective heat transfer acceleration to the conductive heat transfer rate determines the $\overline{Nu}$, which increases with increasing $Ra$ numbers since the convective heat transfer rate rises with more significant fluid motion. The variations in $\overline{Nu}$ were significantly lower as the $Sr$ increased. The $\overline{Sh}$ increased as the $Sr$ increased. The $Sr$ number indicated a rise in the concentration gradient, which suggested a rise in the diffusivity. The concentration gradient and the diffusivity were related to the Sherwood and Soret numbers of the transferred species. The $\overline{Sh}$ also increased with increase in the $Ra$ numbers and decreased as the Darcy number decreased.

In Table 7, the $\overline{Nu}$ and the $\overline{Sh}$ are displayed for the various Dufour numbers ($D_f=0,0.1,0.2$), $\varphi =0.0,0.02,0.04$, and $Br=-2,0,2$ from a quantitative perspective. It can be seen that the $\overline{Nu}$ increased as the $\varphi$ increased as the particles were more thermally conductive. The addition of nanoparticles improved the nanofluid’s overall thermal conductivity. It was observed that when the buoyancy ratio varied from negative to positive for $D_f=0$ and $\varphi =0.0$, the ($\overline{Nu}$) increased by $82.03\%$. Furthermore, in the case of $D_f=0.1$ and $D_f=0.2$ for $\varphi =0.0$, the ($\overline{Nu}$) was enhanced by $81.37\%$ and $80.59\%$, respectively. The acceleration of growth decreased as the $\varphi$ increased. As $\varphi$ increased from 0 to 0.04 for $D_f=0$ and $Br=2$, the $\overline{Nu}$ was reduced by $5.38\%$. At $Br=2$, the decreasing rates were $7.90\%$ and $8.79\%$, respectively, when $D_f=0.1$ and $D_f=0.2$. On the other hand, the same scenario occurred for the $\overline{Sh}$. The $\overline{Sh}$ increased with enhancement in the buoyancy ratio. This could be explained by the fact that as the buoyancy forces grew more powerful than the diffusion forces, the convection mechanism became more dominant in the mass transfer. It can be observed that the $\overline{Sh}$ declined as the $\varphi$ increased, but when $D_f=0.2$ and $Br=-2$, the $\overline{Sh}$ rose. The fluid movement faced more restriction with increased $\varphi$, and the mass transfer rate plummeted concurrently.

The $\overline{Nu}$ against various Lewis numbers for different $Ra$ numbers and Darcy numbers with some constant parameters is shown in Figure 13. It was observed that the $\overline{Nu}$ decreased as the Lewis number increased. However, as the $Da$ numbers decreased and the $Ra$ numbers increased, the $\overline{Nu}$ increased. In Figure 13a, the decreasing rate of the $\overline{Nu}$ was $0.41\%$ for $Le=$ 1 to 5. For higher values of the $Le$ number, less heat will transfer across the cavity due to increasing convective motion. The $Le$ number increased, meaning more convection than conduction occurred in the enclosure. As $Le$ increased, the $\overline{Nu}$ dropped as less energy could be transported from one area to another by conduction than if only convection were in effect. Similarly, the rate of $\overline{Nu}$ decreased by $3.12\%$ and $6.81\%$ in Figure 13b,c, respectively. It was observed that the reduction rate increased as the $Da$ and $Ra$ numbers increased. The $\overline{Nu}$ increased with increase in the $Ra$ numbers, and the $\overline{Nu}$ increased as the $Da$ numbers decreased due to enhanced convective heat transfer. As the $Da$ numbers decreased and the $Ra$ numbers increased, more energy was forced into the system, causing an increase in the fluid motion due to the thermal buoyancy forces. The $\overline{Nu}$ augmented $14.23\%$, when the value of $Da$ = ${10}^{-2}$ and $Ra$ = ${10}^4$ changed to $Da$ = ${10}^{-3}$ and $Ra$ = ${10}^5$ while $Le=1$.

In the last part of this segment, Figure 14 depicts the $\overline{Sh}$ against various Lewis numbers for different $Ra$ and Darcy numbers with some fixed parameters. The $\overline{Sh}$ increased with progress in the $Le$. The $\overline{Sh}$ increased $20\%$ when the $Le$ increased from 1 to 5 in Figure 14a. The increasing rates were $83.17\%$ and $117.90\%$, as shown in Figure 14b,c. As the number of $Ra$ increases, so does the $\overline{Sh}$. As the number of $Da$ decreases, $\overline{Sh}$ increases. The $\overline{Sh}$ was highest for $Le=5$ when the $Da$ numbers were lower and the $Ra$ numbers were higher. The $Sh$ number is the convective and diffusive mass transfer ratio. It indicates the convective heat transfer performance and is usually proportional to the $Le$. Increasing the $Le$ (thermal diffusivity ratio to mass diffusivity) in a nanofluid-filled porous enclosure caused the diffusive mass transport to become relatively slower than the convective mass transport. This caused the $Sh$ to rise as the Lewis number increased. When $Da$ = ${10}^{-2}$ and $Ra$ = ${10}^4$ increased to $Da$ = ${10}^{-4}$ and $Ra$ = ${10}^6$ for $Le=5$, the mean Sherwood number augmented to $174.6\%$.

5.4. Entropy Generation

Entropy is a measurement of the disorder or unpredictability of a technique, and it is produced naturally due to irreversible processes. Table 8 shows the entropy generation with varying parameters while some parameters were kept constant. In general, the entropy generation due to fluid friction $M_F$ increases as the $Ra$ number increases but reduces as the volume fraction of the nanoparticles increases. Fluid friction in natural convection tends to rise as the $Ra$ numbers do. This could be explained by the fact that at higher $Ra$ numbers, the forces generated by buoyancy outweigh the strengths driven by viscosity, resulting in more dominant fluid mobility and stronger convective currents. As the fluid passes over solid surfaces or through porous media, these more significant currents increase the frictional resistance. The fluid friction lowers as the volume percentage of the solid particles in the fluid rises; the solid particles are mostly responsible for this phenomenon. The solid particles obstruct the flow of fluid as the volume fraction rises. As an impact, the overall fluid velocity decreases, and the flow distribution becomes spontaneous due to the reduced fluid friction. As per Table 8, the heat transfer increased for $Ra$ = ${10}^4$ and $Ra$ = ${10}^5$ but decreased for ${10}^6$ as $\varphi$ increased. After $4\%$ of the volume fraction was included, the mass transfer decreased by $16.22\%$, $15.85\%$, and $14.49\%$, while the Rayleigh number was assigned ${10}^4$, ${10}^5$, and ${10}^6$, respectively, for $D_f=0.0$. As the Dufour number increased, the mass transfer also increased. The $M_S$ increased with increase in the $Ra$ numbers but decreased with increasing $\varphi$, while $Ra$ = ${10}^5$ and ${10}^6$. For all cases, the mean Bejan numbers ($B{e}_{avg}$) were always greater than $0.5$, which indicated that irreversibility became dominant due to heat transfer.

The entropy generation due to fluid friction, $\left(M_F\right)$, heat transfer $\left(M_T\right)$ and mass transfer $\left(M_D\right)$ in a u-shaped porous enclosure for various Rayleigh and Darcy numbers while $Sr=0.2$, $D_f$ = $0.2$, $Br=1$, $Pr=6.2$ with $\varphi =0.0$ (solid line) and $0.04$ (dashed line) are shown in Figure 15. The entropy formation from fluid friction rises as the Rayleigh and Darcy numbers, respectively, rise and fall. The maximum entropy is produced near the heated walls and the cold part of the u-shaped porous enclosure when $Da$ = ${10}^{-4}$ and $Ra$ = ${10}^6$ in Figure 15c. Figure 15d–f illustrates the entropy generation due to $M_T$. The $M_T$ increases as the Rayleigh and Darcy numbers increase, respectively. The maximum $M_T$ occurs near the upper cold wall of the enclosure when $Da$ = ${10}^{-4}$ and $Ra$ = ${10}^6$. Natural convection happens when temperature differences in a fluid cause density variations. When a fluid is heated, its density reduces, causing it to become lighter and rise. The more beneficial fluid, on the other hand, grows denser and sinks. This circulation culminates in the production of convection cells, which are patterns of convective movement. The $Ra$ number substantially influences the heat transfer rate in natural convection. When the $Ra$ number grows, the buoyant forces precede the viscous forces. As an outcome, the convective flow becomes more significant and noticeable. With enhanced fluid velocity, faster-moving convection cells correspond to higher Rayleigh numbers. The increased fluid velocity ensures that heat is transported faster throughout the system. The $M_D$ is found to rise when the Rayleigh and Darcy numbers grow and decrease, respectively. After adding the $4\%$ volume fraction, the $M_F$ decreases by $38.46\%$, the $M_T$ increases by $55.85\%$, and the $M_D$ increases by $11.81\%$ for $Da$ = ${10}^{-2}$ and $Ra$ = ${10}^4$. For $Da$ = ${10}^{-3}$ and $Ra$ = ${10}^5$, the $M_F$ decreases by $26.53\%$, the $M_T$ increases by $35.89\%$, and the $M_D$ decreases by $16.05\%$. Similarly, the fluid friction decreases by $22.34\%$, the $M_T$ rises by $21.34\%$, and the $M_D$ reduces by $15.14\%$ when $Ra$ = ${10}^6$ and $Da$ = ${10}^{-4}$.

Figure 16 demonstrates the $M_s$ and the mean $B{e}_{avg}$ number for various $Ra$ numbers and Darcy numbers for $Sr=0.2$, $D_f=0.2$, $Br=1$, $Pr=6.2$ with $\varphi =0.0$ (solid line) and $\varphi =0.04$ (dashed line). The $M_s$ increases as $Ra$ and $Da$ rise and fall. Higher Rayleigh numbers are associated with more muscular buoyancy motion of the fluid, leading to increased irreversibility and energy dissipation within the system. Entropy creation measures a system’s irreversibility or loss of accessible energy. In any real-world process, there are always inefficiencies that lead to the generation of entropy. Entropy creation occurs due to irreversible processes, such as heat transfer across temperature gradients, viscous dissipation, and other dissipative effects. Overall, the combined effects of enhanced viscous dissipation within the fluid and increased entropy generation at the solid–fluid interface result in higher total entropy generation in the system as the Rayleigh number grows. The total entropy generation increases by $27.52\%$ for $Ra$ = ${10}^4$, $Da$ = ${10}^{-2}$, $8.19\%$ for $Da$ = ${10}^{-3}$, $Ra$ = ${10}^5$, and $0.48\%$ for $Da$ = ${10}^{-4}$, $Ra$ = ${10}^6$, after adding the $4\%$ volume fraction. The $B{e}_{avg}$ is almost constant as the Rayleigh and Darcy numbers rise and fall, respectively. From this finding, the $B{e}_{avg}$ is higher than 0.5, which indicates that the irreversibility becomes heat transfer dominant.

6. Response Surface Methodology

A response surface methodology (RSM) is a statistical and mathematical process for modeling and optimizing the connection between numerous independent variables (factors) and a response (dependent variable). It is frequently used in experimental design and analysis to strengthen, enhance, and comprehend intricate systems or processes. RSM uses mathematical models to approximate the link between the independent factors and the response. Higher-order models may also be utilized depending on the system’s complexity, while second-order polynomials are frequently used to express these models. RSM aims to identify the independent variables’ ideal values that produce the best or desired response. Maximizing, decreasing, or setting the reaction to a particular target value is possible. RSM aids in locating and calculating how the independent variables interact with one another. This is essential for comprehending the system’s behavior and determining how different variables affect the response.

Here, $R^2$ = $99.81$, adjusted $R^2$ = $99.43$, and predicted $R^2$ = $93.03\%$. In this RSM section, there were five parameters taken as variables, where ${10}^4\le Ra\le {10}^6$, $0.0\le D_f\le 0.2$, $-1\le Br\le 1$, $0.0\le Sr\le 0.2$ and $0.0\le \varphi \le 0.06$ were considered for this correlation analysis. The other parameters $Pr=6.2$, $Le=5$, $ϵ=0.4$ and $Da={10}^{-2}$ were kept constant for this correlation, as shown in Equation (79).

(79)$\begin{array}{cc}\hfill \overline{Nu}=& 2.05031+6.28\times {10}^{-6}Ra+2.9009D_f+0.044483Br+1.36280Sr+26.27675\varphi +3.13\times {10}^{-6}\hfill \\ \hfill & Ra{D}_f+9.34\times {10}^{-7}RaBr-2.13\times {10}^{-6}RaSr-0.000015Ra\varphi +0.534527D_{f}Br-21.76136\hfill \\ \hfill & D_{f}Sr-62.26668\varphi D_f+1.01756BrSr-2.84153Br\varphi +6.78964Sr\varphi -2.11\times {10}^{-12}R{a}^2+\hfill \\ \hfill & 23.16652{D_f}^2-0.305389B{r}^2+8.19961S{r}^2-201.00616{\varphi }^2\hfill \end{array}$

In Table 9, the model is essential given its standard F-value of $261.81$. A value of F this large might occur due to noise just $0.01\%$ of the time. When the p-value for a model term is less than $0.05$, it is deemed significant [48]. In this case, the essential model terms are $Ra$, $D_f$, $Br$, $Ra\times D_f$,$Ra\times Br$, $Ra\times Sr$, $Ra\times \varphi$, $D_f\times Sr$, $D_f\times \varphi$, $Br\times Sr$, $R{a}^2$, ${D_f}^2$, $B{r}^2$, and ${\varphi }^2$. The model terms are not considered significant if the value exceeds $0.1000$. Model reduction may improve a model if it contains a lot of unnecessary words (except those needed to maintain hierarchy). The lack of fit F-value of $180.98$ showed the significance of the lack of fit. An elevated lack of fit F-value has a $0.01\%$ chance of being brought on by noise. To describe the relation between the response function and the input factors ($Ra$, $D_f$, $Sr$, $Br$, and $\varphi$), consider the following mathematical connection:

(80)$\begin{array}{cc}\hfill \overline{Nu}=& 2.05031+6.28\times {10}^{-6}Ra+2.9009D_f+0.044483Br+3.13\times {10}^{-6}Ra{D}_f+9.34\times {10}^{-7}RaBr\hfill \\ \hfill & -2.13\times {10}^{-6}RaSr-0.000015Ra\varphi -21.76136D_{f}Sr-62.26668\varphi D_f+1.01756BrSr\hfill \\ \hfill & -2.11\times {10}^{-12}R{a}^2+23.16652{D_f}^2-0.305389B{r}^2-201.00616{\varphi }^2\hfill \end{array}$

Response Surface Analysis

Figure 17 shows four different residual plots developed after entering the data into statistical software and performing a variance analysis (ANOVA). The degree of accuracy of fit in regression is measured using the residual graphs. The residuals in Figure 17a show the standard probability plots were very acceptable. Locating any predictable biases in the model is aided by a percent-versus-residual graph. It can reveal over-prediction or under-prediction in a variety of operational circumstances. Over-prediction is the positive percentage of errors, and under-prediction is the harmful percentage. It gives an impression of the regions where the model must be improved or changed to increase accuracy. Figure 17b displays the residual versus the observation plot. The regression model performance can be assessed using this graph, which shows the residuals against the order of the observations. Using the residuals can assist in identifying patterns and trends that can enhance the model or identify data problems. This graph’s main objective is to look for regular patterns or correlations in the residuals caused by the data’s temporal order. The residuals should have a random distribution centered on zero to demonstrate the model’s accuracy and be unaffected by the data recording order. Residuals that change in magnitude as the observations progress indicate a systematic trend or an inaccuracy in the model’s estimates. This pattern implies that the model should consider any sequential or temporal elements that impact $\overline{Nu}$. Looking for outliers or sudden modifications in the residuals at specific statement sites can help detect a clear difference between the expected and the actual values. To enhance the $\overline{Nu}$ forecast’s accuracy, the model may require updating, considering additional variables, or acknowledging temporal dependencies. For all the responses in our simulation, we observed that the most significant residuals for $\overline{Nu}$ were close to $11.5$. The difference between the expected and actual results for $\overline{Nu}$ is seen in Figure 17c. Analyzing the predicted vs. accurate graph will determine the prediction system’s dependability and accuracy. The residuals were found to be close to or on a straight line, demonstrating the regression models’ accuracy in fitting the data. The typical probability plots of the residual distributions are shown to evaluate the observation’s normality. Given that this line is straight, the $\overline{Nu}$ residual distribution is assumed to be predictable. Figure 17d shows an irregular distribution of the residual histograms that does not approximate a symmetrical distribution. The frequency initially rises to 4.03133 and then falls. Between $4.99029$ and $5.94926$, the frequency again rises. After increasing the frequency, it again reduces with increase in $\overline{Nu}$. The maximum frequency was observed between the $\overline{Nu}$ $3.07236$ and $4.03133$.

A mathematical model or outline plot is a two-dimensional visible representation of result variables as a process involving two separate factors. The response surface plotting functions are similar to 3D character plotting but present the response variable as a collection of contour lines in a 2D texture. A model contour plot displays the independent variables on the x- and y-axes, while a series of contours represents the response variable. The model contour and response surface plots for illustrating the impact of various self-dependent factors on $\overline{Nu}$ are shown in Figure 18.

The response function $\overline{Nu}$ with $Ra$ and $\varphi$ for various Dufour numbers is shown in Figure 18 with some fixed parameters. The relationship between $Ra$ and $\varphi$ for the response function $\overline{Nu}$ is analyzed through 2D and 3D graphs, maintaining constant $D_f$. These $2D$ and $3D$ graphs indicate that the $\overline{Nu}$ is maximum at the highest $Ra$ number, and as the $\varphi$ increases, the $\overline{Nu}$ also increases. The system’s improved convective heat transfer indicates superior buoyancy-driven forces over conduction and increased heat exchange interfaces due to higher volume fractions. The $\overline{Nu}$ was slightly higher as the $D_f$ increased.

7. Sensitivity Analysis

Sensitivity analysis in the natural convection of nanofluids refers to systematic research on how modifications to various parameters and variables might impact the heat transfer characteristics of a nanofluid (a fluid containing nanoparticles). Understanding the effects of uncertainties and fluctuations in the input parameters on the system’s output is made possible using this helpful instrument in scientific study and engineering. Researchers and engineers may make wise design and application decisions using sensitivity analysis to better understand the nanofluid’s behavior and to improve its performance, and to pinpoint the essential elements that significantly impact the nanofluid’s natural convection heat transfer mechanism. These variables might be the kind and quantity of nanoparticles, the fluid characteristics (like the viscosity and the thermal conductivity), the size and form of the heated surface, the boundary conditions (like the temperature gradients), and the external variables (like gravity). Sensitivity analysis involves purposefully changing certain factors one at a time while leaving the others constant. Possible real-world values or extreme scenarios may be included in the range of variation. Researchers may, for instance, vary the thermal conductivity of the fluid, the concentration of the nanoparticles, or the temperature differential across the heated surface. A sensitivity analysis examines how parameter changes affect the heat transfer properties, including the temperature distribution in the nanofluid, the convective heat transfer coefficient, and the heat transfer rate. This procedure aids in determining which variables significantly affect how the system behaves. Individual partial derivatives of the dependent variable $\left(Nu\right)$ concerning the independent variables $Ra$, $D_f$, $Br$, $Sr$, and $\varphi$ are calculated in this work. We were able to obtain the intended results by considering the values of the independent components at three levels: $-1$ (low), 0 (medium), and 1 (high). The sensitivity is calculated using the regression Equation (80). The partial derivative of the function that produces the result ($\overline{Nu}$) concerning the independent variables is used to determine the sensitivity. The response function’s sensitivity functions are as follows:

(81)$\begin{array}{cc}\hfill \frac{\partial \overline{Nu}}{\partial Ra}=& 6.28\times {10}^{-6}+3.13\times {10}^{-6}{D}_f+9.34\times {10}^{-7}Br-2.13\times {10}^{-6}Sr-0.000015\varphi \hfill \\ \hfill & -4.22\times {10}^{-12}Ra\hfill \end{array}$

(82)$\begin{array}{c}\hfill \frac{\partial \overline{Nu}}{\partial D_f}=2.9009+3.13\times {10}^{-6}Ra-21.76136Sr-61.26668\varphi +46.33304D_f\end{array}$

(83)$\begin{array}{c}\hfill \frac{\partial \overline{Nu}}{\partial Br}=0.044483+9.34\times {10}^{-7}Ra+1.01756Sr-0.610778Br\end{array}$

(84)$\begin{array}{c}\hfill \frac{\partial \overline{Nu}}{\partial Sr}=-2.13\times {10}^{-6}Ra-21.76136D_f+1.01756Br\end{array}$

(85)$\begin{array}{c}\hfill \frac{\partial \overline{Nu}}{\partial \varphi }=-0.000015Ra-62.26668D_f-402.01232\varphi \end{array}$

A positive sensitivity output indicates that increasing the input variable’s value improves the output operation’s value. A negative sensitivity value also indicates a decrease in the output factor caused by increasing the input factor. The sensitivity analysis for various $Ra$, $D_f$, $Br$, $Sr$ and $\varphi$ quantities on the output response $Nu$ is shown in Table 10 and Figure 19. It was observed that the sensitivity metrics of $Ra$ and $\varphi$ were positive when the volume fraction value was the lowest or medium. However, the sensitivity was negative when the volume fraction was the highest. and positive sensitivity was observed for $Br$ and $Sr$. When the Soret number is raised, the Soret effect becomes more significant than thermal diffusion. In some cases, increasing the Soret number can result in improved separation or concentration gradients of the mixture components, resulting in higher heat and mass transfer rates. This can lead to positive sensitivity, where increasing the Soret number enhances the total transport efficiency. Increasing the buoyancy ratio can lead to stronger natural convection currents in the fluid. This increased fluid motion can improve the heat and mass transfer rates, resulting in positive sensitivity. $D_f$ has a positive sensitivity when the $Sr$ is medium; otherwise, in all cases, the sensitivity of $D_f$ is negative.

8. Conclusions

The study explores the impact of the Dufour and Soret effects on the thermosolutal natural convection of a Newtonian $A{l}_2{O}_3$-$H_{2}O$ nanofluid in a U-shaped porous enclosure. Varying parameters were investigated in the study: the Rayleigh number ($Ra={10}^4,{10}^5,{10}^6$), the Darcy number ($Da={10}^{-2},{10}^{-3},{10}^{-4}$), the Soret number ($Sr=0.0,0.1,0.2$), the Dufour number ($D_f=0.0,0.1,0.2$), the buoyancy ratio ($-2\le Br\le 2$), the Lewis number ($Le=1,3,5$), the volume fraction ($0\le \varphi \le 0.04$), and the porosity $ϵ=(0.2-0.8$). The numerical findings were represented in graphical and tabular form to provide quantitative and physical information regarding the generated solutions. The RSM provided the 2D contour and 3D surface plots that are displayed. The RSM method was also used to investigate the sensitivity of the $\overline{Nu}$ to the input parameters. The following are the principal findings of this study:

• The fluid flow is enhanced as the Rayleigh number and buoyancy force increase. The isotherm and isoconcentration densities close to the hot wall grow when the buoyancy force shifts from a negative to a positive value.

• The u-velocity at mid-x rises between $0\le y\le 0.3$ but decreases between $0.3\le y\le 0.6$ by growing the $Br$. The v-velocity at mid-y enhances close to the walls and reduces in the middle of the cavity for rising $Br$ and $ϵ$. The temperature distribution also decreases by increasing $Br$ and $ϵ$.

• The v-velocity at mid-y decreases near the walls but enhances at the down-indented part of the U-shape with increasing $\varphi$. The temperature distribution rises as the $\varphi$ rises.

• The temperature distribution at the down-indented portion of the U-shape and the v-velocity at mid-y decrease as the $Ra$ number rises. But they increase close to the walls as the $Ra$ number rises.

• The local $Nu$ and the local $Sh$ increase as the $Ra$ number rises but reduce as the volume fraction increases.

• The $\overline{Nu}$ and $\overline{Sh}$ decrease as the Darcy number falls, but both rise as the $Ra$ number increases. The $\overline{Nu}$ and the $\overline{Sh}$ rise as the buoyancy ratio shifts from a negative to a positive value.

• As the $Le$ number reduces, the $\overline{Nu}$ reduces, but the $\overline{Sh}$ rises.

• The fluid friction $M_F$ rises as the $Ra$ number rises, but reduces as the volume fraction rises. The heat transfer $M_T$ increases for $Ra$ = ${10}^4$ and ${10}^5$, but decreases for $Ra$ = ${10}^6$ as the volume fraction rises.

• The $B{e}_{avg}$ is always more than $0.5$, which indicates the irreversibility becomes heat transfer dominant.

• The generated correlation equation from the RSM approach demonstrates the correlation between the output and the input parameters.

• The $Br$ and $Sr$ have a positive sensitivity, but $Ra$ and $\varphi$ have a negative sensitivity only for higher values of $\varphi$.

Author Contributions

Conceptualization, M.M.I. and M.M.M.; methodology, M.M.I. and M.M.M.; software, M.M.I. and M.M.M.; validation, M.M.I. and M.M.M.; formal analysis, M.M.I. and M.F.H.; investigation, M.M.I., M.F.H. and M.M.M.; resources, M.M.I., M.F.H. and M.M.M.; data curation, M.M.I. and M.M.M.; writing—original draft preparation, M.M.I. and M.F.H.; writing—review and editing, M.F.H. and M.M.M.; visualization, M.M.I., M.F.H. and M.M.M.; supervision, M.M.M.; project administration, M.M.M.; funding acquisition, M.M.M. All authors have read and agreed to the published version of the manuscript.

Data Availability Statement

The data are available by request.

Conflicts of Interest

The authors declare no conflict of interest.

Abbreviations

The following abbreviations are used in this manuscript:

Abbreviations
ANOVA	Analysis of variance
BFEDM	Brinkman–Forchheimer extended Darcy model
FEM	Finite element method
FVM	Finite volume method
LBM	Lattice Boltzmann method
MHD	Magnetohydrodynamics
MRT	Multiple-relaxation-time
RSM	Response surface methodology
Nomenclature
English symbols
$Be$	Local Bejan number
$B{e}_{avg}$	Mean Bejan number
$Br$	Buoyancy ratio
C	Concentration of the fluid	kg m${}^3$
$C_p$	Specific heat	J kg${}^{-1}$ K${}^{-1}$
D	Mass diffusivity	m${}^2$ s${}^{-1}$
$Da$	Darcy number
$D_f$	Dufour number
g	Gravitational acceleration	m s${}^{-2}$
H	Height of the cavity	m
L	Length of the cavity
$Le$	Lewis number
k	Thermal conductivity	J m${}^{-1}$ s${}^{-1}$ K${}^{-1}$
$Nu$	Local Nusselt number
$\overline{Nu}$	Mean Nusselt number
$\overline{P}$	Dimensional pressure
P	Dimensionless pressure
$Pr$	Prandtl number
$Ra$	Rayleigh number
M	Entropy generation	J m${}^{-3}$ s${}^{-1}$ K${}^{-1}$
$M_F$	Entropy due to fluid friction
$M_T$	Entropy due to heat transfer
$M_D$	Entropy due to mass transfer
$M_S$	Total entropy
$Sr$	Soret number
$Sh$	Sherwood number
$\overline{Sh}$	Mean Sherwood number
T	Temperature	K
$T_c$	Temperature of the cold wall	(K)
$T_h$	Temperature of the heated wall
t	Dimensional time	s
$\overline{t}$	Dimensionless time
$\overline{u}$,$\overline{v}$	Dimensional velocity components along horizontal and	m s${}^{-1}$
	vertical directions
$u,v$	Dimensionless velocity components
$\overline{x},\overline{y}$	Dimensional Cartesian coordinates	m
$x,y$	Dimensionless Cartesian coordinates
Greek letters
$\alpha$	Thermal diffusivity	m${}^2$ s${}^{-1}$
$\beta$	Thermal expansion coefficient	K${}^{-1}$
$\theta$	Dimensionless temperature function
$\mu$	Dynamic viscosity	kg m${}^{-1}$ s${}^{-1}$
$\nu$	Kinematic viscosity	m${}^2$ s${}^{-1}$
$\rho$	Fluid density	kg m${}^{-3}$
$\varphi$	Volume fraction
$ϵ$	Porosity

Footnotes

Disclaimer/Publisher’s Note: The statements, opinions and data contained in all publications are solely those of the individual author(s) and contributor(s) and not of MDPI and/or the editor(s). MDPI and/or the editor(s) disclaim responsibility for any injury to people or property resulting from any ideas, methods, instructions or products referred to in the content.

Figures and Tables

Figure 1. (a) Model diagram with the coordinate system, (b) mesh composition.

Figure 2. Qualitative comparison of the present MRT-LBM results (b) isotherms and (d) isoconcentrations with the results of Ren and Chan [37] (a) isotherms and (c) isoconcentrations while [Forumla omitted. See PDF.], [Forumla omitted. See PDF.], [Forumla omitted. See PDF.], [Forumla omitted. See PDF.], [Forumla omitted. See PDF.], [Forumla omitted. See PDF.], [Forumla omitted. See PDF.], [Forumla omitted. See PDF.] and [Forumla omitted. See PDF.]. (Reproduced with permission from [37], Elsevier, 2016).

Figure 3. Streamlines (top), isotherms (middle), and isoconcentrations (bottom) for [Forumla omitted. See PDF.] and [Forumla omitted. See PDF.] numbers (a,d,g) [Forumla omitted. See PDF.], [Forumla omitted. See PDF.], (b,e,h) [Forumla omitted. See PDF.], [Forumla omitted. See PDF.], (c,f,i) [Forumla omitted. See PDF.], [Forumla omitted. See PDF.] while [Forumla omitted. See PDF.], [Forumla omitted. See PDF.], [Forumla omitted. See PDF.], [Forumla omitted. See PDF.] and [Forumla omitted. See PDF.] = 0.0 (solid line), 0.04 (long-dash line).

Figure 4. Streamlines for [Forumla omitted. See PDF.], [Forumla omitted. See PDF.] (a) [Forumla omitted. See PDF.], (b) [Forumla omitted. See PDF.] and (c) [Forumla omitted. See PDF.] and for [Forumla omitted. See PDF.], [Forumla omitted. See PDF.] (d) [Forumla omitted. See PDF.], (e) [Forumla omitted. See PDF.] and (f) [Forumla omitted. See PDF.] and for [Forumla omitted. See PDF.], [Forumla omitted. See PDF.] (g) [Forumla omitted. See PDF.], (h) [Forumla omitted. See PDF.] and (i) [Forumla omitted. See PDF.] while [Forumla omitted. See PDF.], [Forumla omitted. See PDF.], [Forumla omitted. See PDF.] and [Forumla omitted. See PDF.] = 0.0 (solid line), 0.04 (long-dash line).

Figure 5. Isotherms for [Forumla omitted. See PDF.], [Forumla omitted. See PDF.] (a) [Forumla omitted. See PDF.], (b) [Forumla omitted. See PDF.] and (c) [Forumla omitted. See PDF.] and for [Forumla omitted. See PDF.], [Forumla omitted. See PDF.] (d) [Forumla omitted. See PDF.], (e) [Forumla omitted. See PDF.] and (f) [Forumla omitted. See PDF.] and for [Forumla omitted. See PDF.], [Forumla omitted. See PDF.] (g) [Forumla omitted. See PDF.], (h) [Forumla omitted. See PDF.] and (i) [Forumla omitted. See PDF.] while [Forumla omitted. See PDF.], [Forumla omitted. See PDF.], [Forumla omitted. See PDF.] and [Forumla omitted. See PDF.] = 0.0 (solid line), 0.04 (long-dash line).

Figure 6. Isoconcentrations for [Forumla omitted. See PDF.], [Forumla omitted. See PDF.] (a) [Forumla omitted. See PDF.], (b) [Forumla omitted. See PDF.] and (c) [Forumla omitted. See PDF.] and for [Forumla omitted. See PDF.], [Forumla omitted. See PDF.] (d) [Forumla omitted. See PDF.], (e) [Forumla omitted. See PDF.] and (f) [Forumla omitted. See PDF.] and for [Forumla omitted. See PDF.], [Forumla omitted. See PDF.] (g) [Forumla omitted. See PDF.], (h) [Forumla omitted. See PDF.] and (i) [Forumla omitted. See PDF.] while [Forumla omitted. See PDF.], [Forumla omitted. See PDF.], [Forumla omitted. See PDF.] and [Forumla omitted. See PDF.] = 0.0 (solid line), 0.04 (long-dash line).

Figure 7. Analysis of velocity and temperature profiles: (a) u-velocity at mid x, (b) v-velocity at mid y, (c) temperature distribution at mid y for various [Forumla omitted. See PDF.] when [Forumla omitted. See PDF.], [Forumla omitted. See PDF.], [Forumla omitted. See PDF.], [Forumla omitted. See PDF.], [Forumla omitted. See PDF.] and [Forumla omitted. See PDF.].

Figure 8. Variations in velocity and temperature profiles: (a) u-velocity at mid x, (b) v-velocity at mid y, (c) temperature distribution at mid y for various porosity, [Forumla omitted. See PDF.], when [Forumla omitted. See PDF.], [Forumla omitted. See PDF.], [Forumla omitted. See PDF.], [Forumla omitted. See PDF.], [Forumla omitted. See PDF.] and [Forumla omitted. See PDF.].

Figure 9. Changes in velocity and temperature profiles: (a) u-velocity at mid x, (b) v-velocity at mid y, (c) temperature distribution at mid y for various [Forumla omitted. See PDF.] when [Forumla omitted. See PDF.], [Forumla omitted. See PDF.], [Forumla omitted. See PDF.], [Forumla omitted. See PDF.], [Forumla omitted. See PDF.], [Forumla omitted. See PDF.] and [Forumla omitted. See PDF.].

Figure 10. Variations in velocity and temperature profiles: (a) u-velocity at mid x, (b) v-velocity at mid y, (c) temperature distribution at mid y for various [Forumla omitted. See PDF.] when [Forumla omitted. See PDF.], [Forumla omitted. See PDF.], [Forumla omitted. See PDF.], [Forumla omitted. See PDF.], [Forumla omitted. See PDF.], [Forumla omitted. See PDF.] and [Forumla omitted. See PDF.].

Figure 11. Local Nusselt number ([Forumla omitted. See PDF.]) for (a) left wall ([Forumla omitted. See PDF.]) and (b) right wall ([Forumla omitted. See PDF.]) and local Sherwood number ([Forumla omitted. See PDF.]) for (c) left wall ([Forumla omitted. See PDF.]) and (d) right wall ([Forumla omitted. See PDF.]) for different Rayleigh numbers, ([Forumla omitted. See PDF.]) while [Forumla omitted. See PDF.], [Forumla omitted. See PDF.], [Forumla omitted. See PDF.], [Forumla omitted. See PDF.], [Forumla omitted. See PDF.], [Forumla omitted. See PDF.] and [Forumla omitted. See PDF.].

Figure 12. Local Nusselt number ([Forumla omitted. See PDF.]) for (a) left wall ([Forumla omitted. See PDF.]) and (b) right wall ([Forumla omitted. See PDF.]) and local Sherwood number ([Forumla omitted. See PDF.]) for (c) left wall ([Forumla omitted. See PDF.]) and (d) right wall ([Forumla omitted. See PDF.]) for different volume fractions, [Forumla omitted. See PDF.] while [Forumla omitted. See PDF.] [Forumla omitted. See PDF.], [Forumla omitted. See PDF.], [Forumla omitted. See PDF.], [Forumla omitted. See PDF.], [Forumla omitted. See PDF.] and [Forumla omitted. See PDF.].

Figure 13. Mean Nusselt number ([Forumla omitted. See PDF.]) against various Lewis numbers ([Forumla omitted. See PDF.]) for (a) [Forumla omitted. See PDF.], [Forumla omitted. See PDF.], (b) [Forumla omitted. See PDF.], [Forumla omitted. See PDF.] and (c) [Forumla omitted. See PDF.], [Forumla omitted. See PDF.] while [Forumla omitted. See PDF.], [Forumla omitted. See PDF.], [Forumla omitted. See PDF.], [Forumla omitted. See PDF.], [Forumla omitted. See PDF.] and [Forumla omitted. See PDF.].

Figure 14. Mean Sherwood number ([Forumla omitted. See PDF.]) against various Lewis numbers ([Forumla omitted. See PDF.]) for (a) [Forumla omitted. See PDF.], [Forumla omitted. See PDF.], (b) [Forumla omitted. See PDF.], [Forumla omitted. See PDF.] and (c) [Forumla omitted. See PDF.], [Forumla omitted. See PDF.] while [Forumla omitted. See PDF.], [Forumla omitted. See PDF.], [Forumla omitted. See PDF.], [Forumla omitted. See PDF.], [Forumla omitted. See PDF.] and [Forumla omitted. See PDF.].

Figure 15. Entropy formation contours showing fluid friction, [Forumla omitted. See PDF.] (top), heat transfer, [Forumla omitted. See PDF.] (middle), and mass transfer, [Forumla omitted. See PDF.] (bottom) for (a,d,g) [Forumla omitted. See PDF.] and [Forumla omitted. See PDF.], (b,e,h) [Forumla omitted. See PDF.] and [Forumla omitted. See PDF.], (c,f,i) [Forumla omitted. See PDF.] and [Forumla omitted. See PDF.] while [Forumla omitted. See PDF.], [Forumla omitted. See PDF.], [Forumla omitted. See PDF.], [Forumla omitted. See PDF.] with [Forumla omitted. See PDF.] (solid line) and [Forumla omitted. See PDF.] (dashed line).

Figure 16. Contours of total entropy generation, [Forumla omitted. See PDF.] (top), and Bejan number, [Forumla omitted. See PDF.] (bottom) for (a,d) [Forumla omitted. See PDF.] and [Forumla omitted. See PDF.] (b,e) [Forumla omitted. See PDF.] and [Forumla omitted. See PDF.] (c,f) [Forumla omitted. See PDF.] and [Forumla omitted. See PDF.] while [Forumla omitted. See PDF.], [Forumla omitted. See PDF.], [Forumla omitted. See PDF.], [Forumla omitted. See PDF.] with [Forumla omitted. See PDF.] (solid line) and [Forumla omitted. See PDF.] (dashed line).

Figure 17. Residual graph for response function [Forumla omitted. See PDF.]. (a) Normal probability graph, (b) residual vs. observation graph, (c) predicted vs. actual graph, and (d) histogram of residuals.

Figure 18. The different forms of [Forumla omitted. See PDF.] for [Forumla omitted. See PDF.], [Forumla omitted. See PDF.] and [Forumla omitted. See PDF.] as a function of the effective parameters ([Forumla omitted. See PDF.]) when (a,b) [Forumla omitted. See PDF.], (c,d) [Forumla omitted. See PDF.], and (e,f) [Forumla omitted. See PDF.] while [Forumla omitted. See PDF.], [Forumla omitted. See PDF.] and [Forumla omitted. See PDF.].

Figure 19. Sensitivity analysis of [Forumla omitted. See PDF.] at different coded values −1 (low), 0 (medium), 1 (high) of volume fraction while [Forumla omitted. See PDF.], [Forumla omitted. See PDF.], [Forumla omitted. See PDF.] and [Forumla omitted. See PDF.].

Table 1

Some relevant studies with regard to the present investigation.

Ref.	NumericalMethods	FlowDomain	StudyMethod	Parametersand Ranges	Grid Size
Esfe et al. [23]	FVM,	U-shaped	Free	${10}^3\le Ra\le {10}^5$,	2100 cells
	Ansys	enclosure	convection,	$0.0\le \varphi \le 0.03$,
	Fluent,		porous media,	$Da=0,60,$
			laminar	$A{l}_2{O}_3/H_{2}O$
				nanofluid
Ali et al. [20]	FEM	U-shaped	Natural	${10}^3\le Ra\le {10}^6$,	$300\times 300$
		cavity	convection	$0.6\le n\le 1.4$,
				$0\le \varphi \le 0.1$,
				$0.2\le AR\le 0.6$,
				$F{e}_3{O}_4/H_{2}O$
				nanofluid
Liu et al. [10]	LBM	Rectangular	Double-	$-0.2\le D_f\le 0.2$,
		enclosure	diffusive	$-0.5\le Sr\le 0.5$,	$200\times 400$
			convection,	$Pr=1$, $Ra={10}^5$
			laminar
Sajjadi et al. [24]	LBM	Square	MHD natural	$Da={10}^{-2}-{10}^{-1}$,
		cavity	convection,	${10}^3\le Ra\le {10}^5$,	$160\times 160$
			porous media	$0.4\le ϵ\le 0.9$,
				$0.0\le \varphi \le 0.003$,
				$0\le Ha\le 50$,
				MWCNT-
				$F{e}_3{O}_4/H_{2}O$,
				Hybrid nanofluid
Xu et al. [11]	LBM	Square	Double-diffusive	$0.0\le Sr\le 0.8$,	$200\times 200$
		enclosure	natural convection	$0.0\le D_f\le 0.8$,
				$-10\le Br\le 10$,
				$Pr=1$, $Ra={10}^5$
Saha et al. [25]	Second	Trapezoidal	MHD double-	$Ra={10}^4,{10}^5,{10}^6$,	$81\times 81$
	order	enclosure	diffusive mixed	$Ha=0,40,80,120$,
	FDM		convection	$\varphi =0-0.1$, $\gamma =0-\pi$,
				$ϵ=0,0.3,0.7,1$.
Nithyadevi  [26]	Control	Square	Double-diffusive	${10}^3\le R{a}_T\le {10}^6$,	$51\times 51$
	volume	enclosure	natural convection	$0.8\le N\le 1.2$,
	method			$0\le R\le 2$,
				$1\le Sc\le 5$,
				$0\le Sr\le 1$,
				$0\le D_f\le 1$
He et al. [27]	LBM	Square	Natural	${10}^3\le Ra\le {10}^5$,	$100\times 100$
		cavity	convection	$0.0\le \varphi \le 0.08$,
				$0\le \zeta \le \pi$,
				${10}^{-2}\le Kr\le {10}^2$

Table 2

Physical properties of the base fluid and $A{l}_2{O}_3$ nanoparticles [30,31,32].

Physical Properties	Base Fluid (${\mathit{H}}_2\mathit{O}$)	${\mathit{Al}}_2{\mathit{O}}_3$
$c_p$ [J kg${}^{-1}$ K${}^{-1}]$	4179	765
$\rho$ [kg m${}^{-3}]$	997.1	3970
k [W m${}^{-1}$ k${}^{-1}]$	0.613	40
$\beta \times {10}^{-5}$ [K${}^{-1}]$	21	0.85
$\mu$ [kg m${}^{-1}$ $s^{-1}]$	$9.1\times {10}^{-4}$	-

Table 3

Grid independent test in terms of the mean Nusselt number ($\overline{Nu}$), and the mean Sherwood number ($\overline{Sh}$), for $Pr=6.2$ across different Soret numbers ($Sr$) while $Ra={10}^6$, $Da={10}^{-4}$, $Le=5$, $\varphi =0.04$, $Br=1$, $D_f=0.2$ and $ϵ=0.4$.

	Mean Nusselt Number ($\overline{\mathit{Nu}}$)			Mean Sherwood Number ($\overline{\mathit{Sh}}$)
Lattice Size	$\mathit{Sr}=\mathbf{0.0}$	$\mathit{Sr}=\mathbf{0.1}$	$\mathit{Sr}=\mathbf{0.2}$	$\mathit{Sr}=\mathbf{0.0}$	$\mathit{Sr}=\mathbf{0.1}$	$\mathit{Sr}=\mathbf{0.2}$
$128\times 128$	4.1052	3.9712	3.9689	7.9802	8.1335	8.4913
$256\times 256$	3.9041	3.8843	3.9689	7.7580	7.9398	8.3437
$512\times 512$	3.8729	3.8336	3.9689	7.6921	7.8811	8.2904

Table 4

Code validation by comparing data of [46,47]. The mean Nusselt number, $\overline{Nu}$, with the BFEDM for $ϵ=0.4$, $\varphi =0.0$, and $Pr=1.0$.

$\mathit{Da}$	$\mathit{Ra}$	Nithiarasu et al. [46]	Guo and Zhao [47]	Present
${10}^{-2}$	${10}^3$	$1.010$	$1.008$	$1.0069$
	${10}^4$	$1.408$	$1.367$	$1.3393$
	${10}^5$	$2.983$	$2.998$	$2.9521$
${10}^{-4}$	${10}^5$	1.067	1.066	1.0579
	${10}^6$	2.550	2.603	2.6263

Table 5

Comparison of the present MRT-LBM results in terms of the mean Nusselt number ($\overline{Nu}$) and the mean Sherwood number ($\overline{Sh}$) with the results of Xu et al. [11] and Ren et al. [37], where $A=H/L=2.0$, $Pr=1.0$, $Le=2.0$, $Sr=0.1$, $Dr=0.1$, $Br=-2$, $\varphi =0.0$, $ϵ=1$ and $Da={10}^9$.

	Mean Nusselt Number ($\overline{\mathit{Nu}}$)			Mean Sherwood Number ($\overline{\mathit{Sh}}$)
$\mathit{Ra}$	Xu et al. [11]	Ren et al. [37]	Present	Xu et al. [11]	Ren et al. [37]	Present
${10}^3$	1.2027	1.2020	1.2055	1.4654	1.4683	1.4879
${10}^4$	1.7443	1.7407	1.7179	2.5766	2.5943	2.5374
${10}^5$	2.8976	2.9153	2.9238	4.6108	4.6993	4.7109
${10}^6$	5.0852	5.1808	5.1688	8.3294	8.5134	8.4593

Table 6

The mean Nusselt number ($\overline{Nu}$), and mean Sherwood number ($\overline{Sh}$), for $Pr=6.2$ across different Rayleigh numbers ($Ra$), Darcy numbers ($Da$), and Soret numbers ($Sr$) while $Le=5$, $\varphi =0.04$, $Br=1$, $D_f=0.2$ and $ϵ=0.4$.

		Mean Nusselt Number ($\overline{\mathit{Nu}}$)			Mean Sherwood Number ($\overline{\mathit{Sh}}$)
$\mathit{Ra}$	$\mathit{Da}$	$\mathit{Sr}=\mathbf{0.0}$	$\mathit{Sr}=\mathbf{0.1}$	$\mathit{Sr}=\mathbf{0.2}$	$\mathit{Sr}=\mathbf{0.0}$	$\mathit{Sr}=\mathbf{0.1}$	$\mathit{Sr}=\mathbf{0.2}$
	${10}^{-1}$	3.1901	3.1754	3.1754	2.8634	3.1279	3.4280
${10}^4$	${10}^{-2}$	3.1735	3.1720	3.1658	2.7396	3.0246	3.2989
	${10}^{-3}$	3.1533	3.1515	3.1544	2.2806	2.5394	2.6936
	${10}^{-4}$	3.1499	3.1512	3.1520	2.1979	2.4648	2.7470
	${10}^{-1}$	4.2372	4.2656	4.2710	6.3332	6.5770	6.9152
${10}^5$	${10}^{-2}$	4.0901	4.1275	4.1257	6.1286	6.3911	6.7514
	${10}^{-3}$	3.3910	3.3860	3.4887	5.0495	5.3723	5.7572
	${10}^{-4}$	3.1482	3.1456	3.1498	2.4609	2.7390	3.0472
	${10}^{-1}$	8.3345	8.3160	8.3493	11.6521	12.2481	12.8763
${10}^6$	${10}^{-2}$	8.1715	8.1556	8.1821	11.5961	12.1203	12.7423
	${10}^{-3}$	7.1137	7.1243	7.1427	10.8865	11.3945	11.9208
	${10}^{-4}$	3.9041	3.8843	3.9689	7.7580	7.9398	8.3437

Table 7

The mean Nusselt number ($\overline{Nu}$), and mean Sherwood number ($\overline{Sh}$), for $Pr=6.2$ across different Dufour numbers ($D_f$), buoyancy ratios ($Br$) and volume fractions ($\varphi$) while $Le=5$, $ϵ=0.4$, $Ra={10}^6$, $Da={10}^{-4}$ and $Sr=0.2$.

		Mean Nusselt Number ($\overline{\mathit{Nu}}$)			Mean Sherwood Number ($\overline{\mathit{Sh}}$)
${\mathit{D}}_{\mathit{f}}$	$\mathit{Br}$	$\mathit{\varphi }=\mathbf{0.0}$	$\mathit{\varphi }=\mathbf{0.02}$	$\mathit{\varphi }=\mathbf{0.04}$	$\mathit{\varphi }=\mathbf{0.0}$	$\mathit{\varphi }=\mathbf{0.02}$	$\mathit{\varphi }=\mathbf{0.04}$
	$-2$	2.2550	2.6084	2.8082	4.8906	4.8280	4.6701
0	0	3.0187	3.0996	3.0270	7.3294	6.8599	5.5523
	2	4.1048	4.0648	3.8840	11.0281	10.4693	9.6477
	$-2$	2.4319	2.7955	2.9959	4.9429	4.87796	4.6419
$0.1$	0	3.2595	3.3191	3.2530	7.2247	6.9384	6.1021
	2	4.4107	4.2971	4.0623	11.0223	10.4048	9.5648
	$-2$	2.6068	2.9882	3.2209	5.0754	5.3167	5.5715
$0.2$	0	3.4705	3.4584	3.3663	7.2541	6.8339	5.5646
	2	4.7076	4.5443	4.2937	11.0216	10.3454	9.5161

Table 8

Total entropy for $Pr=6.2$ across different Dufour numbers ($D_f$), Rayleigh numbers ($Ra$) and volume fractions ($\varphi$) while $Le=5$, $Br=2$, $Sr=0.2$, $ϵ=0.4$, $Da={10}^{-2}$ when $Ra={10}^4$, $Da={10}^{-3}$ when $Ra={10}^5$ and $Da={10}^{-4}$ when $Ra={10}^6$.

${\mathit{D}}_{\mathit{f}}$	$\mathit{Ra}$	$\mathit{\varphi }$	${\left({\mathit{M}}_{\mathit{F}}\right)}_{\mathit{t}}$	${\left({\mathit{M}}_{\mathit{T}}\right)}_{\mathit{t}}$	${\left({\mathit{M}}_{\mathit{D}}\right)}_{\mathit{t}}$	${\left({\mathit{M}}_{\mathit{S}}\right)}_{\mathit{t}}$	${\mathit{Be}}_{\mathit{avg}}$
	${10}^4$	0.0	0.000024	4.639946	3.769352	8.409322	0.999997
		0.04	0.000017	5.638370	3.157865	8.796251	0.999998
$0.0$	${10}^5$	0.0	0.000068	6.265106	6.671517	12.936691	0.999995
		0.04	0.000059	6.659704	5.614256	12.274019	0.999995
	${10}^6$	0.0	0.000131	8.444187	10.408137	18.852455	0.999993
		0.04	0.000114	8.048875	8.899719	16.948708	0.999993
	${10}^4$	0.0	0.000024	4.548592	3.797055	8.345671	0.999997
		0.04	0.000017	5.585220	3.170225	8.755461	0.999998
$0.1$	${10}^5$	0.0	0.000070	6.067520	6.733653	12.801243	0.999995
		0.04	0.000061	6.507360	5.601581	12.109002	0.999995
	${10}^6$	0.0	0.000136	8.102488	10.455931	18.558554	0.999993
		0.04	0.000117	7.812779	8.935463	16.748359	0.999993
	${10}^4$	0.0	0.000024	4.514655	3.811956	8.326635	0.999997
		0.04	0.000017	5.560058	3.182364	8.742438	0.999998
$0.2$	${10}^5$	0.0	0.000073	6.024801	6.794440	12.819314	0.999994
		0.04	0.000062	6.447220	5.640145	12.087426	0.999995
	${10}^6$	0.0	0.000141	8.059910	10.534990	18.595040	0.999992
		0.04	0.000119	7.715879	8.984502	16.700500	0.999993

Table 9

Variance analysis for $\overline{Nu}$.

Source	Totality ofSquares	Degree ofFreedom	MeanSquare	F-Value	p-Value
Model	83.28	20	4.16	261.81	<$0.0001$	Significant
$Ra$	42.77	1	42.77	2689.02	<$0.0001$
$D_f$	2.550	1	2.55	160.34	<$0.0001$
$Br$	4.380	1	4.38	275.46	<$0.0001$
$Sr$	0.0002	1	0.0002	0.0134	0.9102
$\varphi$	0.0066	1	0.0066	0.4141	0.5344
$Ra\times D_f$	0.2062	1	0.2062	12.97	0.0048
$Ra\times Br$	3.090	1	3.09	194.17	<$0.0001$
$Ra\times Sr$	0.1365	1	0.1365	8.58	0.0151
$Ra\times \varphi$	0.5344	1	0.5344	33.6	0.0002
$D_f\times Br$	0.0282	1	0.0282	1.77	0.2125
$D_f\times Sr$	0.2643	1	0.2643	16.62	0.0022
$D_f\times \varphi$	0.1650	1	0.165	10.38	0.0092
$Br\times Sr$	0.1207	1	0.1207	7.59	0.0203
$Br\times \varphi$	0.0731	1	0.0731	4.6	0.0577
$Sr\times \varphi$	0.0042	1	0.0042	0.2645	0.6183
$R{a}^2$	0.1342	1	0.1342	8.44	0.0157
${D_f}^2$	0.2387	1	0.2387	15.01	0.0031
$B{r}^2$	0.2370	1	0.237	14.9	0.0032
$S{r}^2$	0.0302	1	0.0302	1.9	0.1981
${\varphi }^2$	0.1164	1	0.1164	7.32	0.0221
Residual	0.1590	10	0.0159
Lack of Fit	0.1582	5	0.0316	180.98	<$0.0001$	Significant
Pure Error	0.0009	5	0.0002
Cor Total	83.440	30

Table 10

Sensitivity analysis for $\overline{Nu}$.

$\mathit{Ra}$	${\mathit{D}}_{\mathit{f}}$	$\mathit{Br}$	$\mathit{Sr}$	$\mathit{\varphi }$	$\frac{\mathit{\partial }\overline{\mathit{Nu}}}{\mathit{\partial }\mathit{Ra}}$	$\frac{\mathit{\partial }\overline{\mathit{Nu}}}{\mathit{\partial }{\mathit{D}}_{\mathit{f}}}$	$\frac{\mathit{\partial }\overline{\mathit{Nu}}}{\mathit{\partial }\mathit{Br}}$	$\frac{\mathit{\partial }\overline{\mathit{Nu}}}{\mathit{\partial }\mathit{Sr}}$	$\frac{\mathit{\partial }\overline{\mathit{Nu}}}{\mathit{\partial }\mathit{\varphi }}$
1	−1	−1	0	−1	$1.72\times {10}^{-5}$	17.83454	0.6553	20.7438	464.2789
				0	$2.20\times {10}^{-6}$	−43.43214	0.6553	20.7438	62.26667
				1	−$1.28\times {10}^{-5}$	−104.6988	0.6553	20.7438	−339.7457
1	−1	−1	1	−1	$1.51\times {10}^{-5}$	−3.92682	1.6728	20.7438	464.2789
				0	$8.60\times {10}^{-8}$	−65.19350	1.6728	20.7438	62.26667
				1	−$1.49\times {10}^{-5}$	−126.4602	1.6728	20.7438	−339.7457
1	−1	0	0	−1	$1.81\times {10}^{-5}$	17.83454	0.0445	21.7614	464.2789
				0	$3.15\times {10}^{-6}$	−43.43214	0.0445	21.7614	62.26667
				1	−$1.19\times {10}^{-5}$	−104.6988	0.0445	21.7614	−339.7457
1	−1	0	1	−1	$1.60\times {10}^{-5}$	−3.92682	1.0620	21.7614	464.2789
				0	$1.02\times {10}^{-6}$	−65.19350	1.0620	21.7614	62.26667
				1	−$1.40\times {10}^{-5}$	−126.4602	1.0620	21.7614	−339.7457

References

1. Xuan, Y.; Li, Q. Heat transfer enhancement of nanofluids. Int. J. Heat Fluid Flow; 2000; 21, pp. 58-64. [DOI: https://dx.doi.org/10.1016/S0142-727X(99)00067-3]

2. Anee, M.J.; Siddiqa, S.; Hasan, M.F.; Molla, M.M. Lattice Boltzmann simulation of natural convection of ethylene glycol-alumina nanofluid in a C-shaped enclosure with MFD viscosity through a parallel computing platform and quantitative parametric assessment. Phys. Scr.; 2023; 98, 095203. [DOI: https://dx.doi.org/10.1088/1402-4896/ace704]

3. Wang, X.; Song, Y.; Li, C.; Zhang, Y.; Ali, H.M.; Sharma, S.; Li, R.; Yang, M.; Gao, T.; Liu, M. et al. Nanofluids application in machining: A comprehensive review. Int. J. Adv. Manuf. Technol.; 2023; [DOI: https://dx.doi.org/10.1007/s00170-022-10767-2]

4. Can, A.; Selimefendigil, F.; Öztop, H.F. A review on soft computing and nanofluid applications for battery thermal management. J. Energy Storage; 2022; 53, 105214. [DOI: https://dx.doi.org/10.1016/j.est.2022.105214]

5. Kulkarni, S.; Chavali, S.; Dikshit, S. A review on analysis of Vapour Compression Refrigeration System (VCRS) for its performance using different ecofriendly refrigerants and nanofluids. Mater. Today Proceed.; 2023; 72, pp. 878-883. [DOI: https://dx.doi.org/10.1016/j.matpr.2022.09.085]

6. Rahmani, M.; Nejad, A.S.; Barzoki, M.F.; Kasaeian, A.; Sameti, M. Simulation of solar absorption refrigeration cycle with CuO nanofluid for summer cooling of a residential building. Therm. Sci. Eng. Prog.; 2022; 34, 101419. [DOI: https://dx.doi.org/10.1016/j.tsep.2022.101419]

7. Gürdal, M.; Arslan, K.; Gedik, E.; Minea, A.A. Effects of using nanofluid, applying a magnetic field, and placing turbulators in channels on the convective heat transfer: A comprehensive review. Renew. Sustain. Energy Rev.; 2022; 162, 112453. [DOI: https://dx.doi.org/10.1016/j.rser.2022.112453]

8. Mehta, B.; Subhedar, D.; Panchal, H.; Said, Z. Synthesis, stability, thermophysical properties and heat transfer applications of nanofluid–A review. J. Mol. Liq.; 2022; 364, 120034. [DOI: https://dx.doi.org/10.1016/j.molliq.2022.120034]

9. Malika, M.; Ashokkumar, M.; Sonawane, S.S. Mathematical and numerical investigations of nanofluid applications in the industrial heat exchangers. Applications of Nanofluids in Chemical and Bio-medical Process Industry; Elsevier: Amsterdam, The Netherlands, 2022; pp. 53-78.

10. Liu, Q.; Feng, X.B.; Xu, X.T.; He, Y.L. Multiple-relaxation-time lattice Boltzmann model for double-diffusive convection with Dufour and Soret effects. Int. J. Heat Mass Transf.; 2019; 139, pp. 713-719. [DOI: https://dx.doi.org/10.1016/j.ijheatmasstransfer.2019.05.026]

11. Xu, H.; Luo, Z.; Lou, Q.; Zhang, S.; Wang, J. Lattice Boltzmann simulations of the double-diffusive natural convection and oscillation characteristics in an enclosure with Soret and Dufour effects. Int. J. Therm. Sci.; 2019; 136, pp. 159-171. [DOI: https://dx.doi.org/10.1016/j.ijthermalsci.2018.10.015]

12. Magyari, E.; Postelnicu, A. Double-diffusive natural convection flows with thermosolutal symmetry in porous media in the presence of the Soret–Dufour effects. Transp. Porous Media; 2011; 88, pp. 149-167. [DOI: https://dx.doi.org/10.1007/s11242-011-9731-z]

13. Ghachem, K.; Kolsi, L.; Khan, S.U.; Abbas, T.; Maatki, C.; Saeed, M. Soret and Dufour aspect of viscoelastic fluid due to moving cylinder with viscous dissipation and convective boundary conditions. J. Indian Chem. Soc.; 2023; 100, 100913. [DOI: https://dx.doi.org/10.1016/j.jics.2023.100913]

14. Islam, M.M.; Molla, M.M.; Siddiqa, S.; Sheremet, M.A. MRT-LB simulation and response surface analysis of natural convection of non-Newtonian ferrofluid in an enclosure with non-uniformly heated radiator through GPU computing. Eng. Anal. Bound. Elem.; 2023; 157, pp. 92-118. [DOI: https://dx.doi.org/10.1016/j.enganabound.2023.09.003]

15. Zhang, J.; Gao, S.; Xiong, W.; Ye, L.; Liu, H.; Zhu, W.; Mu, Y.; Niu, W. Physical and Numerical Simulation of Tight Gas Flow at the Microscale. Energies; 2023; 16, 5937. [DOI: https://dx.doi.org/10.3390/en16165937]

16. Wanner, J.; Birke, K.P. Investigation of the Influence of Electrode Surface Structures on Wettability after Electrolyte Filling Based on Experiments and a Lattice Boltzmann Simulation. Energies; 2023; 16, 5640. [DOI: https://dx.doi.org/10.3390/en16155640]

17. Islam, M.M.; Hasan, M.F.; Molla, M.M. Analysis of Heat Transfer Characteristics of MHD Ferrofluid by the Implicit Finite Difference Method at Temperature-Dependent Viscosity Along a Vertical Thin Cylinder. Iran. J. Sci. Technol. Trans. Mech. Eng.; 2023; [DOI: https://dx.doi.org/10.1007/s40997-023-00656-8]

18. Amiri, I.; Rashed, A.N.Z. Power enhancement of the U-shape cavity microring resonator through gap and material characterizations. J. Opt. Commun.; 2019; pp. 1-5. [DOI: https://dx.doi.org/10.1515/joc-2019-0108]

19. Wu, H.; Kulatilake, P.H.; Zhao, G.; Liang, W.; Wang, E. A comprehensive study of fracture evolution of brittle rock containing an inverted U-shaped cavity under uniaxial compression. Comput. Geotech.; 2019; 116, 103219. [DOI: https://dx.doi.org/10.1016/j.compgeo.2019.103219]

20. Ali, F.H.; Hamzah, H.K.; Egab, K.; Arıcı, M.; Shahsavar, A. Non-Newtonian nanofluid natural convection in a U-shaped cavity under magnetic field. Int. J. Mech. Sci.; 2020; 186, 105887. [DOI: https://dx.doi.org/10.1016/j.ijmecsci.2020.105887]

21. Ma, Y.; Mohebbi, R.; Rashidi, M.; Yang, Z. Simulation of nanofluid natural convection in a U-shaped cavity equipped by a heating obstacle: Effect of cavity’s aspect ratio. J. Taiwan Inst. Chem. Eng.; 2018; 93, pp. 263-276. [DOI: https://dx.doi.org/10.1016/j.jtice.2018.07.026]

22. Nabwey, H.A.; Rashad, A.; Khan, W.A.; Alshber, S.I. Effectiveness of magnetize flow on nanofluid via unsteady natural convection inside an inclined U-shaped cavity with discrete heating. Alex. Eng. J.; 2022; 61, pp. 8653-8666. [DOI: https://dx.doi.org/10.1016/j.aej.2022.02.010]

23. Esfe, M.H.; Rostamian, H.; Toghraie, D.; Hekmatifar, M.; Abad, A.T.K. Numerical study of heat transfer of U-shaped enclosure containing nanofluids in a porous medium using two-phase mixture method. Case Stud. Therm. Eng.; 2022; 38, 102150. [DOI: https://dx.doi.org/10.1016/j.csite.2022.102150]

24. Sajjadi, H.; Delouei, A.A.; Izadi, M.; Mohebbi, R. Investigation of MHD natural convection in a porous media by double MRT lattice Boltzmann method utilizing MWCNT–Fe₃O₄/water hybrid nanofluid. Int. J. Heat Mass Transf.; 2019; 132, pp. 1087-1104. [DOI: https://dx.doi.org/10.1016/j.ijheatmasstransfer.2018.12.060]

25. Saha, B.C.; Mahapatra, T.; Pal, D. Heat and mass flow in lid-driven magnetohydrodynamic double diffusive mixed convection of nanofluid within a trapezoidal enclosure. J. Nanofluids; 2019; 8, pp. 817-829. [DOI: https://dx.doi.org/10.1166/jon.2019.1623]

26. Nithyadevi, N.; Yang, R.J. Double diffusive natural convection in a partially heated enclosure with Soret and Dufour effects. Int. J. Heat Fluid Flow; 2009; 30, pp. 902-910. [DOI: https://dx.doi.org/10.1016/j.ijheatfluidflow.2009.04.001]

27. He, B.; Lu, S.; Gao, D.; Chen, W.; Li, X. Lattice Boltzmann simulation of double diffusive natural convection of nanofluids in an enclosure with heat conducting partitions and sinusoidal boundary conditions. Int. J. Mech. Sci.; 2019; 161, 105003. [DOI: https://dx.doi.org/10.1016/j.ijmecsci.2019.07.002]

28. El-Shorbagy, M.; Eslami, F.; Ibrahim, M.; Barnoon, P.; Xia, W.F.; Toghraie, D. Numerical investigation of mixed convection of nanofluid flow in a trapezoidal channel with different aspect ratios in the presence of porous medium. Case Stud. Therm. Eng.; 2021; 25, 100977. [DOI: https://dx.doi.org/10.1016/j.csite.2021.100977]

29. Barnoon, P.; Toghraie, D.; Rostami, S. Optimization of heating-cooling generators with porous components/cryogenic conductors on natural convection in a porous enclosure: Using different two-phase models and single-phase model and using different designs. Int. Commun. Heat Mass Transf.; 2020; 111, 104472. [DOI: https://dx.doi.org/10.1016/j.icheatmasstransfer.2019.104472]

30. Bondareva, N.S.; Sheremet, M.A.; Oztop, H.F.; Abu-Hamdeh, N. Heatline visualization of natural convection in a thick walled open cavity filled with a nanofluid. Int. J. Heat Mass Transf.; 2017; 109, pp. 175-186. [DOI: https://dx.doi.org/10.1016/j.ijheatmasstransfer.2017.01.124]

31. Chamkha, A.J.; Miroshnichenko, I.V.; Sheremet, M.A. Numerical analysis of unsteady conjugate natural convection of hybrid water-based nanofluid in a semicircular cavity. J. Therm. Sci. Eng. Appl.; 2017; 9, 041004. [DOI: https://dx.doi.org/10.1115/1.4036203]

32. Bondarenko, D.S.; Sheremet, M.A.; Oztop, H.F.; Ali, M.E. Natural convection of Al₂O₃/H₂O nanofluid in a cavity with a heat-generating element. Heatline visualization. Int. J. Heat Mass Transf.; 2019; 130, pp. 564-574. [DOI: https://dx.doi.org/10.1016/j.ijheatmasstransfer.2018.10.091]

33. Khanafer, K.; Vafai, K.; Lightstone, M. Buoyancy-driven heat transfer enhancement in a two-dimensional enclosure utilizing nanofluids. Int. J. Heat Mass Transf.; 2003; 46, pp. 3639-3653. [DOI: https://dx.doi.org/10.1016/S0017-9310(03)00156-X]

34. Cianfrini, M.; Corcione, M.; Quintino, A. Natural convection in square enclosures differentially heated at sides using alumina-water nanofluids with temperature-dependent physical properties. Therm. Sci.; 2015; 19, pp. 591-608. [DOI: https://dx.doi.org/10.2298/TSCI120328111C]

35. Corcione, M.; Cianfrini, M.; Quintino, A. Two-phase mixture modeling of natural convection of nanofluids with temperature-dependent properties. Int. J. Therm. Sci.; 2013; 71, pp. 182-195. [DOI: https://dx.doi.org/10.1016/j.ijthermalsci.2013.04.005]

36. Chen, X.; Yu, P.; Winoto, S.; Low, H. Free convection in a porous wavy cavity based on the Darcy-Brinkman-Forchheimer extended model. Numer. Heat Transf. Part A Appl.; 2007; 52, pp. 377-397. [DOI: https://dx.doi.org/10.1080/10407780701301595]

37. Ren, Q.; Chan, C.L. Numerical study of double-diffusive convection in a vertical cavity with Soret and Dufour effects by lattice Boltzmann method on GPU. Int. J. Heat Mass Transf.; 2016; 93, pp. 538-553. [DOI: https://dx.doi.org/10.1016/j.ijheatmasstransfer.2015.10.031]

38. Rahman, A.; Nag, P.; Molla, M.M. Non-Newtonian effects on MHD thermosolutal free convection and entropy production of nanofluids in a rectangular enclosure using the GPU-based mesoscopic simulation. Waves Random Complex Media; 2022; [DOI: https://dx.doi.org/10.1080/17455030.2022.2119303]

39. Liu, Q.; He, Y.L.; Li, Q.; Tao, W.Q. A multiple-relaxation-time lattice Boltzmann model for convection heat transfer in porous media. Int. J. Heat Mass Transf.; 2014; 73, pp. 761-775. [DOI: https://dx.doi.org/10.1016/j.ijheatmasstransfer.2014.02.047]

40. Trouette, B. Lattice Boltzmann simulations of a time-dependent natural convection problem. Comput. Math. Appl.; 2013; 66, pp. 1360-1371. [DOI: https://dx.doi.org/10.1016/j.camwa.2013.07.024]

41. Molla, M.M.; Haque, M.J.; Khan, M.A.I.; Saha, S.C. Gpu accelerated multiple-relaxation-time lattice boltzmann simulation of convective flows in a porous media. Front. Mech. Eng.; 2018; 4, 15. [DOI: https://dx.doi.org/10.3389/fmech.2018.00015]

42. Wang, J.; Wang, D.; Lallemand, P.; Luo, L.S. Lattice Boltzmann simulations of thermal convective flows in two dimensions. Comput. Math. Appl.; 2013; 65, pp. 262-286. [DOI: https://dx.doi.org/10.1016/j.camwa.2012.07.001]

43. Mezrhab, A.; Moussaoui, M.A.; Jami, M.; Naji, H.; Bouzidi, M. Double MRT thermal lattice Boltzmann method for simulating convective flows. Phys. Lett. A; 2010; 374, pp. 3499-3507. [DOI: https://dx.doi.org/10.1016/j.physleta.2010.06.059]

44. Ilis, G.G.; Mobedi, M.; Sunden, B. Effect of aspect ratio on entropy generation in a rectangular cavity with differentially heated vertical walls. Int. Commun. Heat Mass Transf.; 2008; 35, pp. 696-703. [DOI: https://dx.doi.org/10.1016/j.icheatmasstransfer.2008.02.002]

45. Magherbi, M.; Abbassi, H.; Brahim, A.B. Entropy generation at the onset of natural convection. Int. J. Heat Mass Transf.; 2003; 46, pp. 3441-3450. [DOI: https://dx.doi.org/10.1016/S0017-9310(03)00133-9]

46. Nithiarasu, P.; Seetharamu, K.N.; Sundararajan, T. Natural convective heat transfer in a fluid saturated variable porosity medium. Int. J. Heat Mass Transf.; 1997; 40, pp. 3955-3967. [DOI: https://dx.doi.org/10.1016/S0017-9310(97)00008-2]

47. Guo, Z.; Zhao, T.S. A lattice Boltzmann model for convection heat transfer in porous media. Num. Heat. Transf. Part B; 2005; 47, pp. 157-177. [DOI: https://dx.doi.org/10.1080/10407790590883405]

48. Ali, L.; Kumar, P.; Iqbal, Z.; Alhazmi, S.E.; Areekara, S.; Alqarni, M.; Mathew, A.; Apsari, R. The optimization of heat transfer in thermally convective micropolar-based nanofluid flow by the influence of nanoparticle’s diameter and nanolayer via stretching sheet: Sensitivity analysis approach. J. Non-Equilib. Thermodyn.; 2023; [DOI: https://dx.doi.org/10.1515/jnet-2022-0064]