Double-Diffusive Natural Convection with Cross-Diffusion Effects in an Anisotropic Porous Enclosure Using ISPH Method

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, inclina‐ tion angle of permeability and Rayleigh numbers on the temperature and concentra‐ tion 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.


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][2][3][4][5][6][7][8][9][10][11][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][15][16][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 formula-tions 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][35][36][37][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][37][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 T h and mass C h and the right wall has cool temperature T c and mass C c and both of the top and bottom walls are adiabatic. While the second case considers natural convection in square annulus, in which the inner surface wall is considered to have a cool temperature T c and mass C c and the outer surface is exposed to a hot temperature T h and mass C h . 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. 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 T h and mass C h and the right wall has cool temperature T c and mass C c . However, the horizontal walls are adiabatic. While for the second case, in square annulus the inner surface wall is considered to have a cool temperature T c and mass C c and the outer surface is exposed to a hot temperature T h and mass C h . In the both cases, the fluid and the porous medium are assumed to be in local thermodynamic equilibrium. The porous medium is anisotropic in permeability as well as in Forchheimer coefficient. The non-Darcy model is used to model this phenomenon. Also, the usual Boussinesq approximation is invoked to model the density variation in the buoyancy terms.

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:

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: Anisotropic Porous media Adiabatic Adiabatic Anisotropic Porous media where, V is the (Darcy) velocity vector, P is the fluid pressure, T ′ is the temperature, and ε is the porosity of the medium, σ is the ratio of heat capacities and C ′ is the concentration of species. D TC is Dufour diffusivity and D CT is Soret diffusivity.
The second order tensor of permeability, K, Forchheimer coefficient Ē, thermal diffusivity, k, and mass diffusivity, D , can be written in the following coordinate system as: cos sin sin cos , sin cos cos sin Substituting K * (Eq. (6)) in Eq. (5) and by applying inverse, the permeability tensor can be reduced to: cos sin 1 sin cos 1 , 1 sin cos cos sin Substituting E * . (Eq. (10)) in Eq. (7), the Forchheimer coefficient tensor can be reduced to: The second order of the thermal diffusivity, k can be written as: Substituting k * (Eq. (14)) in Eq. (13), the thermal conductivity tensor can be reduced to: In similar manner, the second order tensor of mass diffusivity, D , can be written as: Substituting D * (Eq. (17)) in Eq. (16), the mass diffusivity tensor can be reduced to: 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: 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: Double-diffusive natural convection with Cross-Diffusion Effects using ISPH method http://dx.doi.org/10.5772/60879 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:

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.

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:

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.

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 ( U * ) from the momentum equation without including the pressure terms. Thus, the following is obtained: 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, γ, (0 ≤ γ ≤ 1) can be decided from pre-analysis calculation as Asai et al. [35]. In this study, the particle size d 0 is taken as 2.0 cm and then the relaxation coefficient is decided as γ = 0.25.

Step 3.Corrected velocity
In the third step, the real velocity values are obtained using the following correction: Double-diffusive natural convection with Cross-Diffusion Effects using ISPH method http://dx.doi.org/10.5772/60879 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:

SPH Formulations
A spatial discretization using scattered particles, which is based on the SPH, is summarized. First, a physical scalar function ϕ(X i , τ) at a sampling point X i can be represented by the following integral form: where, W is a weight function called by smoothing kernel function in the SPH literature. In the smoothing kernel function, R ij = | X i − X j | and h 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.
( , ) ( , ) ( , ), where, the subscripts i and j indicate positions of labeled particle, ρ j and m j mean density and representative mass related to particle j, respectively. Note that, the triangle bracket ϕ i means SPH approximation of a function ϕ. The gradient of the scalar function can be assumed by using the above defined SPH approximation as follows: Also, the other expression for the gradient can be represented by: 2 2 ( ) ( , ).
In this study, quintic spline function is utilized as a kernel function ( ) where β d is 7 / 478πh 2 and 3 / 358πh 3 , respectively, in two-and three dimension space.

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: 2 2 ( ) ( , ), 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, η is a parameter to avoid a zero dominator, and its value is usually given by η 2 = 0.0001 h 2 . For the case of υ i = υ j and ρ i = ρ j , the Laplacian term is simplified as: . , 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:

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 Da = 10 −2 at buoyancy ratio N = 0, has been introduced in table 1. In this table, the Rayleigh numbers was taken as Ra=10 4 with Darcy parameter, Da = 10 −2 . The current FVM and ISPH results agree well with the data from Nithiarasu et al. [31]. Note, the current comparison was reported before by Aly and Sameh [44].  Comparison of isotherms lines between FVM and ISPH methods under the effect of Rayleigh number and permeability ratio at Darcy number Da = 10 −2 , porosity parameter ε = 0.5, buoyancy ratio N = 0 and an inclination angle of permeability θ = 0° is introduced in figure 2. In this figure, the isotherm lines are affected clearly by both of the Rayleigh number and permeability ratio. As the Rayleigh number increases, the heat conduction increases and the isotherms lines are strongly compressed in the bottom part of the left wall and in the top part of the right wall.
As the permeability ratio increases from 0.1 to 10, then the heat conduction decreases.

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.

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, T h and T c , where, T h is the hot wall and T c is cool wall and constant concentrations, C h and C c , where, C h > C c . However, the horizontal walls are adiabatic.
The effects of combined Soret and Dufour numbers on the isothermal lines, concentration lines and horizontal and vertical velocities contours at inclination angles θ = 0° and 45° at permea-

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.

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, and , where, is the hot wall and is cool wall and constant concentrations, and , where, > .
However, the horizontal walls are adiabatic.
The effects of combined Soret and Dufour numbers on the isothermal lines, concentration lines and horizontal and vertical velocities contours at inclination angles 0° and 45° at   The effects of Darcy number from 10 -2 to 10 -4 at inclination angles θ = 0°, permeability ratio K * = 0.1, porosity parameter ε = 0.5, Rayleigh number Ra = 10 4 , buoyancy ratio N = 1.0, Soret number Sr = 1.0, Dufour number Du = 0.12 and Lewis number Le = 1.0 has been shown in figure  5. For the value of the buoyancy ratio 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 Da * increases from 10 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.  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.    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.

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.

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.

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 heat and mass transfer and the flow characteristics inside the square cavity depended strongly on the Darcy number and Rayleigh number.