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

Abdelraheem M. Aly and Mitsuteru Asai

Additional information is available at the end of the chapter

http://dx.doi.org/10.5772/60879

#### **Abstract**

[15] Polidori G., Taïar R., Fohanno S., Mai T.H. Lodini A., 2006. Skin-friction drag analysis from the forced convection modeling in simplified underwater swimming. Journal of

[16] Vennel, R., Pease, D., Wilson, B., 2006.Wave drag on human swimmers. Journal of Bi‐

[17] Rushall, B. S., Holt, L. E., Sprigings, E. J., Cappaert, J. M., 1994. A re-evaluation of

[18] Sheehan, D. P., Laughrin, D. M., 1992. Device for qualitative measurements of hydro‐

[19] Zaïdi, H., Taïar, R., Fohanno, S., Polidori, G., 2009.An evaluation of turbulence mod‐ els in CFD simulations of underwater swimming. Series in Biomechanics 24, 5–14. [20] Lyttle, A.D., 1999. Hydrodynamics of the Human Body During the Freestyle Tumble Turn. Published doctoral dissertation, The University of Western Australia, Ned‐

[21] Bixler, B., Pease, D., Fairhurst, F., 2007. The accuracy of computational fluid dynam‐ ics analysis of the passive drag of a male swimmer. Sports Biomechanics 6, 81–98. [22] Huijing, P.A., Toussaint, H.M., Clarys, J.P., Groot, G., de., Hollander, A.P., Vervoorn, K., Mackay, R., Savelberg, H.H.C.M., 1988. Active drag related to body dimensions. In: B.E. Ungerechts, K. Reischle, K. Wilke (Eds.), Swimming Science V, 31–37. Cham‐

[23] Schleihauf, R.E., 1979. A Hydrodynamic Analysis of Swimming Propulsion, Swim‐ ming III. International Series of Sports Sciences 8, pp. 70–109. Baltimore: University

[24] Berger, M.A.M., Groot, G., Hollander, A.P., 1995. Hydrodynamic drag and lift forces

[25] Clarys, J.P., 1979. Human morphology and hydrodynamics. In: J. Terauds and E. W. Bedingfield (Eds.), International Series on Sports Science 8; Swimming III, pp. 3–41,

[26] Popa, C.V., Arfaoui A., Fohanno, S., Taïar, R., 2011 Polidori, G., 2014. Influence of a postural change of the swimmer's head in hydrodynamic performances using three dimensional CFD. Computer Methods in Biomechanics and Biomedical Engineering, 17(4), 344–351. CFD simulations of underwater swimming. Series in Biomechanics 24,

on human hand arm models. Journal of Biomechanics, 28, 125–133.

forces in swimming. Journal of Swimming Research 10, 6–30.

dynamic drag on swimmers. Journal of Swimming Research 8, 30–34.

Biomechanics 39, 2535–2541.

omechanics 39, 664–671.

64 Mass Transfer - Advancement in Process Modelling

lands, Australia.

Park Press.

5–14.

paign, IL: Human Kinetics Books.

Baltimore, MD: University Park Press.

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.

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

© 2015 The Author(s). Licensee InTech. This chapter is distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/3.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

#### **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 engi‐ neering 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 signifi‐ cantly [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 tempera‐ ture 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].

**1. Introduction**

66 Mass Transfer - Advancement in Process Modelling

transfer increases but the velocity decreases.

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 engi‐ neering 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 signifi‐ cantly [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 tempera‐ ture 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

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‐

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 invari‐ ance 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 compres‐ sible 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 Hydro‐ dynamics (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 dimen‐ sions. 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 *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 the second case considers natural convection in square annulus, in which the inner surface wall is considered to have a cool temperature *Tc* and mass *Cc* and the outer surface is exposed to a hot temperature *Th* and mass *Ch* . 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.

#### **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 *Th* and mass *Ch* and the right wall has cool temperature *Tc* and mass *Cc*. 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 *Tc* and mass *Cc* and the outer surface is exposed to a hot temperature *Th* and mass *Ch* . 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.

#### **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:

$$\nabla \cdot \mathbf{V} = 0 \tag{1}$$

Double-diffusive natural convection with Cross-Diffusion Effects using ISPH method http://dx.doi.org/10.5772/60879 69

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 *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 the second case considers natural convection in square annulus, in which the inner surface wall is considered to have a cool temperature *Tc* and mass *Cc* and the outer surface is exposed to a hot temperature *Th* and mass *Ch* . 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

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 *Th* and mass *Ch* and the right wall has cool temperature *Tc* and mass *Cc*. 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 *Tc* and mass *Cc* and the outer surface is exposed to a hot temperature *Th* and mass *Ch* . 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

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:

Ñ = . 0 **V** (1)

solutions by the finite volume method.

68 Mass Transfer - Advancement in Process Modelling

density variation in the buoyancy terms.

**3. Mathematical analysis**

**2. Problem description**

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

(b) Second case: square annulus filled with anisotropic porous media

Fig. 1 Physical model for the porous enclosure (a) First case: square cavity filled with anisotropic porous media and (b) Second case: square annulus filled with anisotropic **Figure 1.** Physical model for the porous enclosure (a) First case: square cavity filled with anisotropic porous media and (b) Second case: square annulus filled with anisotropic porous media.

$$\begin{split} \frac{\partial \mathbf{V}}{\partial t} + \left(\frac{\mathbf{V}}{\varepsilon}\right) \cdot \nabla \mathbf{V} &= \\ \frac{\varepsilon}{\rho} \Biggl[ -\nabla P + \frac{\mu}{\varepsilon} \left( \frac{\partial^2 \mathbf{V}}{\partial \boldsymbol{\alpha}^2} + \frac{\partial^2 \mathbf{V}}{\partial \boldsymbol{y}^2} \right) - \frac{\mu}{\overline{K}} \mathbf{V} - \overline{\mathbf{E}} \rho \cdot \mathbf{V} \big| \mathbf{V} \big| + \rho \left\{ \beta\_\tau \left( \mathbf{T}^\prime - \overline{\mathbf{T}}\_\cdot^\prime \right) + \beta\_c \left( \mathbf{C}^\prime - \mathbf{C}\_\cdot^\prime \right) \right\} \end{split} \tag{2}$$

�� �� � �� �

� � ��� � �

porous media.

$$
\sigma \frac{\overline{\mathcal{C}}T'}{\overline{\mathcal{C}}t} + \mathbf{V} \cdot \nabla T' = \frac{1}{\rho \mathbf{C}\_p} \left( \nabla \cdot \left( \overline{k} \nabla T' \right) \right) + D\_{\text{TC}} \left( \nabla \cdot \left( \nabla C' \right) \right), \tag{3}
$$

$$
\varepsilon \frac{\partial \mathbf{C}'}{\partial t} + \mathbf{V} \cdot \nabla \mathbf{C}' = \left(\nabla \cdot \left(\overline{D} \nabla \mathbf{C}'\right)\right) + D\_{CT} \left(\nabla \cdot \left(\nabla T'\right)\right), \tag{4}
$$

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. *DTC* is Dufour diffusivity and *DCT* is Soret diffusivity.

The second order tensor of permeability, *K*¯, Forchheimer coefficient *E*¯, thermal diffusivity, *<sup>k</sup>* ¯, and mass diffusivity, *D*¯, can be written in the following coordinate system as:

$$\overline{K} = \begin{bmatrix} K\_1 \cos^2 \theta + K\_2 \sin^2 \theta & \left( K\_1 - K\_2 \right) \sin \theta \cos \theta \\\\ \left( K\_1 - K\_2 \right) \sin \theta \cos \theta & K\_2 \cos^2 \theta + K\_1 \sin^2 \theta \end{bmatrix}'\tag{5}$$

$$\mathbf{K}^\* = \frac{\mathbf{K}\_1}{\mathbf{K}\_2} \mathbf{\prime} \tag{6}$$

$$\overline{E} = \begin{bmatrix} E\_1 \cos^2 \theta + E\_2 \sin^2 \theta & \left( E\_1 - E\_2 \right) \sin \theta \cos \theta \\\\ \left( E\_1 - E\_2 \right) \sin \theta \cos \theta & E\_2 \cos^2 \theta + E\_1 \sin^2 \theta \end{bmatrix}'\tag{7}$$

$$E\_1 = \frac{F\_1}{\sqrt{K\_1}}, E\_2 = \frac{F\_2}{\sqrt{K\_2}},\tag{8}$$

$$F^\* = \frac{F\_1}{F\_2} \,\prime \tag{9}$$

$$E^\* = \frac{E\_1}{E\_2} = \frac{F^\*}{\sqrt{K^\*}}\tag{10}$$

Substituting *K* \* (Eq. (6)) in Eq. (5) and by applying inverse, the permeability tensor can be reduced to:

$$\overline{K^{-1}} = \frac{1}{K\_1} \begin{bmatrix} \cos^2\theta + K^\*\sin^2\theta & \left(1 - K^\*\right)\sin\theta\cos\theta\\ \left(1 - K^\*\right)\sin\theta\cos\theta & K^\*\cos^2\theta + \sin^2\theta \end{bmatrix} \tag{11}$$

Substituting *E* \* . (Eq. (10)) in Eq. (7), the Forchheimer coefficient tensor can be reduced to: Double-diffusive natural convection with Cross-Diffusion Effects using ISPH method http://dx.doi.org/10.5772/60879 71

$$\overline{E} = \frac{F\_1}{\sqrt{K\_1}} \begin{bmatrix} \cos^2\theta + \frac{\sqrt{K^\*}}{F^\*}\sin^2\theta & \left(1 - \frac{\sqrt{K^\*}}{F^\*}\right)\sin\theta\cos\theta\\ \left(1 - \frac{\sqrt{K^\*}}{F^\*}\right)\sin\theta\cos\theta & \frac{\sqrt{K^\*}}{F^\*}\cos^2\theta + \sin^2\theta \end{bmatrix} \tag{12}$$

The second order of the thermal diffusivity, *k* ¯ can be written as:

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

The second order tensor of permeability, *K*¯, Forchheimer coefficient *E*¯, thermal diffusivity, *<sup>k</sup>*

( )

( )

 q

cos sin sin cos , sin cos cos sin

1 2

, , *F F E E*

\* 1 2 , *<sup>F</sup> <sup>F</sup>*

1 2

\*

(Eq. (6)) in Eq. (5) and by applying inverse, the permeability tensor can be

( )

. (Eq. (10)) in Eq. (7), the Forchheimer coefficient tensor can be reduced to:

 qq

> q

 q

\* <sup>2</sup>

2 \*2 \*

 q

é ù + - <sup>=</sup> ê ú

qq

\* \*2 2 <sup>1</sup> cos sin 1 sin cos <sup>1</sup> , 1 sin cos cos sin *K K*


 q

cos sin sin cos , sin cos cos sin

2 2

2 2

 q

*K K* = = (8)

*<sup>F</sup>* <sup>=</sup> (9)

= = (10)

 qq

 q

*<sup>K</sup>* <sup>=</sup> (6)

 qq

and mass diffusivity, *D*¯, can be written in the following coordinate system as:

*K K KK*

1 2 2 1

\* 1 2 , *<sup>K</sup> <sup>K</sup>*

1 2 2 1

é ù + - <sup>=</sup> ê ú - + ë û

*E E E E*

qq

1 2

\* 1

( )

q

*K K K*

1


*K*

*<sup>E</sup> <sup>F</sup> <sup>E</sup> E K*

 q

é ù + - <sup>=</sup> ê ú - + ë û

*KK K K*

qq

 q

2 2 1 2 12

2 2 1 2 12

q

*E E EE*

¯,

(5)

(7)

(11)

species. *DTC* is Dufour diffusivity and *DCT* is Soret diffusivity.

( )

( )

q

*K*

70 Mass Transfer - Advancement in Process Modelling

*E*

Substituting *K* \*

Substituting *E* \*

reduced to:

$$
\overline{k} = \begin{bmatrix} k\_x & 0 \\ 0 & k\_y \end{bmatrix}'\tag{13}
$$

$$k^\* = \frac{k\_y}{k\_x} \,\prime\tag{14}$$

Substituting *k* \* (Eq. (14)) in Eq. (13), the thermal conductivity tensor can be reduced to:

$$
\overline{k} = k\_x \begin{bmatrix} 1 & 0 \\ 0 & k^\* \end{bmatrix} \tag{15}
$$

In similar manner, the second order tensor of mass diffusivity, *D*¯, can be written as:

$$
\vec{D} = \begin{bmatrix} D\_x & 0 \\ 0 & D\_y \end{bmatrix}'\tag{16}
$$

$$D^\* = \frac{D\_y}{D\_x} \, \, \, \tag{17}$$

Substituting *D* \* (Eq. (17)) in Eq. (16), the mass diffusivity tensor can be reduced to:

$$
\overline{D} = D\_\times \begin{bmatrix} 1 & 0 \\ 0 & D^\* \end{bmatrix} \tag{18}
$$

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

$$\mathbf{C} = \frac{\mathbf{c}^\prime - \mathbf{c}\_l^\prime}{\Delta \mathbf{c}^\prime}, \text{ N} = \frac{\mathcal{B}\_c \left(\mathbf{c}\_h^\prime - \mathbf{c}\_l^\prime\right)}{\mathcal{B}\_T \left(T\_h^\prime - T\_l^\prime\right)}, \text{ L}\_c = \frac{\mathbf{a}\_\chi}{D\_\chi} \tag{19}$$

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

$$
\frac{\partial \mathcal{U}}{\partial X} + \frac{\partial V}{\partial Y} = 0,
\tag{20}
$$

$$\begin{aligned} &\frac{1}{\varepsilon}\frac{\partial II}{\partial \tau} + \frac{\text{II}}{\varepsilon^2}\frac{\partial II}{\partial X} + \frac{V}{\varepsilon^2}\frac{\partial V}{\partial Y} = \\ &-\frac{\partial P}{\partial X} + \frac{\text{Pr}}{\varepsilon}\left[\frac{\partial^2 U}{\partial X^2} + \frac{\partial^2 U}{\partial Y^2}\right] - \frac{\text{Pr}}{D\_s}\left[\text{II}\left(\cos^2\theta + K^\circ\sin^2\theta\right) + \left(\left(1 - K^\circ\right)\sin\theta\cos\theta\right)\right] - \\ &-\frac{F\_1}{\sqrt{D\_s}}\left[\text{II}\left(\cos^2\theta + \frac{\sqrt{K^\circ}}{F}\sin^2\theta\right) + V\left(\left(1 - \frac{\sqrt{K^\circ}}{F}\right)\sin\theta\cos\theta\right)\right]\sqrt{\mathcal{U}^2 + V^2} .\end{aligned} \tag{21}$$

$$\begin{aligned} &\frac{1}{c}\frac{\partial V}{\partial \tau} + \frac{U}{\epsilon^2}\frac{\partial V}{\partial X} + \frac{V}{\epsilon^2}\frac{\partial V}{\partial Y} = -\frac{\partial P}{\partial Y} + \frac{Pr}{\epsilon}\left[\frac{\partial^2 V}{\partial X^2} + \frac{\partial^2 V}{\partial Y^2}\right] - \\ &-\frac{Pr}{Da}\left[\mathcal{U}\Big(\Big(1-K^\prime\Big)\sin\theta\cos\theta\Big) + V\Big(K^\prime\cos^2\theta + \sin^2\theta\Big)\right] - \\ &-\frac{F\_1}{\sqrt{D\_s}}\left[\mathcal{U}\Big(\Big(1-\frac{\sqrt{K^\prime}}{F^\prime}\Big)\sin\theta\cos\theta\Big) + V\Big(\Big(1-\frac{\sqrt{K^\prime}}{F^\prime}\Big)\cos^2\theta\sin^2\theta\Big)\right]\sqrt{\mathcal{U}^2 + V^2} + Ra\Pr\{T+N\mathcal{C}\}, \end{aligned} \tag{22}$$

$$
\sigma \frac{\partial T}{\partial \tau} + \mathcal{U} \frac{\partial T}{\partial \mathbf{X}} + \mathcal{V} \frac{\partial T}{\partial \mathbf{Y}} = \frac{\partial^2 T}{\partial \mathbf{X}^2} + \mathbf{k}^\* \frac{\partial^2 T}{\partial \mathbf{Y}^2} + \mathcal{D}u \left(\frac{\partial^2 \mathbb{C}}{\partial \mathbf{X}^2} + \frac{\partial^2 \mathbb{C}}{\partial \mathbf{Y}^2}\right) \tag{23}
$$

$$\mathcal{L}\frac{\partial \mathbb{C}}{\partial \tau} + \mathcal{U}\frac{\partial \mathbb{C}}{\partial \mathcal{X}} + \mathcal{V}\frac{\partial \mathbb{C}}{\partial \mathcal{Y}} = \frac{1}{Le} \left( \frac{\partial^2 \mathbb{C}}{\partial \mathcal{X}^2} + \mathcal{D}^\* \frac{\partial^2 \mathbb{C}}{\partial \mathcal{Y}^2} \right) + \mathcal{S}r \left( \frac{\partial^2 T}{\partial \mathcal{X}^2} + \frac{\partial^2 T}{\partial \mathcal{Y}^2} \right) \tag{24}$$

where, *Pr* is the Prandtl number, *Ra* is the Rayleigh number and *Le* is Lewis number. *Du* <sup>=</sup> *DTC*Δ*<sup>C</sup>* ′ *<sup>α</sup>x*Δ*<sup>T</sup>* ′ is Dufour number and *Sr* <sup>=</sup> *DCT* <sup>Δ</sup>*<sup>T</sup>* ′ *<sup>α</sup>x*Δ*<sup>C</sup>* ′ is Soret number.

#### **3.1. Boundary treatment**

( ) ( )

*e Th l x*

( )( ) ( )

q q

2 \*2 \*

 q


¶ ¶ (20)

 qq

( )

(21)

(22)

(23)

(24)

*<sup>α</sup>x*Δ*<sup>T</sup>* ′

*ch l l x*

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

0, *U V X Y* ¶ ¶ + =

Pr Pr cos sin 1 sin cos

q

¶ ¶¶ é ù é ù - + + - + +- ê ú - ¶ Î ¶ ¶ ê ú ë û ë û

cos sin 1 sin cos ,

1 sin cos 1 cos sin ,

q q

2 2 22

2 2 22

2 2 22

*<sup>F</sup> K K <sup>U</sup> <sup>V</sup> U V Ra Pr T NC*

 q

\* 2 2 22 , *T T TT T CC U V k Du X Y X Y XY*

\*

where, *Pr* is the Prandtl number, *Ra* is the Rayleigh number and *Le* is Lewis number. *Du* <sup>=</sup> *DTC*Δ*<sup>C</sup>* ′

¶ ¶ ¶¶ ¶ ¶ ¶ æ ö + + =+ + + ç ÷ ¶¶ ¶ ¶ ¶ ¶¶ è ø

*C C C C C TT* <sup>1</sup> *U V D Sr X Y Le X Y XY*

*<sup>α</sup>x*Δ*<sup>C</sup>* ′ is Soret number.

¶ ¶ ¶ ¶ ¶ ¶¶ æ öæ ö ++= + + + ç ÷ç ÷ ¶¶ ¶ ¶ ¶ ¶¶ è øè ø

\* \* 1 2 2 2 2 \* \*

2 2

 q


é ù æ öæ ö æ ö - + +- ê ú <sup>ç</sup> ÷ç ÷ ç ÷ <sup>+</sup> ê ú <sup>ç</sup> ç ÷ ë û <sup>è</sup> øè ø è ø

*<sup>F</sup> K K <sup>U</sup> <sup>V</sup> U V*

*a*

 q

*D F F*

(( ) ) ( )

1 sin cos cos sin

é ù - - + +- ê ú ë û

\* \*2 2

\* \*

é ù ææö æö ö æ ö

\* \*

1 2 2 22

2 2 2 2

*V U V V V P Pr V V X YY X Y*

¶ ¶ ¶ ¶ ¶¶ é ù + + =- + + - ê ú ζ ¶ ¶ ¶ Î Î Î ¶ ¶ ë û

*D F F*

q q

qq

*Pr U K V K*

s

e

t

is Dufour number and *Sr* <sup>=</sup> *DCT* <sup>Δ</sup>*<sup>T</sup>* ′

t

*P UU UK K*

' ' ' ' , , '' '

b

follows:

1

1

*a*

*Da*

t

t

*a*

2 2 2 2

*U UU VV*

¶¶¶ ++=

72 Mass Transfer - Advancement in Process Modelling

ζ ¶ ¶ Î Î

2 2

*X D X Y*

q

*X Y*

b

*c c c c <sup>a</sup> CN L c D T T*

The dimensionless boundary conditions are:

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

$$\begin{aligned} \text{If } U = 0; \; V = 0; \; T = 1, \text{C} = 1 & \text{for} \quad X = 0, \; 0 \le Y \le 1\\ \text{If } U = 0; \; V = 0; \; T = 0, \text{C} = 0 & \text{for} \quad X = 1, \; 0 \le Y \le 1\\ \text{If } U = 0; \; V = 0; \; \frac{\partial T}{\partial X} = 0, \; \frac{\partial \mathcal{C}}{\partial X} = 0 & \text{for} \quad Y = 0, \; 0 \le X \le 1\\ \text{If } U = 0; \; V = 0; \; \frac{\partial T}{\partial X} = 0, \; \frac{\partial \mathcal{C}}{\partial X} = 0 & \text{for} \quad Y = 1, \; 0 \le X \le 1 \end{aligned} \tag{25}$$

**b.** Second case (square annulus filled with anisotropic porous media)

$$\begin{aligned} II &= 0; \ V = 0; \ T = 1, \mathbb{C} = 1 & \text{for} \quad X = 0, \ 0 \le Y \le 1\\ II &= 0; \ V = 0; \ T = 1, \mathbb{C} = 1\\ II &= 0; \ V = 0; \ T = 1, \mathbb{C} = 1\\ II &= 0; \ V = 0; \ T = 1, \mathbb{C} = 1\\ II &= 0; \ V = 0; \ T = 1, \mathbb{C} = 1\\ II &= 0; \ V = 0; \ T = 0, \mathbb{C} = 0\\ II &= 0; \ V = 0; \ T = 0, \mathbb{C} = 0\\ II &= 0; \ V = 0; \ T = 0, \mathbb{C} = 0\\ II &= 0; \ V = 0; \ T = 0, \mathbb{C} = 0\\ II &= 0; \ V = 0; \ T = 0, \mathbb{C} = 0\\ \end{aligned} \tag{26}$$
 
$$\begin{aligned} \text{(27)}\\ II &= 0; \ V = 0; \ T = 0, \mathbb{C} = 0\\ \end{aligned} \qquad \text{for } \quad \begin{cases} L - d \mid \end{cases} / / 2 \le Y \le \left(L + d\right)/2, \text{ } X = \left(L + d\right)/2, \text{ } Y = \left(L - d\right)/2, \text{ } Z\\ \end{aligned} \tag{27}$$

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

$$N\text{\AA} = -\left.\frac{\partial T}{\partial \mathbf{X}}\right|\_{\mathbf{X}=0} \text{\AA} \tag{27}$$

$$Sh = -\left.\frac{\partial \mathbf{C}}{\partial \mathbf{X}}\right|\_{\mathbf{X}=0}.\tag{28}$$

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:

$$\overline{N\mu} = -\frac{1}{L} \frac{\mathbf{i}}{\mathbf{j}} \left( \frac{\partial T}{\partial X} \right) dY. \tag{29}$$

$$\overline{Sh} = -\frac{1}{L} \frac{\mathbf{i}}{\mathbf{j}} \left( \frac{\partial \mathbf{C}}{\partial \mathbf{X}} \right) d\mathbf{Y} \tag{30}$$

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

Horizontal hot walls:

$$\left. \mathcal{N} \Lambda \mathfrak{u} - -\frac{\partial \mathcal{T}}{\partial \mathcal{Y}} \right|\_{\mathcal{T} = \mathcal{T}\_h} \text{ \( \mathcal{S} \mathfrak{h} = -\frac{\partial \mathcal{C}}{\partial \mathcal{Y}} \Big|\_{\mathcal{T} = \mathcal{T}\_h} \text{ \)}\tag{31}$$

Vertical hot walls:

$$\delta N \mu = -\left. \frac{\partial T}{\partial \mathbf{X}} \right|\_{T = T\_h} \text{ , } Sh = -\left. \frac{\partial \mathbf{C}}{\partial \mathbf{X}} \right|\_{T = T\_h} \text{ , } \tag{32}$$

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:

$$\overline{Nu\_L} = -\frac{1}{L} \Big|\_{0}^{L} \Big( \frac{\partial T}{\partial X} \Big) dY,\\ \overline{Sh\_L} = -\frac{1}{L} \Big|\_{0}^{L} \Big( \frac{\partial \mathcal{C}}{\partial X} \Big) dY,\tag{33}$$

Average Nusselt and Sherwood numbers at the right hot surface:

$$\overline{Nu\_{\mathbb{R}}} = -\frac{1}{L} \prod\_{0}^{\mathbb{L}} \left( \frac{\partial T}{\partial \mathbf{X}} \right) dY\_{\prime} \overline{Sh\_{\mathbb{R}}} = -\frac{1}{L} \prod\_{0}^{\mathbb{L}} \left( \frac{\partial \mathbb{C}}{\partial \mathbf{X}} \right) dY\_{\prime} \tag{34}$$

Average Nusselt and Sherwood numbers at the top hot surface:

$$\overline{Nu\_{\rm T}} = -\frac{1}{L} \stackrel{\leftarrow}{\int} \left( \frac{\partial T}{\partial \mathcal{Y}} \right) dX\_{\prime} \overline{Sh\_{\rm R}} = -\frac{1}{L} \stackrel{\leftarrow}{\int} \left( \frac{\partial \mathcal{C}}{\partial \mathcal{Y}} \right) dX\_{\prime} \tag{35}$$

Average Nusselt and Sherwood numbers at the bottom hot surface:

$$\overline{Nu\_{\rm B}} = -\frac{1}{L} \overset{\mathcal{L}}{\underset{0}{\rm f}} \Big| \left( \frac{\partial T}{\partial Y} \right) dX,\\ \overline{Sh\_{\rm B}} = -\frac{1}{L} \overset{\mathcal{L}}{\underset{0}{\rm f}} \Big| \left( \frac{\partial \mathcal{C}}{\partial Y} \right) dX,\tag{36}$$

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

$$\overline{Nu\_{\text{Total}}} = \frac{1}{4} \overline{(Nu\_{\text{L}} + \overline{Nu\_{\text{R}}} + \overline{Nu\_{\text{T}}} + \overline{Nu\_{\text{B}}})} \tag{37}$$

$$\overline{Sh\_{Total}} = \frac{1}{4} (\overline{Sh\_L} + \overline{Sh\_R} + \overline{Sh\_T} + \overline{Sh\_B}) \tag{38}$$

#### **4. Numerical method**

0 <sup>1</sup> *<sup>L</sup> <sup>C</sup> Sh dY L X* æ ö ¶ = - ç ÷

as follows:

Horizontal hot walls:

74 Mass Transfer - Advancement in Process Modelling

Vertical hot walls:

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:

*L L*

*R R*

*T R*

*B B*

Average Nusselt and Sherwood numbers at the bottom hot surface:

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

 , , *h h T T T T*

 , , *h h T T T T*

For the square annulus case, the average Nusselt number and average Sherwood number are

0 0 1 1 , , *L L*

0 0 1 1 , , *L L*

0 0 1 1 , , *L L*

0 0 1 1 , , *L L*

*T C Nu dX Sh dX LY LY*

*T C Nu dX Sh dX LY LY*

*T C Nu dY Sh dY LX LX*

*T C Nu dY Sh dY LX LX*

*X X* = =

*Y Y* = =

*T C Nu Sh*

*T C Nu Sh*

è ø ¶ ò (30)

¶ ¶ <sup>=</sup> - = - ¶ ¶ (31)

¶ ¶ <sup>=</sup> - = - ¶ ¶ (32)

æ ö ¶ ¶æ ö <sup>=</sup> - = ç ÷ - ç ÷ è ø ¶ ¶è ø ò ò (33)

æ ö ¶ ¶æ ö <sup>=</sup> - = ç ÷ - ç ÷ è ø ¶ ¶è ø ò ò (34)

æ ö ¶ ¶æ ö <sup>=</sup> - = ç ÷ - ç ÷ è ø ¶ ¶è ø ò ò (35)

æ ö ¶ ¶æ ö <sup>=</sup> - = ç ÷ - ç ÷ è ø ¶ ¶è ø ò ò (36)

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 orthog‐ onal 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 follow‐ ing criteria of convergence was fulfilled:

$$\sum\_{i,j} \left| \mathcal{Z}\_{i,j}^{\text{new}} - \mathcal{Z}\_{i,j}^{\text{old}} \right| \le 10^{-6} \,\text{\textdegree} \,\text{\textdegree} \,\text{\textdegree} \,\text{\textdegree} \,\text{\textdegree} \,\text{\textdegree} \,\text{\textdegree} \,\text{\textdegree} \,\text{\textdegree} \,\text{\textdegree} \,\text{\textdegree} \,\text{\textdegree} \,\text{\textdegree} \,\text{\textdegree} \,\text{\textdegree} \,\text{\textdegree} \,\text{\textdegree} \,\text{\textdegree} \,\text{\textdegree} \,\text{\textdegree} \,\text{\textdegree} \,\text{\textdegree} \,\text{\textdegree} \,\text{\textdegree} \,\text{\textdegree} \,\text{\textdegree} \,\text{\textdegree} \,\text{\textdegree} \,\text{\textdegree} \,\text{\textdegree} \,\text{\textdegree} \,\text{\textdegree} \,\text{\textdegree} \,\text{\textdegree} \,\text{\textdegree} \,\text{\textdegree} \,\text{\textdegree} \,\text{\textdegree} \,\text{\textdegree} \,\text{\textdegree} \,\text{\textdegree} \,\text{\textdegree} \,\text{\textdegree} \,\text{\textdegree} \,\text{\textdegree} \,\text{\textdegree} \,\text{\textdegree} \,\text{\textdegree} \,\text{\textdegree} \,\text{\textdegree} \,\text{\textdegree} \,\text{\textdegree} \,\text{\textdegree} \,\text{\textdegree} \,\text{\textdegree} \,\text{\textdegree} \,\text{\textdegree} \,\text{\text{\textdegree}} \,\text{\text{\textdegree}} \,\text{\text{\textdegree}} \,\text{\text{\textdegree}} \,\text{\text{\textdegree}} \,\text{\text{\textdegree}} \,\text{\text{\textdegree}} \,\text{\text{\text{\text{\text$$

#### **4.2. ISPH method**

The ISPH algorithm is implemented in a semi-implicit form in order to solve the incompres‐ sible 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:

$$\begin{aligned} &\frac{1}{\pi} \left[ \frac{\mathbf{U}^{n+1} - \mathbf{U}^{n}}{\Delta t} \right] = -\left( \frac{\partial \mathbf{P}}{\partial \mathbf{X}} \right)^{n+1} + \frac{\Pr}{\pi} \left( \frac{\partial^{2} \mathbf{U}}{\partial \mathbf{X}^{2}} + \frac{\partial^{2} \mathbf{U}}{\partial \mathbf{Y}^{2}} \right)^{n} - \\ &- \left( \frac{\Pr}{Da} (\cos^{2} \theta + \mathbf{K}^{\prime} \sin^{2} \theta) + \frac{F\_{1}}{\sqrt{Da}} \left( \cos^{2} \theta + \frac{\sqrt{\mathbf{K}^{\prime}}}{\mathbf{F}^{\prime}} \sin^{2} \theta \right) \sqrt{\mathbf{U}^{2} + \mathbf{V}^{2}} \right) \mathbf{U}^{n+1} - \\ &- \left( \left( \frac{\Pr}{Da} (1 - \mathbf{K}^{\prime}) \sin \theta \cos \theta \right) + \frac{F\_{1}}{\sqrt{Da}} \left( \left( 1 - \frac{\sqrt{\mathbf{K}^{\prime}}}{\mathbf{F}^{\prime}} \right) \sin \theta \cos \theta \right) \sqrt{\mathbf{U}^{2} + \mathbf{V}^{2}} \right) \mathbf{V}^{n} \end{aligned} \tag{40}$$

$$\begin{split} &\frac{1}{\mathsf{er}} \left[ \frac{V^{n+1} - V^{n}}{\Delta t} \right] = -\left( \frac{\partial P}{\partial Y} \right)^{n+1} + \frac{\Pr}{\mathsf{er}} \left( \frac{\partial^{2} U}{\partial X^{2}} + \frac{\partial^{2} U}{\partial Y^{2}} \right)^{n} - \\ & - \left( \frac{\Pr}{Da} \left( \left( 1 - K^{\prime} \right) \sin \theta \cos \theta \right) + \frac{F\_{1}}{\sqrt{Da}} \left( 1 - \frac{\sqrt{K^{\prime}}}{F} \sin^{2} \theta \cos \theta \right) \sqrt{\mathcal{U}^{2} + V^{2}} \right) \mathcal{U}^{n} - \\ & - \left( \frac{\Pr}{Da} \left( K^{\prime} \cos^{2} \theta + \sin^{2} \theta \right) + \frac{F\_{1}}{\sqrt{Da}} \left( \frac{\sqrt{K^{\prime}}}{F} \cos^{2} \theta + \sin^{2} \theta \right) \sqrt{\mathcal{U}^{2} + V^{2}} \right) \mathcal{V}^{n+1} + \\ & + \mathsf{Pr} \operatorname{Ra} \{ T^{n} + \mathcal{N} \mathcal{C}^{n} \} \end{split} \tag{41}$$

For simplicity,

$$\begin{aligned} \text{Porys} &= \frac{Pr}{Da} \Big(K' \cos^2 \theta + \sin^2 \theta\right) + \frac{F\_1}{\sqrt{Da}} \left(\frac{\sqrt{K'}}{F} \cos^2 \theta + \sin^2 \theta\right) \Big(\sqrt{L^2 + V^2}\Big)^n, \\\ \text{Porys} &= \left(\frac{Pr}{Da} \Big(1 - K'\right) \sin \theta \cos \theta\right) + \frac{F\_1}{\sqrt{Da}} \left(\Big(1 - \frac{\sqrt{K'}}{F}\Big) \sin \theta \cos \theta\right) \Big(\sqrt{L^2 + V^2}\Big)^n, \end{aligned} \tag{42}$$

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

$$\begin{split} & \frac{(1 + \epsilon \Delta t \text{ Porsx}) \text{LI}^{n+1} - \left(1 + \epsilon \Delta t \text{ Porsx}\right) \text{LI}^{\*} + \left(1 + \epsilon \Delta t \text{ Porsx}\right) \text{LI}^{\*} - \text{II}^{n}}{\epsilon \Delta t} = \\ & = -\left(\frac{\partial P}{\partial X}\right)^{n+1} + \frac{Pr}{\epsilon} \left(\frac{\partial^{2} \mathcal{U}}{\partial X^{2}} + \frac{\partial^{2} \mathcal{U}}{\partial Y^{2}}\right)^{n} - \text{Porsx}\,\,\,\,\,\text{V}^{n}, \end{split} \tag{43}$$

Double-diffusive natural convection with Cross-Diffusion Effects using ISPH method http://dx.doi.org/10.5772/60879 77

$$\begin{aligned} & \left( \frac{1 + \epsilon \Delta t \text{ Porsy}}{\epsilon} \right) V^{n+1} - \left( 1 + \epsilon \Delta t \text{ Porsy} \right) V^\* + \left( 1 + \epsilon \Delta t \text{ Porsy} \right) V^\* - V^n = \\ & \frac{\epsilon \Delta t}{\epsilon} = -\left( \frac{\partial P}{\partial Y} \right)^{n+1} + \frac{\text{Pr}}{\epsilon} \left( \frac{\partial^2 V}{\partial X^2} + \frac{\partial^2 V}{\partial Y^2} \right)^n - \text{Porsz} \text{ } L^n + \text{Pr} \, \text{Ra} \left( T^n + N \text{ } \mathbb{C}^n \right) \end{aligned} \tag{44}$$

#### **Step 1.**Predict the velocity

*4.2.1. Temporal discretization*

76 Mass Transfer - Advancement in Process Modelling

( )

(( ) )

1 Pr

+ +

( )

q

( )

*Ra T NC*

*n n*

Pr

For simplicity,

+ +

qq

q q

1 Pr

+ +

1 1 2 2

éù æ ö - ¶ ¶¶ æ ö ê ú =- + + - ç ÷ ç ÷ ÎD ¶ Î è ø ¶ ¶ ëû è ø

*Da Da F*

*Da Da F*

1 1 2 2

*V V P UU t Y X Y*

éù æ ö - ¶ ¶¶ æ ö ê ú =- + + - ç ÷ ç ÷ ÎD ¶ Î è ø ¶ ¶ ëû è ø

*Da Da F*

 q

*Da Da F*

q q

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

+

1 2 2

*X X Y*

ò

æ ö ¶ ¶¶ æ ö =- + + - ç ÷ ç ÷ è ø ¶ ¶ ¶ è ø

+

ò

*n n*

 q

q

*Da Da F*

2 2

*P Pr U U <sup>V</sup>*

*Da Da F*

q q

*<sup>n</sup> <sup>n</sup> n n*

Pr 1 sin cos 1 sin cos

Pr cos sin cos sin

*UU P UU t X X Y*

*<sup>n</sup> <sup>n</sup> n n*

Pr cos sin cos sin

Pr 1 sin cos 1 sin cos

( )

follows:

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

2 2

2 \*2 1 2 2 22 1

 q

æ ö æ ö - + + + +- ç ÷ ç ÷ è ø è ø æ ö æ ö æ ö æ ö -- + - ç ÷ ç ÷ ç ÷ ç ÷ <sup>+</sup> ç ÷ ç ÷ ç ÷ è ø è ø è ø è ø

*<sup>F</sup> <sup>K</sup> <sup>K</sup> U VV*

2 2

*<sup>F</sup> <sup>K</sup> <sup>K</sup> U VU*

\* \*2 2 1 2 2 22 1 \*

( ) ( )

*Pr <sup>F</sup> <sup>K</sup> <sup>K</sup> U V*

*t U t U t UU t*

òò

+ D - +D + +D - <sup>=</sup> <sup>D</sup>

Porsy cos sin cos sin ,

æ ö æ ö æ ö =- + - ç ÷ ç ÷ ç ÷ <sup>+</sup> ç ÷ ç ÷ è ø è ø è ø

( ) ( ) ( ) 1\* \*

ò

1 Porsx 1 Porsx 1 Porsx

= ++ + + ç ÷

*Pr <sup>F</sup> <sup>K</sup> <sup>K</sup> U V*

Porsr 1 sin cos 1 sin cos ,

\* \*2 2 1 2 2 22 \*

( ) ( )

Porsr ,

*n n*

*n*

\* \* 1 2 2 \*

q

æ ö

è ø

*<sup>F</sup> <sup>K</sup> <sup>K</sup> U VV*

æ ö æ ö -- +- ç ÷ ç ÷ + è ø è ø æ ö æ ö - ++ + + + ç ÷ ç ÷ è ø è ø

\* \* 1 2 2 2 \*

q

*<sup>F</sup> <sup>K</sup> <sup>K</sup> U VU*

\*

 q

q q

q q

> q

> > q

q q *n*

+

*n*

*n*

*n*

+

*n*

*n*

(40)

(41)

(42)

(43)

\*

\* \* 1 2 2 \*

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:

$$\left(1+\epsilon\Delta t\text{ Porsx}\right)\mathcal{U}^\* = \mathcal{U}^n + \epsilon\Delta t \left(\frac{\Pr\left(\hat{\mathcal{O}}^2\mathcal{U} + \frac{\hat{\mathcal{O}}^2\mathcal{U}}{\hat{\mathcal{O}}\mathcal{X}^2}\right)^n - \text{Porsx}\,\mathcal{V}^n\right),\tag{45}$$

$$\left(1+\epsilon\Delta t \text{ Porsy}\right)V^\* = \left.V^\* + \epsilon\Delta t \left(\frac{\Pr\left(\hat{\sigma}^2 V + \hat{\sigma}^2 V\right)^u}{\epsilon}\right)^u + \Pr\text{Ra}\left(T^u + N \text{ C}^u\right) - \text{Porsr I}\,\text{I}^u\right). \tag{46}$$

**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:

$$\frac{\partial^2 P^{n+1}}{\partial X^2} + \frac{\partial^2 P^{n+1}}{\partial Y^2} = \left( \frac{(1 + \epsilon \Delta t \text{ Pasx})}{\epsilon \Delta t} \bigg| \frac{\partial \mathcal{U}}{\partial X} \right) + \frac{\left(1 + \epsilon \Delta t \text{ Possy} \right) \left(\frac{\partial V}{\partial Y} \right)}{\epsilon \Delta t} \bigg| \frac{\partial V}{\partial Y} \bigg) + \gamma \frac{\rho^0 - \epsilon \rho^\*}{\rho^0 \Delta t^2} \tag{47}$$

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:

$$\left(1+\epsilon\Delta t\text{ Porsx}\right)\text{LI}^{n+1} = \left(1+\epsilon\Delta t\text{ Porsx}\right)\text{LI}^{\*} - \epsilon\Delta t\left[\frac{\partial P}{\partial \mathbf{X}}\right]^{n+1}\text{.}\tag{48}$$

$$\left(1+\epsilon\Lambda t \text{ Porsy}\right)V^{n+1} = \left(1+\epsilon\Lambda t \text{ Porsy}\right)V^{\*} - \epsilon\Lambda t \left[\frac{\partial P}{\partial Y}\right]^{n+1} \tag{49}$$

**Step 4.**Thermal flow problems:

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

$$
\sigma \left( \frac{T^{u+1} - T^u}{\Delta \tau} \right) = \left( \left( \frac{\partial^2 T}{\partial X^2} + \mathbf{k}^\* \frac{\partial^2 T}{\partial Y^2} \right) \right)^u + Du \left( \frac{\partial^2 C}{\partial X^2} + \frac{\partial^2 C}{\partial Y^2} \right)^u,\tag{50}
$$

**Step 5.**Concentration flow problems:

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

$$\varepsilon \left( \frac{\mathbb{C}^{n+1} - \mathbb{C}^{n}}{\Delta \tau} \right) = \frac{1}{Le} \left( \left( \frac{\partial^{2} \mathbb{C}}{\partial \mathbf{X}^{2}} + \mathbf{D}^{\*} \frac{\partial^{2} \mathbb{C}}{\partial \mathbf{Y}^{2}} \right) \right)^{\mathrm{u}} + \mathrm{Sr} \left( \frac{\partial^{2} T}{\partial \mathbf{X}^{2}} + \frac{\partial^{2} T}{\partial \mathbf{Y}^{2}} \right)^{\mathrm{u}}.\tag{51}$$

#### *4.2.2. SPH Formulations*

A spatial discretization using scattered particles, which is based on the SPH, is summarized. First, a physical scalar function *ϕ*(*Xi* , *τ*) at a sampling point *Xi* can be represented by the following integral form:

$$\phi\left(\mathbf{X}\_{i'}\boldsymbol{\tau}\right) = \int \mathcal{W}\left(\mathbf{X}\_{i} - \mathbf{X}\_{j'}h\right)\phi\left(\mathbf{X}\_{j'}\boldsymbol{\tau}\right)d\boldsymbol{v} = \int \mathcal{W}\left(\mathbf{R}\_{ij'}h\right)\phi\left(\mathbf{X}\_{j'}\boldsymbol{\tau}\right)d\boldsymbol{v},\tag{52}$$

where, *W* is a weight function called by smoothing kernel function in the SPH literature. In the smoothing kernel function, *Rij* = | *Xi* − *Xj* | 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.

$$\phi(X\_{i'}\tau) \approx \left\langle \phi\_i \right\rangle = \sum\_{j} \frac{m\_j}{\rho\_j} \mathcal{W}(R\_{ij'}h)\phi(X\_{j'}\tau)\_{\prime} \tag{53}$$

where, the subscripts *i* and *j* indicate positions of labeled particle, *ρ<sup>j</sup>* and *mj* 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:

$$\nabla \phi(X\_i) \approx \left< \nabla \phi\_i \right> = \frac{1}{\mathcal{P}\_i} \sum\_j m\_j (\phi\_j - \phi\_i) \nabla \mathcal{W}(R\_{ij}, h). \tag{54}$$

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

$$\nabla \phi(X\_i) \approx \left\langle \nabla \phi\_i \right\rangle = \rho\_i \sum\_j m\_j \left(\frac{\phi\_j}{\rho\_j^2} + \frac{\phi\_i}{\rho\_i^2}\right) \nabla \mathcal{W}(R\_{ij}, h). \tag{55}$$

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

$$\mathcal{W}\left(R\_{\boldsymbol{y}},h\right) = \mathcal{B}\_{d}\begin{cases} \left(3 - \frac{R\_{\boldsymbol{y}}}{h}\right)^{5} - 6\left(2 - \frac{R\_{\boldsymbol{y}}}{h}\right)^{5} + 15\left(1 - \frac{R\_{\boldsymbol{y}}}{h}\right)^{5} & 0 \le R\_{\boldsymbol{y}} < h\\ \left(3 - \frac{R\_{\boldsymbol{y}}}{h}\right)^{5} - 6\left(2 - \frac{R\_{\boldsymbol{y}}}{h}\right)^{5} & h \le R\_{\boldsymbol{y}} < 2h\\ \left(3 - \frac{R\_{\boldsymbol{y}}}{h}\right)^{5} & 2h \le R\_{\boldsymbol{y}} < 3h\\ 0 & r\_{\boldsymbol{y}} \ge 3h \end{cases} \tag{56}$$

where *βd* is 7 / <sup>478</sup>*π<sup>h</sup>* <sup>2</sup> and 3 / <sup>358</sup>*π<sup>h</sup>* <sup>3</sup> , respectively, in two- and three dimension space.

#### *4.2.3. Projection-based ISPH formulations*

**Step 4.**Thermal flow problems:

78 Mass Transfer - Advancement in Process Modelling

s

**Step 5.**Concentration flow problems:

e

First, a physical scalar function *ϕ*(*Xi*

 t

*4.2.2. SPH Formulations*

following integral form:

f

t

t

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

1 22 2 2 \*

<sup>+</sup> æ öæ ö æ ö - ¶ ¶ ¶¶ æ ö ç ÷ç ÷ ç ÷ =+ + + ç ÷ <sup>D</sup> ¶ ¶ ¶¶ è øè ø è ø è ø

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

*<sup>n</sup> <sup>n</sup> n n TT T T CC k Du*

1 2 2 2 2 \*

<sup>+</sup> æ öæ ö æ ö - ¶ ¶ ¶¶ æ ö ç ÷ç ÷ ç ÷ = + ++ ç ÷ <sup>D</sup> ¶ ¶ ¶¶ è øè ø è ø è ø

*<sup>n</sup> <sup>n</sup> n n CC C C TT D Sr*

2 2 2 2 ,

2 2 2 2 <sup>1</sup> .

, *τ*) at a sampling point *Xi*

f

» = å (53)

and *mj*

mean density and

ò ò (52)

f t

Ñ »Ñ = - Ñ å (54)

 t

*Le X Y XY*

A spatial discretization using scattered particles, which is based on the SPH, is summarized.

( ) *X W X X h X dv W R h X dv <sup>i</sup>* , , , , ,, ( )( ) ( )( ) *i j <sup>j</sup> ij <sup>j</sup>*

( ,) ( , ) ( , ), *<sup>j</sup> i i ij j j j m*

r

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

> <sup>1</sup> ( ) ( ) ( , ). *<sup>i</sup> <sup>i</sup> jj i ij j i*

*X m WR h*

r

 ff

*X WR h X*

 t

where, *W* is a weight function called by smoothing kernel function in the SPH literature. In the smoothing kernel function, *Rij* = | *Xi* − *Xj* | 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.

 =- = f

ftf

using the above defined SPH approximation as follows:

f

where, the subscripts *i* and *j* indicate positions of labeled particle, *ρ<sup>j</sup>*

 f (50)

(51)

can be represented by the

*X Y XY*

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:

$$\nabla P(\mathbf{X}\_i) \approx \left\langle \nabla P\_i \right\rangle = \rho\_i \sum\_j m\_j \left(\frac{P\_j}{\rho\_j^2} + \frac{P\_i}{\rho\_i^2}\right) \nabla W(\mathbf{R}\_{ij}, h)\_{\prime} \tag{57}$$

$$\nabla \cdot \mathbf{U} \left( \mathbf{X}\_{i} \right) \approx \langle \nabla \cdot \mathbf{U}\_{i} \rangle = \frac{1}{\rho\_{i}} \sum\_{j} m\_{j} \left( \mathbf{U}\_{j} - \mathbf{U}\_{i} \right) \cdot \nabla W \left( \mathbf{R}\_{\parallel}, h \right)\_{i} \tag{58}$$

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:

$$\nabla \cdot \left( \nu \nabla U \right) \left( X\_i \right) \approx \langle \nabla \cdot \left( \nu \nabla \cdot U\_i \right) \rangle = \sum\_j m\_j \left( \frac{\rho\_i \nu\_i + \rho\_j \nu\_j}{\rho\_i \rho\_j} \frac{R\_{ij} \cdot \nabla W \left( \left| R\_i - R\_j \right| \cdot h \right)}{R\_{ij}^2 + \eta^2} \right) U\_{ij} \tag{59}$$

where, *η* is a parameter to avoid a zero dominator, and its value is usually given by *η* <sup>2</sup> =0.0001 *h* <sup>2</sup> . For the case of *υ<sup>i</sup>* =*υ<sup>j</sup>* and *ρ<sup>i</sup>* =*ρ<sup>j</sup>* , the Laplacian term is simplified as:

$$\langle \nabla \cdot \left( \nu \nabla \cdot \mathbf{U}\_i \right) \rangle = \frac{2 \nu\_i}{\rho\_i} \sum\_j m\_j \left( \frac{R\_{\vec{\eta}} \cdot \nabla W \left( \left| R\_i - R\_j \right| \, \mu \right)}{R\_{\vec{\eta}}^2 + \eta^2} \right) U\_{ij\prime} \tag{60}$$

In this study, Laplacian of velocity is given by:

$$\nabla^2 \mathbf{U} \left( \mathbf{X}\_i \right) \approx \langle \nabla^2 \mathbf{U}\_i \rangle = \frac{2}{\rho\_i} \sum\_j m\_j \left( \frac{R\_{\vec{\eta}'} \cdot \nabla \mathcal{W} \left( \left| R\_i - R\_j \right| \, \_j \text{ft} \right)}{R\_{\vec{\eta}}^2 + \eta^2} \right) \mathbf{U}\_{\vec{\eta}'} \tag{61}$$

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

$$\nabla^2 P(X\_i) \approx \left\langle \nabla^2 P\_i \right\rangle = \frac{2}{\rho\_i} \sum\_j m\_j \left( \frac{P\_{\vec{\eta}} R\_{\vec{\eta}} \cdot \nabla W \left( \left| R\_i - R\_j \right| \, \mu \right)}{R\_{\vec{\eta}}^2 + \eta^2} \right) \tag{62}$$

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:

$$\nabla^2 T \{ X\_i \} \approx \langle \nabla^2 T\_i \rangle = \sum\_j m\_j \left( \frac{\rho\_i + \rho\_j}{\rho\_i \rho\_j} \frac{T\_{ij} \ R\_{ij} \cdot \nabla \mathcal{W} \{ \left| R\_i - R\_j \right|, h \}}{R\_{ij}^2 + \eta^2} \right) \tag{63}$$

$$\nabla^2 \mathbb{C} \{ X\_i \} \approx \langle \nabla^2 \mathbb{C}\_i \rangle = \sum\_j m\_j \left( \frac{\rho\_j + \rho\_j}{\rho\_j \rho\_j} \frac{\mathbb{C}\_{ij} \cdot \mathbb{R}\_{ij} \cdot \nabla W(\left| \mathbb{R}\_i - \mathbb{R}\_j \right| \cdot h)}{\mathbb{R}\_{ij}^2 + \eta^2} \right), \tag{64}$$

#### **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 *Da* =10−<sup>2</sup> at buoyancy ratio *N* =0, has been introduced in table 1. In this table, the Rayleigh numbers was taken as Ra=104 with Darcy parameter, *Da* =10−<sup>2</sup> . 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].


**Table 1.** 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, *Da* =10−<sup>2</sup> with buoyancy ratio *N* =0.

Comparison of isotherms lines between FVM and ISPH methods under the effect of Rayleigh number and permeability ratio at Darcy number *Da* =10−<sup>2</sup> , porosity parameter *ε* =0.5, buoy‐ ancy 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.

### **6. Results and discussion**

where, *η* is a parameter to avoid a zero dominator, and its value is usually given by

2 2

h

2 2 , <sup>2</sup> , *ij i j*

h

2 2

<sup>+</sup> è ø

h

2 2 ( ,) , *i j ij ij i j*

2 2 ( ,) , *i j ij ij i j*

h

*CR WR R h*

*R*

<sup>+</sup> è ø

æ ö + ×Ñ -

h

å (63)

å (64)

*TR WR R h*

*R*

<sup>+</sup> è ø

æ ö + ×Ñ -

å (62)

*PR W R R h*

*R WR R h*

<sup>+</sup> è ø *U U* å *<sup>U</sup>* (61)

*R WR R h*

<sup>+</sup> è ø

*R*

, the Laplacian term is simplified as:

*U U* å (60)

and *ρ<sup>i</sup>* =*ρ<sup>j</sup>*

( ) ( )

æ ö ×Ñ - ç ÷ áÑ Ñ ñ =

( ) ( ) 2 2

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

( ) 2 2

æ ö ×Ñ - ç ÷ Ñ »Ñ =

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

r r

r r

> r r

Ñ » áÑ ñ= ç ÷

r r

Ñ » áÑ ñ= ç ÷

, <sup>2</sup> ( ) . *ij ij i j*

*R*

*j i ij*

The PPE after SPH interpolation is solved by a preconditioned (diagonal scaling) Conjugate

*j i j ij*

*j i j ij*

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

æ ö ×Ñ - ç ÷ Ñ » áÑ ñ =

*i i j ij j i ij*

*R*

*m*

u

r

*X m*

*ii j*

Gradient (PCG) method [46] with a convergence tolerance (=1.0×10-9).

r

*PX P m*

*i ij*

*i ij*

*CX C m*

*TX T m*

r

, <sup>2</sup> . . , *ij i j <sup>i</sup> i j ij j i ij*

*η* <sup>2</sup> =0.0001 *h* <sup>2</sup>

. For the case of *υ<sup>i</sup>* =*υ<sup>j</sup>*

80 Mass Transfer - Advancement in Process Modelling

u

In this study, Laplacian of velocity is given by:

( ) 2 2

( ) 2 2

**5. Numerical validation tests**

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, *Th* and *Tc*, where, *Th* is the hot wall and *Tc* is cool wall and constant concentrations, *Ch* and *Cc*, where, *Ch* > *Cc*. 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‐

Darcy number �� � ����, porosity parameter� � ��5 and inclination angle � � �° in the absence of Soret and Dufour numbers effects. **Figure 2.** Comparison of the isotherms lines between FVM method (—) and ISPH method (----) at buoyancy ratio *N* =0, Darcy number *Da* =10−<sup>2</sup> , porosity parameter *ε* =0.5 and inclination angle *θ* =0° in the absence of Soret and Dufour numbers effects.

Fig. 2 Comparison of the isotherms lines between FVM method (**—**) and ISPH method (‐‐‐‐) at buoyancy ratio ���,

bility ratio *K* \* =0.1, Darcy parameter *Da* =10−<sup>2</sup> , porosity parameter *ε* =0.5, Rayleigh number *Ra* =10<sup>4</sup> , buoyancy ratio *N* =1.0 and Lewis number *Le* =1.0 has been shown in figure 3. **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

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 anticlock‐ wise 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. 17 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, �� 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 � � �° and 45° at

circular cells next to left and right walls of the enclosure. Also, the horizontal and vertical

decreasing Dufour number.

coupled Soret and Dufour numbers at inclination angle � � ��, permeability ratio �<sup>∗</sup> � ���, �� � ����� � � ���, �� � ���, � � ��� and �� � ���. Fig. 4 presents the effects of Soret and Dufour numbers on temperature, solutal concentration, **Figure 3.** Presents isothermal lines, concentration lines, horizontal and vertical velocities contours under the effects of coupled Soret and Dufour numbers at inclination angle *θ* =0° , permeability ratio *K* \* =0.1, *Da* =10−<sup>2</sup> , *ε* =0.5, *Ra* =10<sup>4</sup> , *N* =1.0 and *Le* =1.0.

bility ratio *K* \*

Dufour numbers effects.

**6. RESULTS AND DISCUSSION**

*N* =0, Darcy number *Da* =10−<sup>2</sup>

82 Mass Transfer - Advancement in Process Modelling

*Ra* =10<sup>4</sup>

=0.1, Darcy parameter *Da* =10−<sup>2</sup>

numbers have been discussed in details using ISPH method.

increasing Soret number with decreasing Dufour number.

However, the horizontal walls are adiabatic.

, buoyancy ratio *N* =1.0 and Lewis number *Le* =1.0 has been shown in figure 3.

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

**Figure 2.** Comparison of the isotherms lines between FVM method (—) and ISPH method (----) at buoyancy ratio

(a) �<sup>∗</sup> � ���� �� � ��� �� � ��� �� � ���

(b) �<sup>∗</sup> � ��� �� � ��� �� � ��� �� � ��� Fig. 2 Comparison of the isotherms lines between FVM method (**—**) and ISPH method (‐‐‐‐) at buoyancy ratio ���, Darcy number �� � ����, porosity parameter� � ��5 and inclination angle � � �° in the absence of Soret and Dufour numbers effects.

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 anticlock‐ wise 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

 **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, �� and ��, where, �� is the hot wall and �� is cool wall and constant concentrations, �� and �� , where, ��>�� .

17

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

, porosity parameter *ε* =0.5, Rayleigh number

, porosity parameter *ε* =0.5 and inclination angle *θ* =0° in the absence of Soret and

18 horizontal velocity profiles and vertical velocity profiles. It is noted that, an increase in Soret 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 concen‐ tration 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 *θ* =0° , permeability ratio *K* \* =0.1, porosity parameter *ε* =0.5, Rayleigh number *Ra* =10<sup>4</sup> , 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

Fig. 4 Presents the effects of combined Soret and Dufour numbers on the temperature, concentration, horizontal and vertical velocity profiles, respectively. **Figure 4.** Presents the effects of combined Soret and Dufour numbers on the temperature, concentration, horizontal and vertical velocity profiles, respectively.

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 102 to 104 , 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. The effects of Darcy number from 10-2 to 10-4 at inclination angles� � �� , permeability ratio �<sup>∗</sup> � ���, porosity parameter � � ���, Rayleigh number��� � ��� , buoyancy ratio � � ���, Soret number �� � ���, Dufour number �� � ���� and Lewis number �� � ��� 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 ��<sup>∗</sup> increases from 102

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=103 and 104 . 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 *X* increases. It is observed that, there is a region in the middle part of the cavity, where the temperature and concentration depart from its normal behavior, which can be called a reverse effect region. The temperature affected clearly by changing the 19 to 104 , 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

Double-diffusive natural convection with Cross-Diffusion Effects using ISPH method http://dx.doi.org/10.5772/60879 85

Fig. 5 Presents isothermal lines, concentration lines, horizontal and vertical velocities contours under the effect of Darcy number 10‐2 (**—**),10‐3(‐‐‐‐) and 10‐4(‐‐‐‐). **Figure 5.** Presents isothermal lines, concentration lines, horizontal and vertical velocities contours under the effect of Darcy number 10-2 (—),10-3(----) and 10-4(----).

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

**Figure 4.** Presents the effects of combined Soret and Dufour numbers on the temperature, concentration, horizontal

 The effects of Darcy number from 10-2 to 10-4 at inclination angles� � �� , permeability ratio �<sup>∗</sup> � ���, porosity parameter � � ���, Rayleigh number��� � ��� , buoyancy ratio � � ���, Soret number �� � ���, Dufour number �� � ���� and Lewis number �� � ��� 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 ��<sup>∗</sup> increases from 102

(c) Horizontal velocity profiles (d) Vertical velocity profiles Fig. 4 Presents the effects of combined Soret and Dufour numbers on the temperature, concentration, horizontal and vertical velocity profiles, respectively.


VY=0.5

Da=10-2

,

Le=1.0, Ra=104

(a) Temperature profiles (b) Concentration profiles


C

Da=10-2

,

Le=1.0, Ra=104

Sr=0.1 & Du=0.6 Sr=0.5 & Du=0.3 Sr=1.0 & Du=0.12 Sr=2.0 & Du=0.03

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=103 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

. 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 *X* increases. It is observed that, there is a region in the middle part of the cavity, where the temperature and concentration depart from its normal behavior, which can be called a reverse effect region. The temperature affected clearly by changing the

, the temperature and concentration contours become more parallel to the vertical walls and

19

increases from 102

0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0

X

, N=1.0,

0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0

Sr=0.1 & Du=0.6 Sr=0.5 & Du=0.3 Sr=1.0 & Du=0.12 Sr=2.0 & Du=0.03

Sr=0.1 & Du=0.6 Sr=0.5 & Du=0.3 Sr=1.0 & Du=0.12 Sr=2.0 & Du=0.03

X

, N=1.0,

to 104

, the

from10-4 to 10-2. In figure 5, as the inverse Darcy number *Da* \*

0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0


UX=0.5

Da=10-2

,

Le=1.0, Ra=104

, N=1.0,

X

, N=1.0,

104

to 104


> 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0

Y

T

Da=10-2

,

84 Mass Transfer - Advancement in Process Modelling

Sr=0.1 & Du=0.6 Sr=0.5 & Du=0.3 Sr=1.0 & Du=0.12 Sr=2.0 & Du=0.03

and vertical velocity profiles, respectively.

Le=1.0, Ra=104

Rayleigh number. On the other hand, concentration profiles are affected slightly by changing the Rayleigh number. The horizontal velocity curves are symmetric, where the line *Y* =0.5 represents the symmetrical axis and the variation of horizontal velocity increases as the Rayleigh number increases. Also, the vertical velocity curves are symmetric, where the line *X* =0.5 represents the symmetrical axis and also the variation of vertical velocity increases as the Rayleigh number increases. 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=103 and 10<sup>4</sup> . 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 ܺ increases. It is observed that, there is a region in the middle part of the

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 19

profiles, respectively at two values of Rayleigh number Ra=103 and 104. **Figure 6.** Depicts the effects of Darcy parameter on the temperature, concentration, horizontal and vertical velocity profiles, respectively at two values of Rayleigh number Ra=103 and 104.

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 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

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 effects of combined an inclination angle with permeability ratio parameter have been shown in Figs. 7 and 8.

The time histories for the effects of the Darcy number with two values of Rayleigh number *Ra* =10<sup>3</sup> and 10<sup>4</sup> and combined Soret with Dufour numbers on the average Nusselt number and the average Sherwood number are presented in Figs. 9 and 10. It is observed that, the average Nusselt number Nu has a decreasing trend with the presence of porous media by changing Darcy number from Da=10<sup>−</sup><sup>2</sup> to Da=10<sup>−</sup><sup>4</sup> . Due to the opposing flow of solute concentration, the average Sherwood number has an increasing trend with the presence of porous media by 20

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

**Figure 6.** Depicts the effects of Darcy parameter on the temperature, concentration, horizontal and vertical velocity

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

(c) Horizontal velocity profiles (d) Vertical velocity profiles Fig. 6 Depicts the effects of Darcy parameter on the temperature, concentration, horizontal and vertical velocity profiles, respectively at two values of Rayleigh number Ra=103 and 104.


VY=0.5

(a) Temperature profiles (b) Concentration profiles


C

Du=0.12,

Le=1.0, Sr=1.0

Du=0.12,

Le=1.0, Sr=1.0

0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0

Da=10-4 Ra=103

0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0

X

=0.1, N=1.0,

Da=10 Ra=10 -4

 Da=10-2 Da=10-3

> Da=10-2 Da=10-3

X

,104

=0.1, N=1.0,

3 ,104

 Da=10-2 Da=10-3 Da=10-4

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

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

the average Sherwood number are presented in Figs. 9 and 10. It is observed that, the average Nusselt number Nu has a decreasing trend with the presence of porous media by changing

average Sherwood number has an increasing trend with the presence of porous media by 20

and combined Soret with Dufour numbers on the average Nusselt number and

. Due to the opposing flow of solute concentration, the

been shown in Figs. 7 and 8.

*Ra* =10<sup>3</sup>

and 10<sup>4</sup>

in Figs. 7 and 8.


> 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0

UX=0.5

T

Du=0.12,

Ra=103 ,104

Le=1.0, Sr=1.0

Ra=103 ,104

 Da=10-2 Da=10-3 Da=10-4

Darcy number from Da=10<sup>−</sup><sup>2</sup>

of an inclination angle with permeability ratio parameter.

0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0

X


Y

profiles, respectively at two values of Rayleigh number Ra=103 and 104.

Du=0.12,

Le=1.0, Sr=1.0

=0.1, N=1.0,

=0.1, N=1.0,

86 Mass Transfer - Advancement in Process Modelling

to Da=10<sup>−</sup><sup>4</sup>

inclination angle for two values of permeability ratio �<sup>∗</sup> � ��� (**—**) and �<sup>∗</sup> � �� (‐‐‐‐). **Figure 7.** Presents isothermal lines, concentration lines, horizontal and vertical velocities contours under the effect of an inclination angle for two values of permeability ratio *K* \* =0.1 (—) and *K* \* =10 (----).

Fig. 7 Presents isothermal lines, concentration lines, horizontal and vertical velocities contours under the effect of an

changing Darcy number from Da=10<sup>−</sup><sup>2</sup> to Da=10<sup>−</sup><sup>4</sup> . Also, the effects of Darcy number on both of average Nusselt and Sherwood numbers appears clearly at the high value of Rayleigh number *Ra* =10<sup>4</sup> . 20

and Dufour number decreases.

inclination angle with permeability ratio parameter.

Fig. 7 Presents isothermal lines, concentration lines, horizontal and vertical velocities contours under the effect of an

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 for two values of permeability ratio �<sup>∗</sup> � ��� (**—**) and �<sup>∗</sup> � �� (‐‐‐‐).

20 profiles, respectively at two values of permeability ratio �<sup>∗</sup> � ��� and 10. **Figure 8.** Depicts the effects of an inclination angle on the temperature, concentration, horizontal and vertical velocity profiles, respectively at two values of permeability ratio *K* \* =0.1 and 10. number on both of average Nusselt and Sherwood numbers appears clearly at the high value of Rayleigh number �� � 10�.

Fig. 8 Depicts the effects of an inclination angle on the temperature, concentration, horizontal and vertical velocity

Sherwood numbers at the heated wall of the cavity. **Figure 9.** Time histories for the effects of Darcy number with two values of Rayleigh number on the average Nusselt and Sherwood numbers at the heated wall of the cavity.

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

Fig. 9 Time histories for the effects of Darcy number with two values of Rayleigh number on the average Nusselt and

26

number decreases. While, the average Sherwood number decreases as the Soret number increases and Dufour number decreases. Double-diffusive natural convection with Cross-Diffusion Effects using ISPH method http://dx.doi.org/10.5772/60879 89

Fig. 9 Time histories for the effects of Darcy number with two values of Rayleigh number on the average Nusselt and

In addition, the average Nusselt number increases as the Soret number increases and Dufour

Fig. 8 Depicts the effects of an inclination angle on the temperature, concentration, horizontal and vertical velocity profiles, respectively at two values of permeability ratio �<sup>∗</sup> � 0�1 and 10.

The time histories for the effects of the Darcy number with two values of Rayleigh number �� � 10� and 10� and combined Soret with Dufour numbers on the average Nusselt number and the average Sherwood number are presented in Figs. 9 and 10. It is observed that, the average Nusselt numberNu has a decreasing trend with the presence of porous media by changing Darcy number from �� � 10�� to �� � 10��. Due to the opposing flow of solute concentration, the average Sherwood number has an increasing trend with the presence of porous media by changing Darcy number from �� � 10�� to �� � 10��. Also, the effects of Darcy number on both of average Nusselt and Sherwood numbers appears clearly at the high value of

Rayleigh number �� � 10�.

Sherwood numbers at the heated wall of the cavity.

26 Fig. 10 Time histories for the effects of combined Soret number with Dufour number on the average Nusselt and **Figure 10.** Time histories for the effects of combined Soret number with Dufour number on the average Nusselt and Sherwood numbers at the heated wall of the cavity.

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**

20

(a) Average Nusselt number (b) Average Sherwood number Fig. 9 Time histories for the effects of Darcy number with two values of Rayleigh number on the average Nusselt and Sherwood numbers at the heated wall of the cavity.

**Figure 9.** Time histories for the effects of Darcy number with two values of Rayleigh number on the average Nusselt

0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0

Average Sherwood 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

26

**Figure 8.** Depicts the effects of an inclination angle on the temperature, concentration, horizontal and vertical velocity

(c) Horizontal velocity profiles (d) Vertical velocity profiles Fig. 8 Depicts the effects of an inclination angle on the temperature, concentration, horizontal and vertical velocity profiles, respectively at two values of permeability ratio �<sup>∗</sup> � ��� and 10.


=0.1 and 10.

Da=10-2

K\* =0.1, 10

Le=1.0, Sr=1.0Ra=104

Du=0.12,

Le=1.0, Sr=1.0

(a) Temperature profiles (b) Concentration profiles

Fig. 8 Depicts the effects of an inclination angle on the temperature, concentration, horizontal and vertical velocity profiles, respectively at two values of permeability ratio �<sup>∗</sup> � 0�1 and 10.

The time histories for the effects of the Darcy number with two values of Rayleigh number �� � 10� and 10� and combined Soret with Dufour numbers on the average Nusselt number and the average Sherwood number are presented in Figs. 9 and 10. It is observed that, the average Nusselt numberNu has a decreasing trend with the presence of porous media by changing Darcy number from �� � 10�� to �� � 10��. Due to the opposing flow of solute concentration, the average Sherwood number has an increasing trend with the presence of porous media by changing Darcy number from �� � 10�� to �� � 10��. Also, the effects of Darcy number on both of average Nusselt and Sherwood numbers appears clearly at the high value of

VY=0.5

C


Da=10-2

Le=1.0, Sr=1.0Ra=104

0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0

X

0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0

X

0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0 1.1 1.2

Ra=103 , 104

, Du=0.12, N=1.0,

=0.1, N=1.0,

, Du=0.12, N=1.0,

K\* =0.1, 10  

> 

 Da=10-2 Da=10-3 Da=10-4

Fig. 7 Presents isothermal lines, concentration lines, horizontal and vertical velocities contours under the effect of an

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 for two values of permeability ratio �<sup>∗</sup> � ��� (**—**) and �<sup>∗</sup> � �� (‐‐‐‐).

inclination angle with permeability ratio parameter.

0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0

X


profiles, respectively at two values of permeability ratio *K* \*

=0.1, N=1.0,

Ra=103 , 104

0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0 1.1 1.2

and Sherwood numbers at the heated wall of the cavity.

Y

Da=10-2

K\* =0.1, 10

Le=1.0, Sr=1.0Ra=104

, Du=0.12, N=1.0,

 Da=10-2 Da=10-3 Da=10-4

, Du=0.12, N=1.0,

K\* =0.1, 10

88 Mass Transfer - Advancement in Process Modelling


> 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0

0

1

2

Average Nusselt Number

3

4

5

UX=0.5

T

Da=10-2

Le=1.0, Sr=1.0Ra=104

 

Rayleigh number �� � 10�.

Du=0.12,

Le=1.0, Sr=1.0

and Dufour number decreases.

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 Da=10<sup>−</sup><sup>2</sup> , *ε* =0.5, Ra=10<sup>4</sup> , Sr=1.0, Du=0.12, *N* =1.0, Le=1.0 and *K* \* =0.1. In this figure, as the width ratio is increased, the oval shape becomes thinner and fluid activity at the top and bottom wall decreases. When the annulus length equals d= 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. Fig. 12 Depicts the time histories for the effects of annulus length on the average Nusselt number and average Sherwood number, respectively at Darcy parameter Da=10<sup>−</sup><sup>2</sup> , *ε* =0.5, Ra=10<sup>4</sup> , Sr=1.0, Du=0.12, *N* =1.0, Le=1.0 and *K* \* =0.1. 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.

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, Da=10<sup>−</sup><sup>2</sup> , *ε* =0.5, Ra=10<sup>4</sup> , *N* =1.0, Le=1.0 and *K* \* =0.1. In this figure, as the Soret number increases from 0.1 to 1.0 with decreasing the value of Dufour number from 0.6 to 0.12, the

velocity contours, respectively at Darcy parameter�� � ���� , � � ���, �� � ���, �� � ���, �� � ����,� � ��� , �� � ��� and �<sup>∗</sup> � ���. **Figure 11.** Depicts the effects of annulus length on the isotherms lines, concentration lines and vertical and horizontal velocity contours, respectively at Darcy parameter Da=10−<sup>2</sup> , *ε* =0.5, Ra=10<sup>4</sup> , Sr=1.0, Du=0.12, *N* =1.0, Le=1.0 and *K* \* =0.1.

Fig. 11 Depicts the effects of annulus length on the isotherms lines, concentration lines and vertical and horizontal

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

#### **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 28

Double-diffusive natural convection with Cross-Diffusion Effects using ISPH method http://dx.doi.org/10.5772/60879 91

Sherwood number, respectively at Darcy parameter �� � ����, � � ���, �� � ���, �� � ���, �� � ����, � � ���, �� � ��� and �<sup>∗</sup> � ���. **Figure 12.** Depicts the time histories for the effects of annulus length on the average Nusselt number and average Sher‐ wood number, respectively at Darcy parameter Da=10−<sup>2</sup> , *ε* =0.5, Ra=10<sup>4</sup> , Sr=1.0, Du=0.12, *N* =1.0, Le=1.0 and *K* \* =0.1.

Fig. 13 Depicts the time histories for the effects of combined Soret number with Dufour number

Fig. 12 Depicts the time histories for the effects of annulus length on the average Nusselt number and average

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. on the average Nusselt number and average Sherwood number, respectively at annulus length 0.5, �� � ����, � � ���, �� � ���, � � ���, �� � ��� and �<sup>∗</sup> � ���. In this figure, as the Soret number increases from 0.1 to 1.0 with decreasing the value of Dufour number from 0.6 to 0.12, the average Nusselt number increases. While, the average Sherwood number decreases as the Soret number increases with decreasing the Dufour number.

The following findings are summarized from the present investigation:

average Nusselt number increases. While, the average Sherwood number decreases as the

d=0.1 d=0.25 d=0.5 d=0.75 Fig. 11 Depicts the effects of annulus length on the isotherms lines, concentration lines and vertical and horizontal velocity contours, respectively at Darcy parameter�� � ���� , � � ���, �� � ���, �� � ���, �� � ����,� � ��� , �� � ��� and �<sup>∗</sup> � ���.

**Figure 11.** Depicts the effects of annulus length on the isotherms lines, concentration lines and vertical and horizontal

0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1



0.6


0.6

1.6

0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1

5

0.2

0.7

0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1

0.45

0.7

0.2

0.45

0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1

0.7



2.5

0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1

0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1

0.2

0.7

0.95

0.45

0.2

0.7

0.95

0.45

0.7

0.7

0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1

0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1

0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1

> 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1

velocity contours, respectively at Darcy parameter Da=10−<sup>2</sup>

5

0.95

0.7

0.45

0.2

0.7

0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1

0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 0.45

0.7

0.2

0.45

0.2

0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1

0.7

0.2

0.45

0.7

0.7

0.7

0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1

0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1

0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 0.2

0.2

0.45

0.45

0.45

0.45

0.2

0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1

0.2

0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1

0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1

, Sr=1.0, Du=0.12, *N* =1.0,

0.2

0.2

0.45

0.45

0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1

0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1

0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1

, *ε* =0.5, Ra=10<sup>4</sup>

0.45

0.7

0.2

0.7



0.6

1.6


0.6

0.45

0.45

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

28

Soret number increases with decreasing the Dufour number.


2.5


**7. Conclusion**

0.45

0.7

0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1

0.45

0.7

0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1


2

3

3


Le=1.0 and *K* \*

1.5

4

0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1

0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 -7




0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1

=0.1.


1.5

0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1

1.5

3



1.5

3

0.45

0.7

0.45

0.95

0.95

0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1

0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 0.7

0.95

90 Mass Transfer - Advancement in Process Modelling


Nusselt number and average Sherwood number, respectively at annulus length 0.5,�� � ���� , � � ���, �� � ���, � � ���, �� � ��� and �<sup>∗</sup> � ���. **Figure 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, Da=10−<sup>2</sup> , *ε* =0.5, Ra=10<sup>4</sup> , *N* =1.0, Le=1.0 and *K* \* =0.1.

Fig. 13 Depicts the time histories for the effects of combined Soret number with Dufour number on the average


Graphic results for the temperature, concentration, horizontal velocity and vertical velocity contours and representative velocity, temperature, and concentration profiles at the cavity


#### **Nomenclature** midsection for various parametric conditions were presented and discussed.

**CONCLUSION**


### **Acknowledgements**

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

### **Author details**

**•** 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.

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

Average Nusselt number Average Sherwood number 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,�� � ���� , � � ���, �� � ���, � � ���, �� � ��� and �<sup>∗</sup> � ���.

**Figure 13.** Depicts the time histories for the effects of combined Soret number with Dufour number on the average

0.0 0.5 1.0 1.5 2.0 2.5 3.0 3.5 4.0 4.5 5.0

Average Sherwood number

 Left Right Bottom Top

0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0

Sr&Du = 0.1&0.6, 1.0&0.12

Length (d=0.5)

, *ε* =0.5, Ra=10<sup>4</sup>

,

**•** 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

*C* ′ concentration of species *t* time *C* dimensionless species concentration *T* ′ temperature

midsection for various parametric conditions were presented and discussed. The following findings are summarized from the present investigation:

*Cp* specific heat *T* dimensionless temperature *Da* Darcy parameter *U* , *V* dimensionless velocity components

The permeability ratio parameter increases leads to decrease in the both of heat conduction

 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 heat and mass transfer and the flow characteristics inside the square cavity depended

*Du* Dufour number *V* velocity vector *d* annulus length *x*, *y* Cartesian coordinates *d*<sup>0</sup> particle size *X* , *Y* dimensionless coordinates

*F* Forchheimer coefficient **Greek symbols**

g gravitational acceleration vector *α* thermal diffusivity *K* permeability *β<sup>T</sup>* thermal expansion coefficient *k* thermal conductivity *β<sup>C</sup>* compositional expansion coefficient

29

 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

number increases with decreasing the Dufour number.

0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0

 Left Right Bottom Top

92 Mass Transfer - Advancement in Process Modelling

=0.1.

Sr&Du = 0.1&0.6, 1.0&0.12

Length (d=0.5)

Nusselt number and average Sherwood number, respectively at annulus length 0.5, Da=10−<sup>2</sup>

strongly on the Darcy number and Rayleigh number.

**Nomenclature**

and flow regime.

the rest of the walls.

**CONCLUSION**

*N* =1.0, Le=1.0 and *K* \*

AverageNusselt number

Abdelraheem M. Aly1,2\* and Mitsuteru Asai1

\*Address all correspondence to: abdelreheam.abdallah@sci.svu.edu.eg

1 Department of Civil Engineering, Kyushu University, 744 Motooka, Nishi-ku, Fukuoka, Japan

2 Department of Mathematics, Faculty of Science, South Valley University, Qena, Egypt

### **References**


radiation and Soret and Dufour effects. Journal of Naval Architecture and Marine Engineering, Vol 8, No 1, 2011.

[16] A. J. Chamkha and A. M. Aly, Heat and mass transfer in stagnation-point flow of a polar fluid towards a stretching surface in porous media in the presence of Soret, Du‐ four and chemical reaction effects. Chemical Engineering Communications, vol. 198, Issue 2, pp. 214-234, 2010.

[4] G. Juncu, "Unsteady conjugate forced convection heat/mass transfer in ensembles of‐ Newtonian fluid spheres," International Journal of Heat and Mass Transfer, vol. 53,

[5] J.-S. Leu, J.-Y. Jang, and Y. Chou, "Heat and mass transfer for liquid film evaporation along a vertical plate covered with a thin porous layer," International Journal of Heat

[6] W. Pirompugd, C.-C.Wang, and S.Wongwises, "Finite circular fin method for heat and mass transfer characteristics for plain fin-and-tube heat exchangers under fully and partially wet surface conditions," International Journal of Heat and Mass Trans‐

[7] W. Pirompugd, S. Wongwises, and C.-C. Wang, "Simultaneous heat and mass trans‐ fer characteristics for wavy fin-and-tube heat exchangers under dehumidifying con‐ ditions," International Journal of Heat and Mass Transfer, vol. 49,no. 1-2, pp. 132-

[8] M. Suresh and A. Mani, "Heat and mass transfer studies on R134a bubble absorber in R134a/DMF solution based on phenomenological theory," International Journal of

[9] P. Talukdar, C.R. Iskra, and C. J. Simonson, "Combined heat and mass transfer for laminar flow of moist air in a 3D rectangular duct: CFD simulation and validation with experimental data," International Journal of Heat and Mass Transfer, vol. 51,

[10] L.-Z. Zhang, "Heat and mass transfer in plate-fin sinusoidal passages with vaporpermeable wall materials," International Journal of Heat and Mass Transfer, vol. 51,

[11] L.-Z. Zhang, "Coupled heat and mass transfer in an application scale cross-flow hol‐ low fiber membrane module for air humidification," International Journal of Heat

[12] F.-Y. Zhao, D. Liu, and G.-F. Tang, "Application issues of the streamline, heatline and massline for conjugate heat and mass transfer," International Journal of Heat

[13] A. Mahdy, "MHD non-Darcian free convection from a vertical wavy surface embed‐ ded in porous media in the presence of Soret and Dufour effect," International Com‐

[14] M. A. Mansour, N. F. El-Anssary and A. M. Aly. Effects of Chemical Reaction and Thermal Stratification on MHD Free Convective Heat and Mass Transfer over a Ver‐ tical Stretching Surface Embedded in Porous Media Considering Soret and Dufour

[15] A. J. Chamkha, M. F. Al-Amin, A. M. Aly, Unsteady double-diffusive natural convec‐ tive MHD flow along a vertical cylinder in the presence of chemical reaction, thermal

munications in Heat and Mass Transfer, vol. 36, no. 10, pp. 1067-1074, 2009.

Numbers. Chemical Engineering Journal, vol. 145, pp. 340-345, 2008.

and Mass Transfer, vol. 49, no. 11-12, pp. 1937-1945, 2006.

Heat and Mass Transfer, vol. 53, no. 13-14, pp. 2813-2825, 2010.

and Mass Transfer, vol. 55, no. 21-22, pp. 5861-5869, 2012.

and Mass Transfer, vol. 50, no. 1-2, pp. 320-334, 2007.

no. 13-14, pp. 2780- 2789, 2010.

94 Mass Transfer - Advancement in Process Modelling

fer, vol. 50, no. 3-4, pp. 552-565, 2007.

no. 11- 12, pp. 3091-3102, 2008.

no. 3-4, pp. 618-629, 2008.

143, 2006.


[42] A. M. Aly, Modeling of multi-phase flows and natural convection in a square cavity using an Incompressible Smoothed Particle Hydrodynamics, International Journal of Numerical Methods for Heat & Fluid Flow, vol. 25 (3), 2015.

[29] G. Degan, P. Vasseur, Boundary-layer regime in a vertical porous layer with aniso‐ tropic permeability and boundary effects, International Journal of Heat and Fluid

[30] P. Nithiarasu, K.N. Seetharamu, T. Sundararajan, Natural convective heat transfer in a fluid saturated variable porosity medium, International Journal of Heat and Mass

[31] P. Nithiarasu, K. S. Sujatha, K. Ravindran, T. Sundararajan, K. N. Seetharamu, Non-Darcy natural convection in a hydrodynamically and thermally anisotropic porous medium, Computer Methods in Applied Mechanics and Engineering, Vol. 188, pp.

[32] J. J. Monaghan. Simulating free surface flows with SPH. Journal of Computational

[33] J. P. Morris, P. J. Fox, and Y. Zhu. Modeling low Reynolds number incompressible flows using SPH. Journal of Computational Physics, 136(1):214 - 226, 1997.

[34] S. J. Cummins and M. Rudman. An SPH projection method. Journal of Computation‐

[35] M. Asai, A. M. Aly, Y. Sonoda and Y. Sakai. A Stabilized Incompressible SPH meth‐ od by relaxing the Density invariance condition. Journal of Applied Mathematics,

[36] A. M. Aly, M. Asai, Y. Sonoda. Modelling of surface tension force for free surface flows in ISPH method. International Journal of Numerical Methods for Heat & Fluid

[37] A. M. Aly, M. Asai, Y. Sonoda. Simulation of free falling rigid body into water by a stabilized incompressible SPH method. Ocean Systems Engineering, An Internation‐

[38] A. M. Aly, An Improved Incompressible Smoothed Particle Hydrodynamics to Simu‐

[39] A. Chaniotis, D. Poulikakosand P. Koumoutsakos. Remeshed smoothed particle hy‐ drodynamics for the simulationof viscous and heat conducting flows. Journal of

[40] K. Szewc, J. Pozorski, and A. Tani`ere. Modeling of natural convection with smooth‐ ed particle hydrodynamics: Non-Boussinesq formulation. International Journal of

[41] M. E. Danis, M. Orhan and A. Ecder, ISPH modeling of transient natural convection, International Journal of Computational Fluid Dynamics, 27:1, pp. 15-31, 2013.

late Fluid-Soil-Structure Interactions, Ph.D., Kyushu University, 2012.

Flow, Vol. 18, (3), pp. 334-343, 1997.

96 Mass Transfer - Advancement in Process Modelling

Physics, 110(2):399 - 406, 1994.

al Physics, 152(2):584 - 607, 1999.

2012; 2012: 24. doi:10.1155/2012/139583.

Flow, vol. 23 Iss: 3, pp.479-498, 2013.

al Journal, vol. 1(3): pp. 207-222, 2011.

ComputationalPhysics,182, no. 1: pp. 67-90, 2002.

Heat and Mass Transfer, vol. 54, no. 23-24: pp. 4807-16, 2011.

413-430, 2000.

Transfer, Vol. 40, (16), pp. 3955-3967, 1997.

