Open access peer-reviewed chapter

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

Written By

Abdelraheem M. Aly and Mitsuteru Asai

Submitted: 31 October 2014 Reviewed: 22 May 2015 Published: 22 October 2015

DOI: 10.5772/60879

From the Edited Volume

Mass Transfer - Advancement in Process Modelling

Edited by Marek Solecki

Chapter metrics overview

1,763 Chapter Downloads

View Full Metrics

Abstract

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

Keywords

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

1. Introduction

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

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

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

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

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

Advertisement

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.

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.

Advertisement

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:

.V=0E1
Vt+(Vε)V ==ερ[P+με(2Vx2+2Vy2)μK¯VE¯ρ V|V|+ρ g(βT(TTl')+βC(CCl'))]E2
σTt+VT=1ρCP((k¯T))+DTC((C)),E3
εCt+VC=((D¯C))+DCT((T)),E4

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, k¯, and mass diffusivity, D¯, can be written in the following coordinate system as:

K¯=[K1cos2θ+K2sin2θ(K1K2)sinθ cosθ(K1K2)sinθcosθK2cos2θ+K1sin2θ],E5
K*=K1K2,E6
E¯=[E1cos2θ+E2sin2θ(E1E2)sinθ cosθ(E1E2)sinθ cosθE2cos2θ+E1sin2θ],E7
E1=F1K1,E2=F2K2,E8
F*=F1F2, E9
E*=E1E2=F*K*E10

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

K1¯=1K1[cos2θ+K*sin2θ(1K*)sinθcosθ(1K*)sinθ cosθK*cos2θ+sin2θ],E11

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

E¯=F1K1[cos2θ+K*F*sin2θ(1K*F*)sinθ cosθ(1K*F*)sinθ cosθK*F*cos2θ+sin2θ],E12

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

k¯=[kx00ky],E13
k*=kykx,E14

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

k¯=kx[100k*],E15

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

D¯=[Dx00Dy],E16
D*=DyDx,E17

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

D¯=Dx[100D*],E18

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

C=c'cl'Δc',N=βc(ch'cl')βT(Th'Tl'),Le=axDxE19

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

UX+VY=0,E20
1Uτ+U2UX+V2VY=PX+Pr[2UX2+2UY2]PrDa[U(cos2θ+K*sin2θ)+((1K*)sinθcosθ)]F1Da[U(cos2θ+K*F*sin2θ)+V((1K*F*)sinθcosθ)]U2+V2,E21
1Vτ+U2VX+V2VY=PY+Pr[2VX2+2VY2]PrDa[U((1K*)sinθcosθ)+V(K*cos2θ+sin2θ)]F1Da[U((1K*F*)sinθcosθ)+V((1K*F*)cos2θsin2θ)]U2+V2+RaPr(T+NC),E22
σTτ+UTX+VTY=2TX2+k*2TY2+Du(2CX2+2CY2),E23
εCτ+UCX+VCY=1Le(2CX2+D*2CY2)+Sr(2TX2+2TY2)E24

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

3.1. Boundary treatment

The dimensionless boundary conditions are:

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

U=0;  V=0;  T=1, C=1 for   X=0,   0Y1U=0;  V=0;  T=0, C=0 for   X=1,   0Y1U=0;  V=0;  TX=0,  CX=0 for    Y=0,   0X1U=0;  V=0; TX=0,  CX=0 for    Y=1,   0X1E25
  1. Second case (square annulus filled with anisotropic porous media)

U=0;  V=0;  T=1, C=1                    for   X=0,   0Y1U=0;  V=0;  T=1, C=1                    for   X=1,   0Y1U=0;  V=0;  T=1, C=1                    for    Y=0,   0X1U=0;  V=0; T=1, C=1                     for    Y=1,   0X1U=0;  V=0; T=1, C=1                     for    Y=1,   0X1U=0;  V=0;  T=0, C=0                    for(Ld)/2Y(L+d)/2,   X=(L+d)/2,U=0;  V=0;  T=0, C=0                    for(Ld)/2X(L+d)/2,   Y=(Ld)/2,U=0;  V=0;  T=0, C=0                    for(Ld)/2X(L+d)/2,   Y=(L+d)/2,E26

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

Nu= TX|X=0,E27
Sh= CX|X=0.E28

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:

Nu¯= 1L0L(TX)dY.E29
Sh¯= 1L0L(CX)dYE30

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

Horizontal hot walls:

Nu= TY|T=Th, Sh= CY|T=Th,E31

Vertical hot walls:

Nu= TX|T=Th, Sh= CX|T=Th,E32

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:

NuL¯= 1L0L(TX)dY, ShL¯= 1L0L(CX)dY, E33

Average Nusselt and Sherwood numbers at the right hot surface:

NuR¯= 1L0L(TX)dY, ShR¯= 1L0L(CX)dY, E34

Average Nusselt and Sherwood numbers at the top hot surface:

NuT¯= 1L0L(TY)dX, ShR¯= 1L0L(CY)dX, E35

Average Nusselt and Sherwood numbers at the bottom hot surface:

NuB¯= 1L0L(TY)dX, ShB¯= 1L0L(CY)dX, E36

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

NuTotal¯= 14(NuL¯+NuR¯+NuT¯+NuB¯),E37
ShTotal¯= 14(ShL¯+ShR¯+ShT¯+ShB¯),E38
Advertisement

4. Numerical method

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

4.1. Finite volume method

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

i,j|χi,jnewχi,jold|106,E39

4.2. ISPH method

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

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

4.2.1. Temporal discretization

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

1[Un+1UnΔt]=(PX)n+1+Pr(2UX2+2UY2)n(PrDa(cos2θ+K*sin2θ)+F1Da(cos2θ+K*F*sin2θ)U2+V2)Un+1((PrDa(1K*)sinθcosθ)+F1Da((1K*F*)sinθcosθ)U2+V2)VnE40
1[Vn+1VnΔt]=(PY)n+1+Pr(2UX2+2UY2)n(PrDa((1K*)sinθcosθ)+F1Da(1K*F*sin2θcosθ)U2+V2)Un(PrDa(K*cos2θ+sin2θ)+F1Da(K*F*cos2θ+sin2θ)U2+V2)Vn+1++PrRa(Tn+NCn)E41

For simplicity,

Porsy=PrDa(K*cos2θ+sin2θ)+F1Da(K*F*cos2θ+sin2θ)(U2+V2)n,Porsr=(PrDa(1K*)sinθcosθ)+F1Da((1K*F*)sinθcosθ)(U2+V2)n,E42

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

(1+ϵΔt Porsx)Un+1(1+ϵΔt Porsx)U*+(1+ϵΔt Porsx)U*UnϵΔt==(PX)n+1+Prϵ (2UX2+2UY2)nPorsr Vn,E43
(1+ϵΔt Porsy)Vn+1(1+ϵΔt Porsy)V*+(1+ϵΔt Porsy)V*VnϵΔt==(PY)n+1+Prϵ (2VX2+2VY2)nPorsr Un+PrRa (Tn+N Cn)E44

Step 1.Predict the velocity

The first step of the predictor-corrector scheme is the calculation of an intermediate velocity (U*) from the momentum equation without including the pressure terms. Thus, the following is obtained:

(1+ϵΔt Porsx)U*=Un+ϵΔt(Prϵ (2UX2+2UY2)nPorsr Vn),E45
(1+ϵΔt Porsy)V*= Vn+ϵΔt(Prϵ (2VX2+2VY2)n+PrRa (Tn+N Cn)Porsr Un) .E46

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:

2Pn+1X2+2Pn+1Y2=((1+ϵΔt Porsx)ϵΔt(U*X)+(1+ϵΔt Porsy)ϵΔt(V*Y))+γρ0<ρ*>ρ0Δt2E47

Note, the relaxation coefficient, γ,(0γ1) can be decided from pre-analysis calculation as Asai et al. [35]. In this study, the particle size d0 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:

(1+ϵΔt Porsx)Un+1=(1+ϵΔt Porsx)U*ϵΔt[PX]n+1,E48
(1+ϵΔt Porsy)Vn+1=(1+ϵΔt Porsy)V*ϵΔt[PY]n+1,E49

Step 4.Thermal flow problems:

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

σ(Tn+1TnΔτ)=((2TX2+k*2TY2))n+Du(2CX2+2CY2)n,E50

Step 5.Concentration flow problems:

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

ε(Cn+1CnΔτ)=1Le((2CX2+D*2CY2))n+Sr(2TX2+2TY2)n.E51

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:

ϕ(Xi,τ)=W(XiXj,h)ϕ(Xj,τ)dv=W(Rij,h)ϕ(Xj,τ)dv,E52

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

ϕ(Xi,τ)ϕi=jmjρjW(Rij,h)ϕ(Xj,τ),E53

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

ϕ(Xi)ϕi=1ρijmj(ϕjϕi)W(Rij,h).E54

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

ϕ(Xi)ϕi=ρijmj(ϕjρj2+ϕiρi2)W(Rij,h).E55

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

W(Rij,h)=βd[(3Rijh)56(2Rijh)5+15(1Rijh)50Rij<h(3Rijh)56(2Rijh)5hRij<2h(3Rijh)52hRij<3h0rij3h,E56

where βd is 7/478πh2 and3/358πh3, respectively, in two- and three dimension space.

4.2.3. Projection-based ISPH formulations

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

P(Xi)Pi=ρijmj(Pjρj2+Piρi2)W(Rij,h),E57
.U(Xi).Ui=1ρijmj(UjUi).W(Rij,h),E58

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:

.(υU)(Xi).(υ.Ui)=jmj(ρiυi+ρjυjρiρjRijW(|RiRj|,h)Rij2+η2)Uij,E59

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

.(υ.Ui)=2υiρijmj(RijW(|RiRj|,h)Rij2+η2)Uij,E60

In this study, Laplacian of velocity is given by:

2U(Xi)2Ui=2ρijmj(RijW(|RiRj|,h)Rij2+η2)Uij,E61

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

2P(Xi)2Pi=2ρijmj(PijRijW(|RiRj|,h)Rij2+η2).E62

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:

2T(Xi)2Ti=jmj(ρi+ρjρiρjTijRijW(|RiRj|,h)Rij2+η2),E63
2C(Xi)2Ci=jmj(ρi+ρjρiρjCijRijW(|RiRj|,h)Rij2+η2),E64
Advertisement

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=102 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=102. 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].

θ K* Ra Nithiarasu et al. [31] FVM ISPH
0.1 104 1.587 1.570 1.517
45° 0.1 104 1.573 1.563 1.522
90° 0.1 104 1.579 1.573 1.524
10 104 1.106 1.110 1.132
45° 10 104 1.119 1.130 1.139
90° 10 104 1.106 1.111 1.134

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=102 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=102, porosity parameter ε=0.5, buoyancy ratio N=0 and an inclination angle of permeability θ=0° is introduced in figure 2. In this figure, the isotherm lines are affected clearly by both of the Rayleigh number and permeability ratio. As the Rayleigh number increases, the heat conduction increases and the isotherms lines are strongly compressed in the bottom part of the left wall and in the top part of the right wall. As the permeability ratio increases from 0.1 to 10, then the heat conduction decreases.

Figure 2.

Comparison of the isotherms lines between FVM method (—) and ISPH method (----) at buoyancy ratio N=0, Darcy number Da=102, porosity parameter ε=0.5 and inclination angle θ=0° in the absence of Soret and Dufour numbers effects.

Advertisement

6. Results and discussion

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

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

In this study, the unsteady natural convection in square cavity saturated with anisotropic porous media under the effects of Soret and Dufour numbers has been solved numerically using ISPH method. It is assumed that, the vertical walls are kept at constant temperatures, 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 permeability ratio K*=0.1, Darcy parameter Da=102, porosity parameter ε=0.5, Rayleigh number  Ra=104, buoyancy ratio N=1.0 and Lewis number Le=1.0 has been shown in figure 3.

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

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=102, ε=0.5, Ra=104, N=1.0 and Le=1.0.

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

Figure 4.

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

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=104, buoyancy ratio N=1.0, Soret number Sr=1.0, Dufour number Du=0.12 and Lewis number Le=1.0 has been shown in figure 5. For the value of the buoyancy ratio N considered, the interaction between the thermal and compositional buoyancy effects is significant. It is observed that, the heat and mass conduction in cavity are increase clearly and the flow moves faster at the case of increasing Darcy number from10-4 to 10-2. In figure 5, as the inverse Darcy number Da* increases from 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.

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(----).

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

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.

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

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 (----).

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.

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.

The time histories for the effects of the Darcy number with two values of Rayleigh number Ra=103 and 104 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=102 to Da=104. 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 Da=102 to Da=104. Also, the effects of Darcy number on both of average Nusselt and Sherwood numbers appears clearly at the high value of Rayleigh number Ra=104.

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 and Dufour number decreases.

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.

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

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

The effects of annulus length in the temperature, concentration, average Nusselt number and average Sherwood number have been introduced in Figs. 11 and 12. Fig. 11 Depicts the effects of annulus length on the isotherms lines, concentration lines and vertical and horizontal velocity contours, respectively at Darcy parameter Da=102, ε=0.5, Ra=104, 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=102, ε=0.5, Ra=104, 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.

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=102, ε=0.5, Ra=104, Sr=1.0, Du=0.12, N=1.0, Le=1.0 and K*=0.1.

Figure 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=102, ε=0.5, Ra=104, 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 on the average Nusselt number and average Sherwood number, respectively at annulus length 0.5, Da=102, ε=0.5, Ra=104, 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 average Nusselt number increases. While, the average Sherwood number decreases as the Soret number increases with decreasing the Dufour number.

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=102, ε=0.5, Ra=104, N=1.0, Le=1.0 and K*=0.1.

Advertisement

7. Conclusion

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

The following findings are summarized from the present investigation:

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

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

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

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

Advertisement

Nomenclature

C concentration of species t time
C dimensionless species concentration T temperature
Cp specific heat T dimensionless temperature
Da Darcy parameter U, V dimensionless velocity components
Du Dufour number V velocity vector
d annulus length x, y Cartesian coordinates
d0 particle size X, Y dimensionless coordinates
F Forchheimer coefficient Greek symbols
g gravitational acceleration vector α thermal diffusivity
K permeability βT thermal expansion coefficient
k thermal conductivity βC compositional expansion coefficient
Nu Nusselt number ε porosity
N Buoyancy ratio μ viscosity
Le Lewis number υ kinematic viscosity
P pressure σ ratio of heat capacities
Pr Prandtl number ρ density
Ra Rayleigh number τ dimensionless time
Sh Sherwood number 2 Laplacian operator
Sr Soret number

Advertisement

Acknowledgments

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

References

  1. 1. H. Demir, M. Mobedi, and S. Ülkü, “Effects of porosity on heat and mass transfer in a granular adsorbent bed,” International Communications in Heat and Mass Transfer, vol. 36, no. 4, pp. 372-377, 2009.
  2. 2. S. Han and R. J.Goldstein, “The heat/mass transfer analogy for a simulated turbine blade,” International Journal of Heat and Mass Transfer, vol. 51, no. 21-22, pp. 5209-5225, 2008.
  3. 3. S. Han and R. J. Goldstein, “The heat/mass transfer analogy for a simulated turbine end wall,” International Journal of Heat and Mass Transfer, vol. 51, no. 11-12, pp. 3227-3244, 2008.
  4. 4. G. Juncu, “Unsteady conjugate forced convection heat/mass transfer in ensembles ofNewtonian fluid spheres,” International Journal of Heat and Mass Transfer, vol. 53, no. 13-14, pp. 2780- 2789, 2010.
  5. 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 and Mass Transfer, vol. 49, no. 11-12, pp. 1937-1945, 2006.
  6. 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 Transfer, vol. 50, no. 3-4, pp. 552-565, 2007.
  7. 7. W. Pirompugd, S. Wongwises, and C.-C. Wang, “Simultaneous heat and mass transfer characteristics for wavy fin-and-tube heat exchangers under dehumidifying conditions,” International Journal of Heat and Mass Transfer, vol. 49,no. 1-2, pp. 132- 143, 2006.
  8. 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 Heat and Mass Transfer, vol. 53, no. 13-14, pp. 2813-2825, 2010.
  9. 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, no. 11- 12, pp. 3091-3102, 2008.
  10. 10. L.-Z. Zhang, “Heat and mass transfer in plate-fin sinusoidal passages with vapor-permeable wall materials,” International Journal of Heat and Mass Transfer, vol. 51, no. 3-4, pp. 618-629, 2008.
  11. 11. L.-Z. Zhang, “Coupled heat and mass transfer in an application scale cross-flow hollow fiber membrane module for air humidification,” International Journal of Heat and Mass Transfer, vol. 55, no. 21-22, pp. 5861-5869, 2012.
  12. 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 and Mass Transfer, vol. 50, no. 1-2, pp. 320-334, 2007.
  13. 13. A. Mahdy, “MHD non-Darcian free convection from a vertical wavy surface embedded in porous media in the presence of Soret and Dufour effect,” International Communications in Heat and Mass Transfer, vol. 36, no. 10, pp. 1067-1074, 2009.
  14. 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 Vertical Stretching Surface Embedded in Porous Media Considering Soret and Dufour Numbers. Chemical Engineering Journal, vol. 145, pp. 340-345, 2008.
  15. 15. A. J. Chamkha, M. F. Al-Amin, A. M. Aly, Unsteady double-diffusive natural convective MHD flow along a vertical cylinder in the presence of chemical reaction, thermal radiation and Soret and Dufour effects. Journal of Naval Architecture and Marine Engineering, Vol 8, No 1, 2011.
  16. 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, Dufour and chemical reaction effects. Chemical Engineering Communications, vol. 198, Issue 2, pp. 214-234, 2010.
  17. 17. A. M. Aly, M. A. Mansour, Ali J. Chamkha. Effects of soret and Dufour numbers on free convection over isothermal and adiabatic stretching surfaces embedded in porous media, Journal of Porous Media, vol. 14, pp. 67-72, 2011.
  18. 18. N. Nithyadevi and R.-J. Yang, “Double diffusive natural convection in a partially heated enclosure with Soret and Dufour effects,” International Journal of Heat and Fluid Flow, vol. 30, no. 5, pp. 902-910, 2009.
  19. 19. J. A. Weaver and R. Viskanta, “Natural convection due to horizontal temperature and concentration gradients—2. Species inter diffusion, Soret and Dufour effects,” International Journal of Heat and Mass Transfer, vol. 34, no. 12, pp. 3121-3133, 1991.
  20. 20. S. K. Sinha, T. Sundararajan, V.K. Garg, A variable property analysis of alloy solidification using the anisotropic porous medium approach, Int. J. Heat Mass Transfer Vol. 35, pp. 2865-2877, 1992.
  21. 21. A. R. Chaudhuri, K.N. Seetharamu, T. Sundararajan, Modelling of steam surface condenser using finite element method, Comm. Numer. Meth. Engrg. Vol. 13, pp. 909-921, 1997.
  22. 22. P. Cheng, Heat transfer in geothermal systems, Adu. Heat Transfer. vol. 14, pp. l-105, 1978.
  23. 23. V. Prasad, F. A. Kulacki and M. Keyhani, Natural convection in porous media, J. Fluid Mech. vol. 150, pp. 89- 119, 1985.
  24. 24. C. L. Tien and J. T. Hong, Natural convection in porous media under non-Darcian and non-uniform permeability conditions. In Natural Convection (Edited by S. Kakac, W. Aung and R. Viskanta). Hemisphere, Washington, DC, 1985.
  25. 25. P. Cheng, Wall effects on fluid flow and heat transfer in porous media. Proc. ASME/JSME Heat Transfer Conf., pp. 297-303, 1987.
  26. 26. V. Prasad, G. Lauriat and N. Kladias, Reexamination of Darcy-Brinkman solutions for free convection in porous media, Int. Symp. Convection in Porous Media: Non- Darcy Effects, Proc. 25th Nat. Heat Transfer Conf. vol. 1, pp. 569-580, 1988.
  27. 27. N. Kladias and V. Prasad, Benard convection in porous media: effects of Darcy and Prandtl numbers, Int. Symp. Convection in Porous Media: Non-Darcy Effects, Proc. 25th Nat. Heat Transfer Conf. vol. 1, pp. 593-604, 1988.
  28. 28. D. A. Nield and A. Bejan, Convection in Porous Media. Springer, New York, 1992.
  29. 29. G. Degan, P. Vasseur, Boundary-layer regime in a vertical porous layer with anisotropic permeability and boundary effects, International Journal of Heat and Fluid Flow, Vol. 18, (3), pp. 334-343, 1997.
  30. 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 Transfer, Vol. 40, (16), pp. 3955-3967, 1997.
  31. 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. 413-430, 2000.
  32. 32. J. J. Monaghan. Simulating free surface flows with SPH. Journal of Computational Physics, 110(2):399 - 406, 1994.
  33. 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. 34. S. J. Cummins and M. Rudman. An SPH projection method. Journal of Computational Physics, 152(2):584 - 607, 1999.
  35. 35. M. Asai, A. M. Aly, Y. Sonoda and Y. Sakai. A Stabilized Incompressible SPH method by relaxing the Density invariance condition. Journal of Applied Mathematics, 2012; 2012: 24. doi:10.1155/2012/139583.
  36. 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 Flow, vol. 23 Iss: 3, pp.479-498, 2013.
  37. 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 International Journal, vol. 1(3): pp. 207-222, 2011.
  38. 38. A. M. Aly, An Improved Incompressible Smoothed Particle Hydrodynamics to Simulate Fluid-Soil-Structure Interactions, Ph.D., Kyushu University, 2012.
  39. 39. A. Chaniotis, D. Poulikakosand P. Koumoutsakos. Remeshed smoothed particle hydrodynamics for the simulationof viscous and heat conducting flows. Journal of ComputationalPhysics,182, no. 1: pp. 67-90, 2002.
  40. 40. K. Szewc, J. Pozorski, and A. Tani`ere. Modeling of natural convection with smoothed particle hydrodynamics: Non-Boussinesq formulation. International Journal of Heat and Mass Transfer, vol. 54, no. 23-24: pp. 4807-16, 2011.
  41. 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.
  42. 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.
  43. 43. A. M. Aly, Mitsuteru Asai, Modelling of non-Darcy Flows through porous media using an extended Incompressible Smoothed Particle Hydrodynamics, Numerical Heat Transfer, Part B: Fundamentals, 2014, doi = 10.1080/10407790.2014.955772.
  44. 44. A. M. Aly and S. E. Ahmed, An incompressible smoothed particle hydrodynamics method for natural/mixed convection in a non-Darcy anisotropic porous medium, Int. J. Heat Mass Transfer, vol. 77, pp. 1155-1168, 2014.
  45. 45. S. V. Patankar, Numerical Heat Transfer and Fluid Flow, Hemisphere, McGraw-Hill, Washington DC, 1980.
  46. 46. J. A. Meijerink and H. A. van der Vorst, An iterative solution method for linear systems of which the coefficient matrix is a symmetric M-matrix, l Mathematics of Computation, vol. 31, no. 137, pp. 148-162, 1977.

Written By

Abdelraheem M. Aly and Mitsuteru Asai

Submitted: 31 October 2014 Reviewed: 22 May 2015 Published: 22 October 2015