The deficiency of freshwater resources would cause a serious crisis in the foreseeable future. Searching for an efficient desalination method, as a result, has been surged in the past recent years. Among conventional desalination methods, membrane distillation (MD) is one of the newest technologies that is under development1,2,3,4,5. MD is a thermally driven system which water vapor is condensed by crossing the hydrophobic membrane. This technology could be used in a vast range of applications ranging from food industries to wastewater treatments. This technology benefits from several advantages. Operating in a lower temperature in comparison with other desalination systems such as multiple-effect desalination (MED) system and the ability to operate in highly saline water are two of the major advantages of this technology6,7,8,9,10. Such merits make the MD suitable for integrating with solar flat plate collectors. Several MDs are investigated such as direct contact membrane desalination (DCMD), vacuum membrane desalination (VMD), air gap membrane desalination (AGMD) and swiping gas membrane desalination (SGMD). At the present study, DCMD has been chosen with a focus on the flat sheet membrane configuration11. Countless mathematical models were developed to simulate the DCMD performance12,13,14.
Bandini et al.15 studied the effect of different operating conditions and membrane properties on the efficiency and mass transfer of the DCMD. Their results show that applying polarization coefficient alone, cannot lead to a satisfactory analysis of various resistances. They suggested novel sensitivity factors for the distillate flux based on convective heat transfer rate, membrane heat and mass transfer coefficient.
Al-Obaidani et al.13 conducted an experiment, investigating the hollow fiber DCMD exergy analysis, sensitivity study and economic evaluation in an effort to provide an optimization guideline. Hwang et al.16 carried out an experimental and analytical study on the flat sheet DCMD. They studied the impact of module dimensions and operating parameters on the DCMD performance. They also, obtained mass transfer coefficients for different module designs.
Bahmanyar et al.17 studied the impact of membrane morphology and operating parameters on the DCMD performance. The effect of those parameters on the concentration polarization and temperature was examined. Their results demonstrated that one of the most prominent parameters for permeate flux is the feed inlet temperature. Also, they showed that the feed temperature and feed circulation velocity have a meaningful correlation with concentration polarization.
Wu et al.18 investigated the impact of the membrane thickness on the transmembrane flux for DCMD. They calculated that the optimal membrane thickness increases with the reduction of the heat transfer coefficient and feed inlet temperature. However, the optimal membrane thickness decreases with the drop of the feed salinity and membrane permeability.
Bouchrit et al.19 evaluated the capability of treating hypersaline solution using flat sheet DCMD. They developed a model for mass and heat transfer and studied the effect of operating parameters on the DCMD performance. They determined the optimum operation parameters and constructed a lab-scale DCMD to concentrate the RO discharge for achieving the super-saturation saline with the aim of the salt crystallization. A long-term experiment was also constructed to monitor the permeability and the scaling of the membrane.
Imdakm and Matsuura20 developed a Monte Carlo simulation model to study the effect of membrane physical properties such as membrane pore size distribution, thickness, porosity, and thermal conductivity on the membrane flux. They proposed a DCMD’s porous membrane by a three-dimensional network model of the interconnected cylindrical pores with the distributive effective pore sizes.
Deshpande et al.21 completed a systematic numerical analysis, focusing on exergy destruction inside the hollow fiber DCMD module. They investigated the impact of different design and operation parameters on the exergy efficiency and recovery ratio of the DCMD system. They found that regardless of packing arrangement, the higher inlet feed temperature, the better recovery ratio. Also, they pointed out the important role of membrane properties on the recovery ratio.
Heat transfer is of great importance in a vast range of technologies. The very epitome of this could be seen at the heat exchangers, fuel cells and so forth. DCMD is one of the technologies which is highly dependent on the heat transfer. Multi-objective optimization, therefore, has been employed to obtain the optimized operation conditions22,23,24,25,26. Bejan et al.27,28 calculated the suitable size and shape of the heat exchanger which leads to a minimum thermodynamic loss. Bejan’s theory, constructal law, has the ability to increase the rate of the heat transfer significantly. In broad terms, constructal law can enhance the effectiveness of heat exchangers29. For instance, constructal theory can reduce the system cost by 50% in comparison with its conventional system30.
The present study tries to improve the performance of the flat sheet DCMD by combining the DCMD enhanced model with constructal law. The Nusselt correlation which considers the entry effects (both thermal and hydrodynamic) for rectangular cross-sections31 is also employed to improve the accuracy of the model. The most commonly cited issues in the literature were lack of proposing a constructal design with a thorough analysis of the DCMD design parameters. The following highlight the novelty of this research and its importance:
-
Proposing a novel constructal DCMD design.
-
Conducting a thorough analysis of DCMD design and operational parameters.
-
Conducting an evolutionary-based optimization for constructal and conventional DCMD design.
-
The performance of constructal and conventional DCMD design is compared with each other.
-
Correlations between decision parameters are obtained.
Modeling
The schematic of the flat sheet DCMD is illustrated in Fig.1. The DCMD consists of three main parts namely: feed channel, permeate channel and membrane. The whole flat plate is divided into smaller control volumes and governing equations are solved for each control volume. It is worth mentioning that because of the coupling nature of governing equations these control volumes should be co-directed along the width of the DCMD.
The schematic of the DCMD module.
While hot saltwater flows through the feed channel, cold freshwater streams in the opposite direction at the permeate channel whereas a porous hydrophobic membrane separates these water streams. This temperature difference causes a water vapor pressure gradient which is perpendicular to the membrane direction. The mass transfer through the membrane is proportional with the pressure gradient32:
$$ J_{m} = B\left( {P_{w.f.m} - P_{w.p.m} } \right) $$
(1)
where \(P_{w.f.m}\) and \(P_{w.p.m}\) are the water pressure at the feed and permeate channel, respectively. B is the membrane flux coefficient which can be evaluated by empirical or theoretical correlations. Chung et al.33 proposed an empirical equation to evaluate the mass flux through the membrane. Although the membrane flux coefficient changes with operating conditions, some of existed DCMD models select a constant value and employ it as a calibration parameter.
$$ B = \frac{\pi }{RT}\frac{\epsilon}{\tau \delta }\left[ {\left( {\frac{2}{3}\left( {\frac{8RT}{{\pi M_{w} }}} \right)^{\frac{1}{2}} r_{p}^{3} } \right)^{ - 1} + \left( {\frac{PD}{{p_{a} }}r_{p}^{2} } \right)^{ - 1} } \right]^{ - 1} $$
(2)
where R is the gas constant, \(M_{w}\) is the water molecular weight and T is the membrane temperature. The membrane properties \(r_{p}\), \(\delta\), \(\epsilon\) and \(\tau\) are pore size, thickness, porosity and tortuosity, respectively. The membrane tortuosity and porosity can be calculated as33:
$$ \tau = \frac{1}{\epsilon} $$
(3)
\(p_{a}\) is the air pressure inside the membrane pores and D is the water diffusion coefficient. The total pressure can be assumed constant in both DCMD’s channels module and equal to the atmospheric pressure. The viscous flow, as a result, would be negligible34. PD for the water–air mixture can be written:
$$ PD = 1.895 \times 10^{ - 5} T^{2.072} $$
(4)
where the water vapor pressure at each side of the membrane can be written:
$$ P_{w.p.m} \left( T \right) = P_{sat.w} \left( T \right) $$
(5)
here \(P_{sat.w} \left( T \right)\) indicates the saturation pressure of water at the temperature T. The presence of salt at the feed channel causes the deviation of vapor from pure water. To considering impurities due to the NaCl, the activity coefficient (\(\gamma_{w}\)) can be determined as follows35:
$$ \gamma_{w} = 1 - 0.5X - 10X^{2} $$
(6)
The water vapor pressure at the feed channel can be expressed as:
$$ P_{w.f.m} \left( T \right) = \gamma_{w} \left( {1 - X} \right) P_{sat.w} \left( T \right) $$
(7)
here X is the molar solute concentration. The saturation pressure of water can be calculated by the Antoine equation. \(P_{sat.w}\) is calculated as follows:
$$ P_{sat.w} \left( T \right) = {\exp}\left( {23.1964 - \frac{3816.44}{{T - 46.13}}} \right) $$
(8)
Mass flux throughout the membrane is dependent on the temperature of channels. By applying the energy equation on each control volume, the temperature distribution along the channel can be obtained. Since the membrane module consists of three segments, the energy equation applies to each separately as is demonstrated in Fig.2. In this regard, the energy equation is applied for the feed channel with the assumption of adiabatic walls:
$$ \dot{m}_{f} h_{f.b} |_{z + dz} = \dot{m}_{f} h_{f.b} |_{z} - \left( {J_{m} h_{v.f.m} + q_{m} } \right)dA $$
(9)
where \(\dot{m}_{f}\) is the water mass flow rate, \(h_{f.b}\) is the bulk enthalpy of liquid water and \(h_{v.f.m}\) represents the vapor enthalpy at the membrane surface which flows through the membrane into the permeate channel. Also, \(q_{m}\) is the heat flux through the membrane.
The schematic of control volumes.
An energy balance for the permeate channel yields:
$$ \dot{m}_{p} h_{p.b} |_{z} = \dot{m}_{p} h_{p.b} |_{z + dz} + \left( {J_{m} h_{v.p.m} + q_{m} } \right)dA $$
(10)
where \(\dot{m}_{p}\) is the water mass flow rate at the permeate channel. \(h_{p.b}\) and \(h_{v.p.m}\) represent the bulk liquid and vapor enthalpy at the permeate channel, respectively.
It is assumed that no mass is added or removed inside the membrane and all vapor flows from the feed channel into the permeate channel without any condensation in the membrane pores. Following equations, as a result, are given based on a continuity equation:
$$ \dot{m}_{f} |_{z + dz} = \dot{m}_{f} |_{z} - J_{m} dA $$
(11)
$$ \dot{m}_{p} |_{z} = \dot{m}_{p} |_{z + dz} + J_{m} dA $$
(12)
The vapor enthalpy at the boundary, \(h_{v.m}\) can be considered as the latent heat of vaporization and saturated enthalpy of water at that temperature.
$$ h_{v.f.m} = h_{f.m} - h_{fg.f.m} $$
(13)
$$ h_{v.p.m} = h_{p.m} - h_{fg.p.m} $$
(14)
Substituting previous equations into Eqs.9 and 10 yields:
$$ \dot{m}_{f} dh_{f.b} = - \left[ {J_{m} \left( {h_{fg.f.m} + h_{f.m} - h_{f.b} } \right) + q_{m} } \right]dA $$
(15)
$$ \dot{m}_{p} dh_{p.b} = - \left[ {J_{m} \left( {h_{fg.p.m} + h_{p.m} - h_{p.b} } \right) + q_{m} } \right]dA $$
(16)
here \(dh_{f.b}\) and \(dh_{p.b}\) represent the bulk enthalpy difference in the feed and permeate channel, respectively. The convective heat transfer of the membrane sides can be calculated as follows:
$$ J_{m} \left( {h_{fg.f.m} + h_{f.m} - h_{f.b} } \right) + q_{m} = h_{t.f} \left( {T_{f.b} - T_{f.m} } \right) $$
(17)
$$ J_{m} \left( {h_{fg.p.m} + h_{p.m} - h_{p.b} } \right) + q_{m} = h_{t.p} \left( {T_{p.m} - T_{p.b} } \right) $$
(18)
where \(T_{f.b}\), \(T_{p.b}\), \(T_{f.m}\), \(T_{p.m}\) are the bulk and membrane surface temperature of the feed and permeate channel, respectively. \(h_{t.f}\) and \(h_{t.p}\) are the convective heat transfer coefficient at the feed and permeate side which determined by the appropriate Nusselt number correlation. The total membrane thermal conductivity consists of two parts, the thermal conductivity of the solid part (\(k_{s}\)) and the heat conduction through the gas inside pores with the thermal conductivity of \(k_{g}\). The thermal conductivity of air and water vapor, however, have the same order of magnitude so the thermal conductivity of water vapor is considered as \(k_{g}\). By applying the parallel heat flow assumption through the gas inside pores and solid parts, the membrane thermal conductivity \((k_{m} ) \) can be calculated as follows 36:
$$ k_{m} =\epsilon k_{g} + \left( {1 -\epsilon } \right)k_{s} $$
(19)
The thermal conductivity of water vapor and PTFE membrane, which is employed in this study, are depended on the temperature. These values can be calculated using the following equations:
$$ k_{s} = 0.087 + 6.0 \times 10^{ - 4} T $$
(20)
$$ k_{v} = 2.72 \times 10^{ - 3} + 5.71 \times 10^{ - 5} T $$
(21)
Previous equations are valid between the temperature of 273–373K. It is worth mentioning that the convection heat transfer within the membrane pores is negligible:
$$ q_{m} = k_{m} \frac{{T_{f.m} - T_{p.m} }}{{\delta_{m} }} $$
(22)
where \(\epsilon \) and \(\delta_{m}\) are porosity and membrane thickness, respectively.
For solving the above equations along the channel length, a proper Nusselt number is needed for each section of the channel. Various correlations based on the geometry and flow regime are suggested and used in the previous studies. At the present study, a general relationship for an arbitrary cross-sectional shape, considering the entry effects (both thermally and hydrodynamically), is applied31.
$$ Nu_{\sqrt A } \left( {z^{*} } \right) = \left[ {\left( {\frac{{f\left( {Pr} \right)}}{{\sqrt {z^{*} } }}} \right)^{m} + \left( {\left\{ {0.501\left( {\frac{{fRe_{\sqrt A } }}{{z^{*} }}} \right)^{\frac{1}{3}} } \right\}^{5} + \left\{ {3.86\left( {\frac{{fRe_{\sqrt A } }}{{8\sqrt \pi \varepsilon^{\gamma } }}} \right)} \right\}^{5} } \right)^{m/5} } \right]^{1/m} $$
(23)
Here \(m\) is the blending parameter which is given by:
$$ {\text{m}} = {\text{a}} + {\text{b}} Pr^{1/3} $$
(24)
where a = 2.27 and b = 1.65. The product of friction factor (f) and Reynolds number which is defined based upon the square root of cross-sectional area, is given by:
$$ fRe_{\sqrt A } = \frac{12}{{\sqrt \varepsilon \left( {1 + \varepsilon } \right)\left[ {1 - \frac{192\varepsilon }{{\pi^{5} }}\tanh \left( {\frac{\pi }{2\varepsilon }} \right)} \right]}} $$
(25)
\(f\left( {Pr} \right)\) for the uniform wall heat flux boundary condition is calculated as follows:
$$ f\left( {pr} \right) = \frac{0.886}{{\left[ {1 + \left( {1.909Pr^{1/6} } \right)^{9/2} } \right]^{2/9} }} $$
(26)
where \(\varepsilon\), \(z^{*}\), Pr and \(\gamma\) are the aspect ratio of duct, dimensionless position of thermally developed flows, Prandtl number and shape parameter, respectively. \(\gamma\) is considered − 0.3 due to the low aspect ratio of the duct. The input parameters of the DCMD is represented in Table1.
Validation
Based on the equations proposed in the modelling section, the DCMD model is developed and validated with the experimental data which are reported by Martı́nez-Dı́ez and Vazquez-Gonzalez37. The PTFE membrane (with 80% void fraction, 60μ\(\mathrm{m}\) thickness, 0.2mm nominal pore size) was implemented to perform experimental tests. The membrane module in that study was a counter flow flat sheet with 9 feed channels and 9 permeate channels. The water mass flux (\({J}_{m})\) versus the average bulk temperature of the feed saltwater (\({T}_{b1}\)) are compared in Fig.3 for two recirculation rates (7 and 11cm3/s) at various salt feed water concentrations (0, 0.55, 1.15 and 1.67 molar NaCl). As is represented in Fig.3 the simulation results are in an excellent agreement with experimental data. It is worth mentioning that \({T}_{b1}\) in Fig.3 is the average of the inlet and outlet bulk temperature of the feed channel.
Water mass flux across membrane (a) Inlet flow rate, Q = 7cm3/s (b) Inlet flow rate, Q = 11cm3/s. Symbols indicate experimental data25 and lines show model results.
Constructal design and optimization
Speed-constrained Multi-objective particle swarm optimization (SMPSO) is an improved particle swarm optimization (PSO) characterized by the use of a strategy to limit the velocity of the particles. PSO is a bio-inspired evolutionary-based which imitating social behavior of bird flocking or fish schooling38. SMPSO use a strategy to put a constraint on particle velocity. In this regard, SMPSO produces new particles when the velocity becomes too high. Also, SMPSO benefits from polynomial mutation as a turbulence factor and external archive to store the non-dominated solutions found during the search39.
Constructal law, which was proposed for the first time in 199640, states that for a flow system to persist in time it needs to evolve in a way that it provides easier access to its current41. The very epitome of this notion could be seen everywhere in nature42. Dendritic shapes, which is ubiquities in nature, are an icon of constructal design. These patterns could be seen in blood vessels, river basins and most recently the quantum footprint43,44,45. The present study tries to implement constructal law to enhance DCMD performance. While, the water mass flow rate and the size of channels are kept constant, as is demonstrated in Fig.4, a new configuration is adopted to give the flow architecture the freedom to morph. The hot brine after entering in the feed channel is distributed into two other equally spaced feed channels (1 \(\to\) 2,3). The same design is also considered for permeate channels. This configuration is also known as a canopy-to-canopy design.
Schematic of constructal DCMD structure.
Two different designs including constructal and conventional DCMD are considered for optimization. The schematic diagram of the constructal DCMD is illustrated in Fig.4. The constructal DCMD consists of three parts. Each part of the constructal DCMD is one-third of the conventional DCMD. The inlet volume flow rate is divided into two equal portions while entering part 2 and 3. Nine design parameters including the inlet volume flow rate at the feed side, inlet volume flow rate at the permeate side, inlet feed water salinity, bulk temperature of the feed side, bulk temperature of the Permeate side, membrane length, membrane width, porosity and membrane pore size are considered for both conventional and constructal design. SMPSO algorithm is employed for optimization. The thermal efficiency and water mass flux are considered as two objective functions. The thermal efficiency of the DCMD can be defined as the ratio of latent heat of vaporization to the total heat transfer across the membrane which can written as:
$$ \eta = \frac{{J_{m} h_{fg} }}{{J_{m} h_{fg} + \frac{{k_{m} }}{\delta }\left( {T_{f.m} - T_{p.m} } \right)}} $$
(27)
The amount of the water mass flux (\(J_{m}\)) is calculated based on the previous equations. At the present study, the main goal of the optimization is to maximize both objective functions. The applicable ranges for relevant decision variables in this study are listed in Table2.
