CFD modeling of humidification dehumidification distillation process

CFD modeling of humidification dehumidification distillation process

Desalination 395 (2016) 46–56 Contents lists available at ScienceDirect Desalination journal homepage: www.elsevier.com/locate/desal CFD modeling o...

4MB Sizes 1 Downloads 74 Views

Desalination 395 (2016) 46–56

Contents lists available at ScienceDirect

Desalination journal homepage: www.elsevier.com/locate/desal

CFD modeling of humidification dehumidification distillation process Adnan Saeed a, Mohamed A. Antar a,⁎, Mostafa H. Sharqawy b, Hassan M. Badr a a b

Department of Mechanical Engineering, King Fahd University of Petroleum and Minerals, Dhahran 31261, Saudi Arabia School of Engineering, University of Guelph, Guelph, Ontario N1G 2W1, Canada

H I G H L I G H T S • • • •

Natural convection humidification dehumidification process is modeled in a rectangular enclosure. Maximum distillate rate occurs at an aspect ratio of 1.5 at Rayleigh number = 104 and buoyancy ratio = 1. Nusselt number increases with high convection rates and high buoyancy ratio. Nusselt number based on latent heat exceeds the one based on sensible heating.

a r t i c l e

i n f o

Article history: Received 26 October 2015 Received in revised form 16 March 2016 Accepted 17 March 2016 Available online 31 May 2016 Keywords: Water distillation Humidification cycle Dehumidification Natural convection Mass transfer

a b s t r a c t Exploring new methodologies to meet the growing demands for water and energy sources is an important step towards a sustainable future. HDH is a promising technology applied for both water desalination and industrial water treatment. However, due to the lower Gain Output Ratio (GOR) of HDH systems, its wide-ranging applications are restricted. This issue can be resolved solely by increasing the GOR of the HDH systems. The main principle of HDH systems considered in this study is the distillation of water under atmospheric condition by an air loop saturated with water vapor and driven by natural convection thus eliminating the need for pumping power. In this work, a numerical investigation is carried out for predicting the performance of a simple HDH cycle taking place in a rectangular enclosure (cavity) of different aspect ratios. The cavity has two horizontal adiabatic walls and two vertical isothermal walls; one of them is heated while the other one is cooled. The vertical hot wall is kept wet by maintaining a thin layer of water flowing downstream and the cavity is assumed long enough to justify the assumption of twodimensional flow. A computational model is developed for predicting the velocity, temperature and concentration fields within the cavity as well as calculating the rate of water evaporating from the cavity hot side and condensing on its cold side. The model is based on numerical solution of the conservation equations of mass, momentum, energy and species considering the Boussinesq approximation. The model has been validated and the effect of aspect ratio on the performance of the HDH system is investigated. The results show that the aspect ratio of 1.5 produces the maximum heat and mass transfers. In addition, high convection rates tend to maximize the Nusselt number. The values of Nusselt based on latent heat exceed the Nusselt based on sensible heating. © 2016 Elsevier B.V. All rights reserved.

Nomenclature A aspect ratio C concentration/mass fraction (kg vapor/kg mixture) D diffusion coefficient g gravitational acceleration H height k thermal conductivity L length Le Lewis number Mv molecular mass of water vapor (kg/kmol) Ma molecular mass of air (kg/kmol) N buoyancy ratio

⁎ Corresponding author.

http://dx.doi.org/10.1016/j.desal.2016.03.011 0011-9164/© 2016 Elsevier B.V. All rights reserved.

Nu P Pr Rax X Th Tc

Nusselt number atmospheric pressure Prandtl number Rayleigh number at position x characteristic length hot surface temperature quiescent temperature (fluid temperature far from the surface)

Greek symbols

α ν ρ

thermal diffusivity kinematic viscosity density (kg/m3)

A. Saeed et al. / Desalination 395 (2016) 46–56

β ψ

thermal expansion coefficient stream function

1. Introduction Desalination is a process of removing salts from saline water. 97% of global water is in oceans and the percentage of clean water available for human beings is only 0.014% [1]. Natural resources of fresh water are mitigating abruptly and are not enough for growing industrial and population needs for supplies of clean drinking water. Diminishing water resources forced human race to search new sources. Besides, the quick reduction of aquifers caused water scarcity issues in different places of the world especially in areas lacking renewable resources. Conventional desalination technologies are well established such that they operate and produce large fresh water rates of about 2.5– 10 MIGD (Million Imperial Gallons per day). Few other technologies are under research and development phases. They are expected to lead new approaches to be applied to design systems in order to enhance the effectives of technologies intended for sea water desalination. Humidification–dehumidification (HDH) is a desalination process in which air is used as a working fluid and a carrier gas to obtain desalinated water. In this process, raw water is first evaporated in dry air (known as air humidification) and the generated vapor is then condensed to obtain clean water in an air dehumidifier. The energy required to evaporate and condense sea water can be obtained from solar, thermal, geothermal, and combinations of renewable energy sources. Air is driven in this system by forced or natural convection. If it is driven by natural convection, the cycle is referred as “humidification cycle”. To understand this cycle, natural convection heat transfer in various geometries has been the subject of numerous experimental and computational studies. Because of temperature and concentration differences in the gravity field, natural convection streams are created by density gradients. Most common studies encompass heat transfer in enclosures and over spheres and cylinders. Ostrach [2] solved the boundary layer equation for free convection commencing at a flat upright plate by means of numerical techniques. He reduced the governing equations (continuity, momentum and energy) into two equations with their particular boundary condition. He concluded that such flow relies on Prandtl and Grashof numbers. Kumar and Dalal [3] used vorticity and stream function formulations of Navier-Stokes equations to model the free convection around a heated cylinder placed in an enclosure. They used finite-difference approximations for solving the governing equations numerically. The effect of geometry on flow field and heat transfer has been simulated for 3 aspect ratios. The authors reported that overall heat transfer and flow pattern change with changing aspect ratios. It was found that the constant wall temperature heating mode is more efficient than the constant wall heat flux mode for higher heat transfer. Free convection heat transfer and fluid flow in a porous triangular cavity having 2 conducting solid baffles have been presented recently by Mushatet [4]. Darcy model is used for mathematical relations of mass, momentum and energy. Different parameters such as baffle location, Rayleigh number, baffle height and baffle width are investigated to explain flow and thermal field characteristics. The author reported that increased deviations of streamlines and a multi-cellular flow have been obtained for Ra ≥ 150. They reported that conducing baffles have a significant effect on the flow and heat transfer. In most of the double diffusive natural convection cases, opposite surfaces of the cavity or enclosure are subjected to temperature and concentration gradients. Wee et al. [5] analyzed heat and moisture transfers in a cavity both experimentally and numerically. They considered both vertical and horizontal cavities. They studied simultaneous transfer of moisture and heat by free convection in a cavity having air as the working fluid and an aspect ratio of 7. The values of the Pr and Schmidt numbers used were 0.7 and 0.6, respectively. Both opposing and aiding flows were computed for different ranges of Grashof number using the finite-difference approximation. The rate of moisture and heat

47

transfers in the vertical and horizontal cavities with upward transfer was found to be five to nine times higher than those for downward flow in the horizontal cavity. A similar numerical study is conducted by Han and Kuehn [6] on double diffusive natural convection in a vertical rectangular cavity of an aspect ratio of 4, the concentration gradient as well as temperature were imposed in a horizontal direction using a finite difference algorithm. Various types of flow structure regimes were obtained for aiding and opposing buoyancy conditions as functions of both thermal and solutal Grashof numbers. Lee and Hyun [7] studied combined horizontal temperature and concentration gradients in a rectangular cavity having an aspect ratio of 2. In this case, thermal and solutal buoyancy forces are counteracting. The thermal and solutal velocity boundary layer structures are examined near to the side wall. They reported that in the beginning thermal buoyancy in much of the cavity is dominant as a result of the larger difference of the diffusivities of solute and heat in a large Lewis number. In most cases, the flow behavior strongly depended on the buoyancy ratio. An experimental study by Kamotani et al. [8] on natural convection in rectangular enclosures of low aspect-ratios was conducted with joint horizontal concentration and temperature gradients. The solutal buoyancy force either augments or opposes the thermal buoyancy force. The flow possessed double-diffusive characteristics as a result of a large difference between the solutal and thermal diffusion rates. The solutal boundary layers beside the vertical walls were much thinner than the thermal boundary layer. In the case of the buoyancy ratio below 6, flow had a unicellular pattern with secondary cells close to the perpendicular walls. A numerical study of free convection considering double diffusion in a partial water heated field with Soret and Dufour coefficients was presented by Nithyadevi and Yang [9]. They showed that when the values of thermal Rayleigh number increase, the rate of heat and mass transfers increases too. Another numerical study was conducted by Kuznetsov and Sheremet [10] on double-diffusive conjugate natural convection cubical (three-dimensional) enclosure with finite wall thicknesses filled with air. The study showed an increase of the heat and mass transfer intensity at the Rayleigh number range 104 ≤ Ra ≤ 105. Mezhrab et al. [11] studied double diffusive natural convection governed by the combined effects of radiation and natural convection in a gray gas square cavity. It was concluded that the effect of radiation depends on the flow regime (opposite or assisting), buoyancy ratio ‘N’ and the Rayleigh number ‘Ra’. Recently in 2013, a solar distiller having a working principal based on double diffusive natural convection is presented by Ghachem et al. [12]. The conservation equations were solved for the laminar flow case using the finite-volume technique. Rayleigh number is considered invariant and of the order of 105. For opposed temperature and concentration gradients, the effect of the buoyancy ratio was studied for different aspect ratios and the entropy generation was studied in details. It was concluded that the isotherms, isoconcentration lines and the structure of the flow are significantly affected by buoyancy ratio. Heat and mass irreversibilities for N = 1 were found to have high Bejan number. Complex distribution of Nusselt number was observed for buoyancy ratio of 1. Recent research in double diffusive natural convection in enclosures includes the analytical treatment of Nield and Kuznetsov [13] in a porous medium saturated by a nanofluid. They used Darcy model and presented a similarity solution. They have only considered aiding buoyancy forces. Kifayati [14,15] used finite difference lattice Boltzmann method to study double diffusive natural convection with Soret and Dufour effects. The cavity is filled with power-law non-Newtonian fluid with entropy generation. They reported that while Dufour effect increased heat and mass transfers, Soret parameter only augments mass transfer and entropy generation. To the best of the author's knowledge and from the literature review cited above, evaporation and condensation in conjunction with double diffusive natural convection have not received sufficient attention. In this paper, the phenomenon of heat and mass transfers is applied to

48

A. Saeed et al. / Desalination 395 (2016) 46–56

humidification dehumidification technique mainly used for the desalination of water. A study based on double diffusive natural convection with evaporation and condensation can simulate humidification cycle, i.e. an HDH system based on natural convection. In the present work, an attempt is made to study the effect of aspect ratio, buoyancy ratio and Raleigh number on heat and mass transfers in an enclosure intended for HDH desalination system. 2. Problem formulation 2.1. Basic idea Muller Holst et al. [16] presented an idea of multi-effect HDH in order to enhance the heat recovery in their experimental work. They optimized the desalination method for producing small scale (0.5– 2 m3/day) production for standalone operation of desalination system that is supposed to be a maintenance-free system that can operate either on solar energy or diesel engine waste heat. They proposed that the multi-effect humidification–dehumidification desalination system works on atmospheric pressure for the evaporation and condensation of water in a thermally insulated box. Air is used to transport vapor from the evaporation portion to the condenser portion where the distillate is then collected. It was concluded that higher heat recovery could be achieved with optimized size of HDH system. This process was further improved to allow the circulation of air to be carried out entirely by natural convection. The reason for the highest GOR of this HDH system lies in a fact that less energy is needed in the humidifier due to heat recovery from dehumidifier that results in lowering the energy requirement. Fig. 1 shows a schematic of the system used in an experiment to study multi-effect HDH desalination. 2.2. Physical description Müller-Holst [16,17] showed experimentally that a high performance HDH system can be achieved based on flat evaporation surface of plastic sheets in an enclosed space (insulated box, or a cavity for example). Therefore, the evaporator-condenser chamber is an enclosure

Fig. 2. Schematic of the enclosure.

with temperature differences between the evaporator and condenser walls which develop the natural convection loops inside the enclosure and drive the evaporation and condensation processes. Fig. 2 shows a simplified sketch of a rectangular enclosure with a length—L and a height—H, which simulates an evaporation–condensation cell of HDH system. The left side wall has higher temperature as well as concentration compared to the right wall. Top and bottom walls are considered thermally adiabatic. The boundary conditions are selected to obtain opposite thermal and solutal buoyancies. The left side wall is wetted by a thin liquid water film and maintained at constant temperature. It is considered that liquid film is extremely thin and can be treated as a boundary condition [18]. Vaporization occurs at the hot surface from the liquid vapor interface. Combined heat as well as mass transfer resulting from temperature and concentration gradients is numerically investigated in this enclosure. The combined heat and mass transfers are dependent on fluid properties as well as the geometry of the enclosure whose aspect ratio is given by (H/L). 2.3. Mathematical modeling Laminar and steady two-dimensional governing equations inside an enclosure are considered (Fig. 2) for natural convection heat and mass transfers of water vapor. The model is based on the governing equations of mass, momentum and energy in addition to the concentration equation. The water on the left wall of the enclosure, due to the contact with high temperature wall, is vaporized from the liquid–vapor interface. The vapor is carried to the opposite side to condense at the colder side wall.

Fig. 1. Schematic of multi-effect HDH process [16].

Fig. 3. Variations of the average Nusselt number on the left wall for A = 2.

A. Saeed et al. / Desalination 395 (2016) 46–56

49

Motion of fluid mixture, air-water vapor inside the enclosure can be described by the following conservation equations [19,20]: Continuity equation ∂u ∂v þ ¼ 0: ∂x ∂y

ð1Þ

Navier Stokes equations !   2 2 ∂u ∂u ∂P ∂ u ∂ u þ ¼− þμ ρ u þv ∂x ∂y ∂x ∂x2 ∂y2 !   2 2 ∂v ∂v ∂P ∂ v ∂ v þ ρ u þv ¼− þμ ∂x ∂y ∂y ∂x2 ∂y2 þ ρ∘ gðβT ðT−T ∘ Þ þ βM ðC−C ∘ ÞÞ:

ð2Þ

ð3Þ Fig. 5. Validation of vertical velocity profile validation at enclosure mid-section with Chamkha [27].

Energy equation !   2 2 ∂T ∂T ∂ T ∂ T þ : þv ¼K ρC p u ∂x ∂y ∂x2 ∂y2

ð4Þ

Species transport equation !   2 2 ∂C ∂C ∂ C ∂ C þ þv ¼D : u ∂x ∂y ∂x2 ∂y2

ð5Þ

The fluid is considered to be Newtonian and the flow is incompressible. Due to the low Mach number, viscous dissipation is considered negligible. Considering constant properties (the Lewis number is taken as 1 and Prandtl number is taken as 0.7) except for the density in the body force term, Boussinesq approximation is applied while neglecting radiation, heat generation, and Soret as well as Dufour effects. The Soret effect (concentration gradients driving heat flux) and Dufour effect (temperature gradients driving mass flux) are neglected because of their smaller order of magnitude compared to the effects produced by Fick's and Fourier's laws. At water–air interface, thermodynamic equilibrium is assumed and constant thermophysical properties are assumed except for the density in the body force term, that is considered as linear function of temperature and mass fractions as described by Boussinesq approximation [21]; ρðT; C Þ ¼ ρ∘ ½1−βT ðT−T C Þ−½1−βM ðC−C C Þ:

ð6Þ

The coefficients of thermal and mass fraction expressions, are defined by,

βC ¼

  −1 ∂ρ : ρ ∂C P;C

ð8Þ

The governing transport equations for two-dimensional steady laminar incompressible flows written in Cartesian coordinates describing double diffusive natural convection are normalized using appropriate parameters and then solved numerically as a function of vorticity, stream function, temperature and concentration. The coupled equations are descritized using finite difference approximations. The Newton–Raphson method is used in solving these equations. The parameters used to put the governing equation in a dimensionless form are given as: X¼

x y Y ðNon‐dimensional lengthsÞ L L Ω

ωL2 α

T−T c T h −T c

C

Ψ¼ θ¼

ψ α

ðStream function and vorticityÞ c−cc ch −cc

ðTemperature and concentrationÞ :

Whereas vorticity and stream functions are defined as: ω ¼ΔV ¼

∂v ∂u dψ ¼u − ∂y ∂y dy

dψ ¼v dx

ð7Þ

where L is the characteristic length of the cavity, βC & β T are the coefficients of the concentration and thermal expansions respectively. Th, Tc, Ch, and Cc are the hot and cold wall temperatures and concentrations respectively. α is thermal diffusivity of the fluid.

Fig. 4. Validation of the temperature profile results at enclosure mid-section with Barakos and Mitsoulis [26].

Fig. 6. Validation of horizontal velocity profile at enclosure mid-section with Chamkha [27].

βT ¼

  −1 ∂ρ : ρ ∂T P;C

50

A. Saeed et al. / Desalination 395 (2016) 46–56

Momentum equation !     2 2 ∂Ψ ∂Ω ∂Ψ ∂Ω ∂ Ω ∂ Ω ∂θ ∂C þ RaPr þ − ¼ Pr þ N : ∂Y ∂X ∂X ∂Y ∂X ∂X ∂X 2 ∂Y 2

ð10Þ

Energy equation   2 2 ∂Ψ ∂θ ∂Ψ ∂θ ∂ θ ∂ θ − ¼ 2þ 2: ∂Y ∂X ∂X ∂Y ∂X ∂Y

Fig. 7. Validation of temperature and concentration profiles at enclosure mid-section with the results of Chamkha [26].

ð11Þ

Species Transport equation !   2 2 ∂Ψ ∂C ∂Ψ ∂C 1 ∂ C ∂ C þ − ¼ : Le ∂X 2 ∂Y 2 ∂Y ∂X ∂X ∂Y

ð12Þ

The above quantities are used to define the non-dimensional variables: Ra = Rayleigh number, Le = Lewis number, Pr = Prandtl number and N = buoyancy ratio. 2.4. Boundary conditions ν υ Pr ¼ ; Sc ¼ ; α D β ðc −cc Þ N¼ C h : βT ðT h −T c Þ

α Le ¼ ; D

gβðT h −T c Þ 3 L ; RaL ¼ να

hL Nu ¼ ; k

The model boundary conditions in the dimension less form are given below: 2

The vorticity and stream function formulations are used to simplify the equations by omitting the pressure terms from the momentum equations. The two-dimensional equations for continuity, momentum energy and species transfer considering steady natural convection can be written as: Continuity equation 2

∂ Ψ ∂X 2

X ¼ 0; 0 ≤Y ≤1

:

∂Ψ ∂Ψ ∂ Ψ ¼0 ; ¼ 0 ;Ω ¼ − 2 ; θ ¼ C ¼ 1 ∂X ∂Y ∂X

X ¼ 1; 0 ≤Y ≤1

:

∂Ψ ∂Ψ ∂ Ψ ¼0 ; ¼ 0 ;Ω ¼ − 2 ; θ ¼ C ¼ 0 ∂X ∂Y ∂X

2

2

2

þ

∂ Ψ ∂Y 2

¼ −Ω:

ð9Þ

Y ¼ 0 or Y ¼ 1; 0≤X ≤1 :

∂Ψ ∂Ψ ∂ Ψ ∂θ ∂C ¼ ¼ 0; Ω ¼ − 2 ; ¼ ¼ 0: ∂X ∂Y ∂Y ∂Y ∂Y

Table 1 Effect of variation in buoyancy ratio and Raleigh number on stream lines for A = 1. N=0 Ra = 10

4

Ra = 105

Ra = 106

Stream lines for aspect ratio = 1

N=1

N=2

A. Saeed et al. / Desalination 395 (2016) 46–56

51

Table 2 The effect of changing buoyancy ratio and Raleigh number on the isotherms for aspect ratio A = 1. N=0

N=1

N=2

Ra = 104

Ra = 105

Ra = 106

Isotherms for aspect ratio = 1

The evaporated liquid film into the air stream or condensed humid air vapors on the wet plate is represented by Velocity ‘Vw’ that can be evaluated by [19,20]: In non-dimensional form, Vw ¼

 PrðC H −C C Þ ∂C  : Scð1−C w Þ ∂X X¼0

ð13Þ

As per Dalton's law and the state equation of ideal gas mixture, the interfacial mass fraction of water vapor on the wetted wall can be calculated as: Mv Ma Cw ¼ Mv P −1 þ Ma P sat ðT w Þ

ð14Þ

2.5. Flow, heat and mass transfer parameters Heat and mass transfers between the wetted walls and the fluid motion within the enclosure depend on two factors. (i) Fluid temperature gradient at the wetted walls, resulting in sensible heat transfer. (ii) The mass transfer rate, resulting in latent heat transfer. Hence the total heat flux from wetted wall can be expressed as [23]: qx ¼

 ∂T  qs þ ql ¼ −k  − ∂X X¼0

 ρDhfg ∂C  : 1−C w ∂X X¼0

ð16Þ

Therefore, the average Nusselt number along the wetted surface is defined as Nux ¼

Nus þ Nul

ð17Þ

where Psat is the saturated water vapor pressure on the wetted wall. The following equation is used to calculate the saturation pressure [22]:

where Nus and Nul, are the local Nusselt numbers for sensible and latent heat transfer respectively.

log10 P sat ðT Þ ¼ 28:59051−8:2 logT þ 2:4804  10−3 T−3142:32=T:

Nus ¼ ð15Þ

 −1 ∂θ   1−θb ∂X 

X¼0

ð18Þ

52

A. Saeed et al. / Desalination 395 (2016) 46–56

Table 3 Effect of variation in buoyancy ratio and Raleigh number on stream lines for A = 2. N=0

N=1

N=2

Ra = 104

Ra = 105

Ra = 106

Stream lines for aspect ratio = 2

Nul ¼

 −S 1 ∂C  :  1−θb 1−C w ∂X 

ð19Þ

X¼0

where S indicates the importance of energy transfer through species diffusion relative to thermal diffusion [24]. It is given by

3.1. Grid independence

ρDhfg ðC H −C c Þ : S¼ kðT H −T C Þ

ð20Þ

The dimensionless bulk temperature is defined as:

θb ¼

1 um

0:5 Z

ð21Þ

u:θdX: 0

Similar to mass transfer rate at the interface, the locally averaged Sherwood number on the wetted wall can be calculated as: Sh ¼

 −1 ∂C   ð1−C b Þ ∂X 

:

ð22Þ

X¼0

Whereas,

Cb ¼

1 um

obtaining the temperature and concentration distributions within the flow field. This process is repeated until satisfying the convergence criteria. The stream function (Ψ) and the average Nusselt and Sherwood numbers are calculated. The values for mass transfer are computed and velocity, temperature and concentration profiles are recorded.

0:5 Z

u:CdX:

ð23Þ

0

3. Solution procedure The solution starts by specifying values for the following parameters: Rayleigh number (Ra), Prandtl numbers (Pr), Lewis number (Le), buoyancy ratio (N), diffusion coefficient (D), thermal diffusivity (α), the enclosure hot and cold wall temperatures and the aspect ratio (A). The vorticity and stream function equations are similar to those presented by Rahman et al. [25]. They solved to obtain ζ and ψ distributions, which are then substituted in energy and species equations for

The grid sensitivity tests are conducted for the cases of heat transfer in an enclosure in order to save computational time and to obtain better accuracy. Different uniform mesh sizes were employed ranging from 11∗11 to 99∗99 and Nusselt number was calculated. The Nusselt number is unaffected by increasing the mesh size beyond 51∗51, so this grid was deemed adequate for the Rayleigh number in consideration as shown in Fig. 3. 3.2. Model validation To validate the developed numerical model, it is simplified and reduced for the application to 2 basic cases that are presented in the literature. First case is the problem of heat transfer in differentially heated enclosure and the second case is the problem of double diffusive natural convection in a square cavity. 3.2.1. Case 1: differentially heated enclosure model validation The heat transfer in a differentially heated enclosure is modeled such that the momentum and mass conservation equations along with the energy equation are solved. Comparisons of the results with Barakos and Mitsoulis [26] show a very good agreement as illustrated in Fig. 4. The model is validated for differentially heated enclosure with the results reported earlier for Ra = 1.89 × 105. After this validation, the work is extended and species transport equation was added in the model in case 2. 3.2.2. Case 2: double diffusive natural convection model validation The case is validated with an earlier work on double diffusive natural convection in a rectangular enclosure that was reported by Chamkha [27]. Comparisons of results are illustrated in Fig. 5–7. The comparison

A. Saeed et al. / Desalination 395 (2016) 46–56

53

Fig. 10. Variation of horizontal velocity (U) with buoyancy ratio (N) for A = 1 at the middle section along x axis (X = 0.5). Fig. 8. Variation of vertical velocity (V) with buoyancy ratio (N) for A = 1 at the middle section X = 0.5.

results showed a very good agreement with double diffusive natural convection case for heat and mass transfers. The flow is driven by both temperature and concentration gradients in this case. When Darcy number is set equal to zero, the resulting equations become facsimile to our current work. In order to validate the complete flow regime, both velocity components are validated and results for velocity profiles at enclosure midsection are compared. The profile for dimensionless temperature and concentration coincides since the value for Lewis number, Le = 1, considering our working fluid as a moist air. The concentration and energy equations become similar and their profiles exactly lie on top of each other. 4. Results and discussion In this section, the results are presented for fluid flow, heat and mass transfers for different values of the non-dimensional governing parameters, including the Rayleigh number (104 ≤ Ra ≤ 106) and the buoyancy ratio (0 ≤ N ≤ 3) for various aspect ratios of the enclosure (1 ≤ A ≤ 2). Since the flow, thermal and concentration fields are dependent on the geometrical aspects of the enclosure; we will present the stream lines and isotherms for aspect ratios of 1 and 2. Results presented in Table 1 show the effect of buoyancy ratio on streamlines for N = 0, 1 and 2 for Raleigh number of 104,105 and 106 respectively. Buoyancy ratio N describes the relative effect of species diffusion on thermal diffusion that creates in the density gradient driving the

Fig. 9. Variation of vertical velocity (V) with buoyancy ratio (N) for A = 2 at the middle section X = 0.5.

flow. It is worth to note that when N = 0, the temperature difference results solely in producing buoyancy force that cause the flow and there is no mass diffusion effect. When N N 0, the combined buoyancy force resulting from mass diffusion and thermal diffusion assists the flow. As the buoyancy ratio is increased from N = 0 to N = 3, the streamlines change from circular flow to asymmetric flow but the temperature gradient pattern remains unchanged (Table 2). As the Rayleigh number is increased from 104 to 106, velocity gradients become high near both (right and left) walls, this is due to the fact that temperature gradients near right and left walls are increased, which enhances natural convection. The temperature gradients (with respect to y-axis) near wall are zero, since the top and bottom walls are insulated and there is no heat transfer from these walls as evident from the. A boundary layer at the left hot wall of the enclosure is noticed at Ra = 105. This boundary layer becomes thinner and more visible with the increase of Rayleigh number. The streamlines exhibit a clockwise vortex that is generated at the center. This cell is solely generated due to temperature difference for N = 0. With the increase of Rayleigh number, for Ra = 104, elliptic-shaped velocity contours are seen at the center of the enclosure and convection effects are more prominent in the isotherms. The temperature gradients are higher at the walls and diminish at the center. At Ra = 105, velocity contours at the center are stretched and as a result, two small vortices/cells arise in the main cell. In the boundary layer, convection heat transfer distorts the isotherms and temperature gradients become equal to zero. These vortices do not appear for the case of Ra number less than 105 due to the presence of viscous diffusion force. The straightened isotherms prevent the motion in the normal direction at the center of the enclosure. The dependence of natural convection on Rayleigh number can be better presented from the isothermal contours. When the Rayleigh number is low, parallel isotherms are observed close to the left hot

Fig. 11. Variation of horizontal velocity (U) with buoyancy ratio (N) for A = 2 at the middle section along x axis X = 0.5.

54

A. Saeed et al. / Desalination 395 (2016) 46–56

Fig. 12. Variation of average Nusselt number with buoyancy ratio and Rayleigh number.

Fig. 14. Variation of vertical velocity (V) with Rayleigh number for A = 1 at the middle section along x axis (X = 0.5).

wall. When the isotherms are straight and close to each other, it indicates that heat conduction dominates. The convection phenomenon prevails when the isotherms deform and take sharp turns. For high Rayleigh numbers the isotherm deformation becomes more significant since convection heat transfer prevails at high Rayleigh numbers. Changing aspect ratio from A = 1 to A = 2, similar trends are observed for the temperature gradients, which means that the trends would be similar for the concentration contours. Increasing Rayleigh number from 104 to 105 has no significant effect on these contours for conduction heat transfer close to the walls. However, increasing the aspect ratio would cause eddy formation in the middle section of the domain as evident from the streamline contours in Table 3. The variation of non-dimensional velocity ‘V’ with respect to nondimensional length is shown in Fig. 8 for various buoyancy ratios for A = 1 and Ra = 105. At the enclosure mid-section, it is clear that the magnitude of maximum velocity increases with the increase of buoyancy ratio and maximum velocity points are closer to the walls. The maximum non-dimensional velocity increases from a value of 60 for N = 0 to 130 for N = 3 for the case of A = 1. Fig. 9 represents the vertical velocity component profiles for A = 2 and Ra = 105 at the enclosure mid-section for different buoyancy ratios. It is clear that in the enclosure, the magnitude of maximum velocity increases with the increase of buoyancy ratio and maximum velocity points get closer to the walls due to high gradients and the velocity difference between any two maximum points on the graph mitigates. It is also evident from both Figs. 8 and 9 that at N = 0, for the case of pure natural convection due to thermal gradient only, the vertical velocity values increase when the buoyancy force generated by mass diffusion

adds to the total force for N = 1, 2 and 3.The highest vertical velocity is achieved at N = 3. As the aspect ratio changes from A = 1 to A = 2, the maximum nondimensional velocity increase from approximately 130 to 180 for the case of N = 3, this is due to the fact that the heated wall length changes as the aspect ratio changes (increases convection area) that enhances convection driving forces. Fig. 10 represents the variation of horizontal velocity component with buoyancy ratio for A = 1 and Ra = 105.The flow makes a clockwise rotation in the enclosure and the highest value of U is observed at H = 0.85 and at H = 0.15 with flow in the opposite direction. The fluid is stagnant at the middle of the enclosure and it takes the same path for the high values of buoyancy ratios. Fig. 11 represents the variation of horizontal velocity component with buoyancy ratio for A = 2 and Ra = 105. Unlike A = 1, a secondary vortex is also generated at the enclosure mid-section that also rotates in the same fashion for the case A = 2. Fig. 12 represents the average Nusselt number variation with Rayleigh number and buoyancy ratio. It is found that Nu increases with the increase of Rayleigh number and buoyancy ratio. Since flow is subjected to additional force due to solutal buoyancy force for N greater than 1, hence high convection rates tend to maximize the Nu number on the wall. Fig. 13 represents the variation of average latent Nu with Rayleigh number and buoyancy ratio. Since the flow is unique where both sensible and latent heat transfers are subjected to the same fluid rotating due to natural convection in an enclosure. The trends observed for latent

Fig. 13. Variation of average latent Nusselt number with buoyancy ratio and Rayleigh number.

Fig. 15. Variation of vertical velocity (V) with aspect ratio for same Ra number and buoyancy ratio at the middle section along x axis (X = 0.5).

A. Saeed et al. / Desalination 395 (2016) 46–56

55

Fig. 16. Variation of horizontal velocity (U) with Rayleigh number for A = 1 at the middle section along x axis X = 0.5. Fig. 18. Average Nusselt number variation with enclosure aspect ratio.

heat transfer are similar to those for sensible heat transfer. Values of Nusselt number based on latent heat transfer are higher than the sensible heat values since heat transfer results in phase change that incorporates higher energy transfer. Fig. 14 depicts the vertical velocity component profiles for A = 1 for varying Rayleigh numbers from 104 to 106 at the enclosure mid-section for constant buoyancy ratio value of 1. The figure shows that in an enclosure, the magnitude of maximum velocity increases with the increase of Rayleigh number and maximum velocity points get closer to the walls due to higher gradients. The change in the value of the vertical velocities is due to the increase in convection flow in an enclosure due to increasing Rayleigh number. Similarly variation of vertical velocity with aspect ratio is shown in Fig. 15. For the same Rayleigh number and buoyancy ratio, vertical velocity profile for A = 2 is having a higher value than that for A = 1 at the center of the enclosure. When the enclosure height is doubled (A = 2), retardation to the flow generated is due to the combined effect of thermal and solutal buoyancy forces. It is less experienced by the horizontal walls than that for the case of A = 1. Continuous heat is supplied to the fluid near the hot wall for bigger height and hence due to high density variations, high convection rate leads to higher vertical velocities. Fig. 16 represents the horizontal velocity component profiles for A = 1 with varying Rayleigh number that varies from 104 to 106 at the enclosure mid-vertical section for constant buoyancy ratio value of 1. The magnitude of maximum velocity increases with the increase of Rayleigh number. The increase in the value of the horizontal velocities is due to the increase in convection flow in an enclosure with increasing Rayleigh number as shown in Fig. 17. In Fig. 18, variation of average sensible Nusselt number with respect to buoyancy ratio is shown for different aspect ratios. The average Nusselt number increases with the buoyancy ratio, the Nusselt number. is maximum for A = 1.5 for any value of buoyancy ratio.

Fig. 17. Variation of horizontal velocity (U) with buoyancy ratio (N) for A = 2 at the middle section along x axis X = 0.5.

Fig. 19 shows the effect of enclosure aspect ratio on the average distillate rate per day. It is clear from the figure that the maximum rate is obtained when A = 1.5. In connection to Fig. 18, the dependence of heat and mass transfers on geometry has significant dependency on the enclosure aspect ratios in terms of the output of the system. 5. Conclusion The performance of humidification dehumidification desalination process taking place in a simple rectangular enclosure is numerically investigated. The investigation is based on the solution of the conservation equations of mass, momentum and energy in addition to the species transport equation. The cavity, which has two horizontal adiabatic walls and two vertical isothermal walls is assumed to have different aspect ratios. The effects of Rayleigh number (104 ≤ Ra ≤ 106) and buoyancy ratio (0 ≤ N ≤ 3) are investigated for the different values of the cavity aspect ratios (1 ≤ A ≤ 2). The obtained results indicate a maximum rate of average distillate occurring at an aspect ratio of 1.5 for the case of Rayleigh number, Ra = 104 and buoyancy ratio, N = 1. The results have also shown the dependence of heat and mass transfers on the enclosure aspect ratios in terms of output distilled water. Significant effect on heat and mass transfers is noted for buoyancy ratio, aspect ratio and the Rayleigh number. The Nusselt number is observed to increase with high convection rates and high buoyancy ratio, its component associated with latent heat is higher than the sensible heat component. Acknowledgements The authors would like to thank King Fahd University of Petroleum and Minerals in for funding the research reported in this work through

Fig. 19. Variation of computed distillation rate with cavity aspect ratio.

56

A. Saeed et al. / Desalination 395 (2016) 46–56

National Science, Technology and Innovation Plan office, NSTIP Project # 11-WAT1625-04. The authors would also like to thank Mr. Waqas Akram for his helpful comments during the development of this work. References [1] J. Bendfeld, C. Broker, K. Menne, E. Ortjohann, L. Temme, J. Vob, P.C.M. Carvallo, Design of a PV-powered reverse osmosis plant for desalination of brackish water, Proceedings of 2nd World Conference and Exhibition on Photovoltaic Solar Energy Conversion 1998, pp. 3075–3077. [2] Ostrach, An analysis of laminar Free-Convection Heat Transfer about a flat plate parallel to the direction of the generating body force, Supersedes NACA TN 2635, 1952 1111. [3] De A. Kumar, A. Dalal, A numerical study of natural convection around a square, horizontal, heated cylinder placed in an enclosure, Int. J. Heat Mass Transf. 49 (2006) 4608–4623. [4] K.S. Mushatet, Natural convection heat transfer inside a porous triangular enclosure with baffles, Global J. of Res. in Eng. 10 (1) (2010) 63–75. [5] H.K. Wee, R.B. Keey, M.J. Cunningham, Heat and moisture transfer by natural convection in a rectangular cavity, Int. J. Heat Mass Transf. 32 (9) (1989) 1765–1778. [6] H. Han, T.H. Kuehn, Double diffusive natural convection in a vertical rectangular enclosure- II. Numerical study, Int. J. Heat Mass Transf. 34 (2) (1991) 461–471. [7] J.W. Lee, J.M. Hyun, Double-diffusive convection in a rectangle with opposing horizontal temperature and concentration gradients, Int. J. Heat Mass Transf. 33 (8) (1990) 1619–1632. [8] Y. Kamotani, L.W. Wang, S. Ostrach, H.D. Jiang, Experimental study of natural convection in shallow enclosures with horizontal temperature and concentration gradients, Int. J. Heat Mass Transf. 28 (1) (1985) 165–173. [9] N. Nithyadevi, R. Yang, Double diffusive natural convection in a partially heated enclosure with Soret and Dufour effects, Int. J. Heat Fluid Flow 30 (5) (2009) 902–910. [10] G.V. Kuznetsov, M.A. Sheremet, A numerical simulation of double-diffusive conjugate natural convection in an enclosure, Int. J. Therm. Sci. 50 (10) (2011) 1878–1886. [11] A. Mezrhab, D. Lemonnier, S. Meftah, A. Benbrik, Numerical study of doublediffusion convection coupled to radiation in a square cavity filled with a participating grey gas, J. Phys. D. Appl. Phys. 41 (19) (2008) 1–16. [12] K. Ghachem, L. Kolsi, C. Mâatki, A.K. Hussain, M.N. Borjini, Numerical simulation of three-dimensional double diffusive free convection flow and irreversibility studies in a solar distiller, Int. Commun. in Heat Mass Transfer 39 (6) (2012) 869–876.

[13] D.A.A. Nield, A.V. Kuznetsov, The cheng-minkowycz problem for the double diffusive natural convection boundary layer flow in a porous medium saturated by a nanofluid, Int. J. Heat Mass Transf. 54 (1–3) (2016) 374–378. [14] G.H.R. Kifayati, FDLBM Simulation Of Entropy Generation In Double Diffusive Natural Convection of Power-Law Fluids in an Enclosure with Soret and Dufour Effects, Int. J. Heat Mass Transfer 89 (2015) 267–290. [15] G.H.R. Kefayati, Simulation of double diffusive natural convection and entropy generation of power-law fluids in an inclined porous cavity with Soret and Dufour effects (Part II: Entropy generation), Int. J. Heat Mass Transfer 94 (2016) 582–624. [16] H. Muller-Holst, M. Engelhardt, W. Scholkopf, Small-scale thermal seawater desalination simulation and optimization of system design, Desalination 122 (2) (1999) 255–262. [17] H. Müller-Holst, M. Engelhardt, M. Herve, W. Schölkopf, Solarthermal seawater desalination systems for decentralised use, Renew. Energy 14 (1998) 311–318. [18] H. Tsai, W. Yan, Mixed convection heat and mass transfer in vertical rectangular ducts, Int. J. Heat Mass Transf. 40 (7) (1997) 1621–1631. [19] L.C. Burmeister, Convective Heat Transfer, McGraw-Hill, New York, 1983. [20] R.B. Bird, W.E. Stewart, E.N. Lightfood, Transport Phenomena, Wiley, New York, 2002. [21] J. Boussinesq, Théorie analytique de la chaleur, Gauthier-Villars, Paris, France, 1903. [22] M.A. Kassim, B. Benhamou, S. Harmand, Effect of air humidity at the entrance on heat and mass transfers in a humidifier intended for a desalination system, Appl. Therm. Eng. 31 (11–12) (2011) 1906–1914. [23] J.H. Jang, W.M. Yan, C.C. Huang, Mixed convection heat transfer enhancement through film evaporation in inclined square ducts, Int. J. Heat Mass Transf. 48 (11) (2005) 2117–2125. [24] Y. Wei-Mon, Effects of film evaporation on laminar mixed convection heat and mass transfer in a vertical channel, Int. J. Heat Mass Transf. 35 (12) (1992) 3419–3429. [25] M.M. Rahman, H.F. Öztop, A. Ahsan, M.A. Kalam, Y. Varol, Double-diffusive natural convection in a triangular solar collector, Int. Commun. Heat Mass Transfer 39 (2) (2012) 264–269. [26] G. Barakos, E. Mitsoulis, Natural convection flow in a square cavity revisited : laminar and turbulent models With Wall functions, Int. J. Numer. Methods Fluids 18 (1994) 695–719. [27] A.J. Chamkha, Double-diffusive convection In a porous enclosure with cooperating temperature and concentration gradients and heat generation or absorption effects, Numer. Heat Transfer 41 (A) (2002) 65–87.