Comparison of average Nusselt number between generalized porous medium models for Nithiarasu et al. [31], FVM and ISPH method in the non-Darcy flow regime,

## Abstract

A study on heat and mass transfer behaviour on an anisotropic porous medium embedded in square cavity/annulus is conducted using Incompressibe Smoothed Particle Hydrodynamics (ISPH) method. In the case of square cavity, the left wall has hot temperature Th and mass Ch and the right wall has cool temperature Tc and mass Cc and both of the top and bottom walls are adiabatic. While in the case of square annulus, the inner surface wall is considered to have a cool temperature Tc and mass Cc while the outer surface is exposed to a hot temperature Th and mass Ch. The governing partial differential equations are transformed to non-dimensional governing equations and are solved using ISPH method. The results present the influences of the Dufour and Soret effects on the heat and mass transfer. The effects of various physical parameters such as Darcy parameter, permeability ratio, inclination angle of permeability and Rayleigh numbers on the temperature and concentration profiles together with the local Nusselt and Sherwood numbers are presented graphically. The results from the current ISPH method are well validated and have favorable comparisons with previously published results and solutions by the finite volume method.

### Keywords

- Anisotropic porous media
- Double-diffusive
- Dufour number
- ISPH
- Natural convection
- Non-Darcy flow
- Soret number
- Square annulus

## 1. Introduction

The enormous interest in the double diffusive convection in the recent years has led researchers to an extensive study on this topic due to its applicability in the industry as well as in engineering fields. The various aspects related to the heat and mass transfer have been addressed in the extensive literature [1-12]. The thermal diffusion (Soret) and diffusion-thermo (Dufour) effects in the combined heat and mass transfer affect the flow field in free convection significantly [13]. Thus both Soret and Dufour effects have considerable implications on the flow field. Aly and his co-authors [14-17] introduced several studies related to the effects of thermal diffusion (Soret) and diffusion-thermo (Dufour) in the combined heat and mass transfer. Mansour et al. [14] investigated the effects of chemical reaction, thermal stratification, Soret number and Dufour number on MHD free convective heat and mass transfer of a viscous, incompressible and electrically conducting fluid on a vertical stretching surface embedded in a saturated porous medium. Chamkha et al. [15] studied unsteady double-diffusive natural convective MHD flow along a vertical cylinder in the presence of chemical reaction, thermal radiation and Soret and Dufour effects. Chamkha and Aly [16] solved numerically the steady boundary-layer stagnation-point flow of a polar fluid towards a stretching surface embedded in porous media in the presence of the effects of Soret and Dufour numbers and first-order homogeneous chemical reaction. Aly et al. [17] studied the effects of soret and Dufour numbers on free convection over two different types of flows, namely isothermal and adiabatic stretching surfaces embedded in porous media in the presence of a homogeneous first-order chemical reaction. A study made by Nithyadevi and Yang [18] and Weaver and Viskanta [19] discovered that the Soret and Dufour effects will become very significant when the temperature and concentration gradients are high. The fluid velocity and heat transfer rates increase when Dufour parameter is very high. For Soret effect, as the parameter is high, the mass transfer increases but the velocity decreases.

On the other hand, the natural convection in an anisotropic porous medium is an important area of research due to its wide range of applications including thermal insulation, flow in mushy region of a solidifying alloy [20] and flow past heat exchanger tubes [21]. The non-Darcy effects on natural convection in porous media have also received significant attention as a result of the experiments conducted with several combinations of solids and fluids. These experiments covered wide ranges of governing parameters that indicate that the experimental data for systems other than glass water at low Rayleigh numbers do not agree with the theoretical predictions based on the Darcy flow model. This divergence in the heat transfer results has been reviewed in detail in Cheng [22] and Prasad et al. [23], among others. Thus, extensive efforts are being made to include the inertia and viscous diffusion terms in the flow equations and to examine their effects in order to develop a reasonably accurate mathematical model for convective transport in porous media. Detailed accounts of the research into non-Darcy convection have been reported in Tien and Hong [24], Cheng [25], Prasad et al. [26], and Kiadias and Prasad [27]. Nield and Bejan [28] provided an excellent summary of the subject regarding porous media models. The numerical studies of the natural convection flow in anisotropic porous media were conducted by use of Brinkman equation [29] or Brinkman-Forchheimer equation with permeability tensor [30]. They demonstrated that their formulations were accurate in predicting the flow and heat transfer for various inclinations of the principal permeability direction, permeability ratios, and Darcy numbers. The natural convective flow and heat transfer in a fluid saturated anisotropic porous medium have been investigated using the generalized non-Darcy models Nithiarasu et al. [31].

In recent years, the SPH method had been applied into compressible and incompressible viscous fluid flow problems [32, 33]. The SPH was originally developed in compressible flow, and then some special treatment was required to satisfy the incompressible condition. A proposal for developing an incompressible SPH (ISPH) model has been introduced, which pressure is implicitly calculated by solving a discretized pressure Poisson equation at every time step [34-38]. Cummins and Rudman [34] introduced a new formulation for enforcing incompressibility in Smoothed Particle Hydrodynamics (SPH). The method uses a fractional step with the velocity field integrated forward in time without enforcing incompressibility. The resulting intermediate velocity field is then projected onto a divergence-free space by solving a pressure Poisson equation derived from an approximate pressure projection. Asai et al. [35] introduced the stabilized incompressible SPH method by relaxing the density invariance condition. Aly et al. [36-38] applied the stabilized incompressible SPH method to simulate multi-fluid problems, fluid-structure interaction and fluid-soil-structure interactions. Aly et al. [36] modeled the surface tension force for free surface flows and an eddy viscosity based on the Smagorinsky sub-grid scale model using ISPH (Incompressible smoothed particle hydrodynamics) method. They declared that, the eddy viscosity has clear effects in adjusting the splashes and reduces the deformation of free surface in the interaction between the two fluids. In addition, Aly et al. [37] applied the stabilized incompressible SPH method to simulate free falling of rigid body and water entry/exit of circular cylinder into water tank. Aly [38] discussed in details the simulation of fluid-soil-structure interactions using an improved ISPH method.

Numerical modeling of transient natural convection by using SPH method has also been investigated. Chaniotis et al. [39] proposed a remeshing algorithm based on weakly compressible flow approach and performed a comprehensive study for non-isothermal flows. SPH simulation of flow and energy transport using SPH was performed by Szewc et al. [40]. In their study, natural convection in a square cavity problem with a Boussinesq and a non-Boussinesq formulation was performed. They introduced a new variant of the Smoothed Particle Hydrodynamics (SPH) simulations of the natural convection phenomena. Danis et al. [41] modeled the transient and laminar natural convection in a square cavity using SPH method with a discretization tool on uniform Eulerian grids. Recently, Aly [42] modelled the multi-phase flow and natural convection in a square/cubic cavity using ISPH method in two and three dimensions. Rayleigh-Taylor instability between two and three adjacent fluid layers has been simulated and also the natural convection in a square/cubic cavity has been introduced with a good agreement compared to benchmark tests. Aly and Asai [43] modelled non-Darcy flows through porous media using an extended ISPH method. In their study, unsteady lid-Driven flow, natural convection in non-Darcy porous cavities and natural convection in porous medium-fluid interface are examined separately by using ISPH method. In addition, Aly and Sameh [44] modelled the non-Darcy flows through anisotropic porous media for natural/mixed convection and heat transfer in a cavity using ISPH method.

In this study, we presented a generalized porous medium model based on the ISPH method for natural convection and heat and mass transfer in cavity saturated with anisotropic porous media under the effects of Soret and Dufour numbers. Here, we described the implementation of the projection method procedure for a more general hydrodynamically, thermally and diffusion anisotropic porous medium. A semi-implicit time integration scheme is applied for double-diffusive natural convection with anisotropic porous media under the effects of Soret and Dufour numbers. Here, two different cases of boundary condition for the square cavity have been studied numerically using ISPH method. The first case considers natural convection in a square cavity, in which the left wall has hot temperature

## 2. Problem description

Fig. 1 shows the physical model for the current problem for the natural convection in porous enclosure with two different cases. In the first case, it is assumed that, the left wall has hot temperature

## 3. Mathematical analysis

The continuity, momentum, energy and concentration equations for the unsteady natural convection in the two-dimensional enclosure can be written in dimensional form as follows:

where,

The second order tensor of permeability,

Substituting

Substituting

The second order of the thermal diffusivity,

Substituting

In similar manner, the second order tensor of mass diffusivity,

Substituting

The following dimensionless parameters are used to convert equations (1)-(4) to non-dimensional forms in natural convection case:

Then, the dimensionless continuity, momentum, energy and mass equations are written as follows:

where,

### 3.1. Boundary treatment

The dimensionless boundary conditions are:

First case (square cavity filled with anisotropic porous media):

Second case (square annulus filled with anisotropic porous media)

The local Nusselt number and local Sherwood number can be defined for the hot wall of square cavity as follows:

For the majority of design problems, the knowledge of the average Nusselt number is very useful. The average Nusselt number and average Sherwood number are obtained through the integration:

The local Nusselt number and local Sherwood number can be defined for the square annulus as follows:

Horizontal hot walls:

Vertical hot walls:

For the square annulus case, the average Nusselt number and average Sherwood number are obtained through the following integration:

Average Nusselt and Sherwood numbers at the left hot surface:

Average Nusselt and Sherwood numbers at the right hot surface:

Average Nusselt and Sherwood numbers at the top hot surface:

Average Nusselt and Sherwood numbers at the bottom hot surface:

Then, the total average Nusselt and Sherwood numbers are defined as:

## 4. Numerical method

The dimensionless governing equations and boundary conditions were solved using both of finite volume method and ISPH method. Here, the finite volume method is introduced for the validation tests.

### 4.1. Finite volume method

The finite control volume method used here, was described in Patankar [45]. The combined continuity, momentum, energy and concentration equations have been solved numerically using the SIMPLE algorithm [45]. It is worth mentioning that, collocated, regular and orthogonal grids were used in this implementation. Also, Rhie-Chow interpolation was used. The algebraic equations resulting from this treatment are solved using alternating direct implicit (ADI) procedure. It is found that, the 121×121 grid is sufficiently enough to solve the system of equations. The unknowns dependent variables were calculated iteratively until the following criteria of convergence was fulfilled:

### 4.2. ISPH method

The ISPH algorithm is implemented in a semi-implicit form in order to solve the incompressible viscous flow equations. In this section, the procedure for the solution of porous medium equations is described.

The ISPH method is based on the calculation of an intermediate velocity from a momentum equation where the pressure gradients are omitted. Then, the pressure is evaluated through solving the pressure Poisson equation (PPE). The PPE after SPH interpolation is solved by a preconditioned diagonal scaling Conjugate Gradient PCG method [46] with a convergence tolerance (=1.0x10^{-9}). Finally, the velocity is corrected using the evaluated pressure.

#### 4.2.1. Temporal discretization

The momentum equation can be discretized in time using predictor-corrector scheme. Here, the momentum, energy and concentration equations are described in lagrangian description. In particular, the time discrete momentum equation in its semi-implicit form can be written as follows:

For simplicity,

Then, equations (40) and (41) can be defined as follows:

** Step 1.**Predict the velocity

The first step of the predictor-corrector scheme is the calculation of an intermediate velocity

** Step 2.**Solving the pressure Poisson equation

In the second step, the pressure is calculated using the modified Poisson equation, which ensures that the continuity equation is satisfied, and for the generalized model can be written as follows:

Note, the relaxation coefficient,

** Step 3.**Corrected velocity

In the third step, the real velocity values are obtained using the following correction:

** Step 4.**Thermal flow problems:

In this step, the time discretization of the energy equation is introduced:

** Step 5.**Concentration flow problems:

In this step, the time discretization of the concentration equation is introduced:

#### 4.2.2. SPH Formulations

A spatial discretization using scattered particles, which is based on the SPH, is summarized. First, a physical scalar function

where, * W*is a weight function called by smoothing kernel function in the SPH literature. In the smoothing kernel function,

*are the distance between neighbor particles and smoothing length respectively. For SPH numerical analysis, the integral Eq. (52) is approximated by a summation of contributions from neighbor particles in the support domain.*h

where, the subscripts * i*and

*indicate positions of labeled particle,*j

*, respectively. Note that, the triangle bracket*j

Also, the other expression for the gradient can be represented by:

In this study, quintic spline function is utilized as a kernel function

where

#### 4.2.3. Projection-based ISPH formulations

Here, the projection method for incompressible fluid problem, which is summarized in section 4, is discretized into particle quantities based on the SPH methodology. For this purpose, the gradient of pressure and the divergence of velocity are approximated as follows:

Although Laplacian could be derived directly from the original SPH approximation of a function in Eq. (58), this approach may lead to a loss of resolution. Then, the second derivative of velocity for viscous force and the Laplacian of pressure have been proposed by Morris et al. [33] by approximation expression as follows:

where,

In this study, Laplacian of velocity is given by:

Similarly, the Laplacian of pressure in pressure Poisson equation (PPE) is given by:

The PPE after SPH interpolation is solved by a preconditioned (diagonal scaling) Conjugate Gradient (PCG) method [46] with a convergence tolerance (=1.0×10^{-9}).

The Laplacian operators for both of the temperature and concentration are given as:

## 5. Numerical validation tests

To check the accuracy of the numerical method employed for the solution of the problem under consideration, it was validated (after making the necessary modifications) with the problem of Nithiarasu et al. [31]. Comparison of average Nusselt number between Brinkman and Generalized porous medium models for Nithiarasu et al. [31], FVM and ISPH methods in the non-Darcy flow regime ^{4} with Darcy parameter,

^{*} | |||||

0° | 0.1 | 10^{4} | 1.587 | 1.570 | 1.517 |

45° | 0.1 | 10^{4} | 1.573 | 1.563 | 1.522 |

90° | 0.1 | 10^{4} | 1.579 | 1.573 | 1.524 |

0° | 10 | 10^{4} | 1.106 | 1.110 | 1.132 |

45° | 10 | 10^{4} | 1.119 | 1.130 | 1.139 |

90° | 10 | 10^{4} | 1.106 | 1.111 | 1.134 |

Comparison of isotherms lines between FVM and ISPH methods under the effect of Rayleigh number and permeability ratio at Darcy number

## 6. Results and discussion

In this section, two different cases of boundary conditions for unsteady natural convection flow in square cavity saturated with anisotropic porous media under the effects of Soret and Dufour numbers have been discussed in details using ISPH method.

### 6.1. Natural convection in square cavity saturated with anisotropic porous media

In this study, the unsteady natural convection in square cavity saturated with anisotropic porous media under the effects of Soret and Dufour numbers has been solved numerically using ISPH method. It is assumed that, the vertical walls are kept at constant temperatures,

The effects of combined Soret and Dufour numbers on the isothermal lines, concentration lines and horizontal and vertical velocities contours at inclination angles

It is observed that, in general, the isotherms lines gather beside the vertical walls, indicating thermals boundary layers beside the bottom of the left wall and top of the right wall of the cavity. As the Soret number increases with decrease in Dufour number, the isothermal lines are slightly increase the thermals boundary layers beside the bottom of the left wall and top of the right wall of the cavity. Regarding the solutal concentration lines, these contours are parallel to each other within the core of the cavity. Moreover, as the Soret number increases with decrease in Dufour number leads to the solutal concentration lines are highly affected and their contours are paralleled within the core of the cavity. The contours of horizontal velocity component formed in the shape of two vertically-extended clockwise and anticlockwise circular cells next to the bottom and top walls of the enclosure. On the contrary, the contours of the vertical velocity component formed in the shape of two horizontally-extended clockwise and anticlockwise circular cells next to left and right walls of the enclosure. Also, the horizontal and vertical velocity components contours become slightly stretch in case of increasing Soret number with decreasing Dufour number.

Fig. 4 presents the effects of Soret and Dufour numbers on temperature, solutal concentration, horizontal velocity profiles and vertical velocity profiles. It is noted that, an increase in Soret number accompanied by a decrease in Dufour number leads to increase both of fluid concentration and horizontal and vertical velocity profiles and decrease in fluid temperature.

The effects of Darcy number from 10^{-2} to 10^{-4} at inclination angles * N*considered, the interaction between the thermal and compositional buoyancy effects is significant. It is observed that, the heat and mass conduction in cavity are increase clearly and the flow moves faster at the case of increasing Darcy number from10

^{-4}to 10

^{-2}. In figure 5, as the inverse Darcy number

^{2}to 10

^{4}, the temperature and concentration contours become more parallel to the vertical walls and the horizontal and vertical values are decrease. Two parallel recirculations for both of horizontal and vertical velocities are formed near to the boundaries. It can be concluded that the main contributions of the presence of the porous medium for a buoyancy ratio of unity are predicted to be a flow retardation effect and a suppression of the overall heat transfer in the enclosure.

Fig. 6 depicts the effects of Darcy parameter on temperature, solutal concentration, horizontal velocity profiles and vertical velocity profiles at two values of Rayleigh number Ra=10^{3} and 10^{4}. At the presence of the porous medium and Darcy number changes from 10^{-2} to 10^{-4}, both of the horizontal and vertical velocity profiles are decrease. The temperature and concentration profiles are decrease linearly as

In order to demonstrate the variation in velocity, temperature and concentration fields under the effects of combined an inclination angle with permeability ratio parameter. The horizontal and vertical velocity, temperature and concentration profiles and contours through the cavity under the effects of combined an inclination angle with permeability ratio parameter have been shown in Figs. 7 and 8.

In these figures, the permeability ratio parameter showed clear effects in isothermal lines, and horizontal and vertical velocities contours and also in temperature and horizontal and vertical velocity profiles compared to the effects of an inclination angle. As the permeability ratio parameter increases from 0.1 to 10, the temperature profiles are increase until their profiles become linear, this indicates the limit of no-flow. Also, the net fluid velocity in the vicinity of the heated wall is decrease. The concentration profiles are affected slightly by changing both of an inclination angle with permeability ratio parameter.

The time histories for the effects of the Darcy number with two values of Rayleigh number

In addition, the average Nusselt number increases as the Soret number increases and Dufour number decreases. While, the average Sherwood number decreases as the Soret number increases and Dufour number decreases.

### 6.2. Natural convection in a square annulus filled with anisotropic porous media

In this section, we studied heat and mass behavior on a square annulus saturated with anisotropic porous media.

The effects of annulus length in the temperature, concentration, average Nusselt number and average Sherwood number have been introduced in Figs. 11 and 12. Fig. 11 Depicts the effects of annulus length on the isotherms lines, concentration lines and vertical and horizontal velocity contours, respectively at Darcy parameter

Fig. 13 Depicts the time histories for the effects of combined Soret number with Dufour number on the average Nusselt number and average Sherwood number, respectively at annulus length 0.5,

## 7. Conclusion

The unsteady Double-diffusive natural convection in an anisotropic porous square cavity/annulus has been investigated numerically using the stabilized ISPH method. In the current ISPH algorithm, a semi-implicit velocity correction procedure was successfully used. The ISPH solution is validated by direct comparisons with previously published work and FVM solution on special cases of the problem, and the results show good agreements with these references. Graphic results for the temperature, concentration, horizontal velocity and vertical velocity contours and representative velocity, temperature, and concentration profiles at the cavity midsection for various parametric conditions were presented and discussed.

The following findings are summarized from the present investigation:

The permeability ratio parameter increases leads to decrease in the both of heat conduction and flow regime.

For the square annulus, the bottom wall has higher average Nusselt and Sherwood numbers than the rest of the walls. The top wall has lower average Nusselt and Sherwood numbers than the rest of the walls.

When the annulus length equals 0.75, there is almost no fluid activity at the top and bottom walls and from our observation due to restriction in space and movement of the fluid, the isotherms and concentration lines are almost equal under the effect of any parameter.

As the Soret number increases with decreasing the value of Dufour number, the average Nusselt number increases. While, the average Sherwood number decreases as the Soret number increases with decreasing the Dufour number.

## Nomenclature

concentration of species | time | ||

dimensionless species concentration | temperature | ||

specific heat | dimensionless temperature | ||

Darcy parameter | dimensionless velocity components | ||

Dufour number | velocity vector | ||

annulus length | Cartesian coordinates | ||

particle size | dimensionless coordinates | ||

Forchheimer coefficient | |||

gravitational acceleration vector | thermal diffusivity | ||

permeability | thermal expansion coefficient | ||

thermal conductivity | compositional expansion coefficient | ||

Nusselt number | porosity | ||

Buoyancy ratio | viscosity | ||

Lewis number | kinematic viscosity | ||

pressure | ratio of heat capacities | ||

Prandtl number | density | ||

Rayleigh number | dimensionless time | ||

Sherwood number | Laplacian operator | ||

Soret number |

## Acknowledgments

This work was supported by JSPS KAKENHI Grant Numbers 26282106 and 13F03051.