Interactive comment on “ Modeling evaporation processes in a saline soil from saturation to oven dry conditions ”

First of all, I wish to express my sincere appreciation for the author’s courage and competence to model such a complex experimental system. However, in the current state, the numerical simulations seem to be incompletely described which makes it quite hard to understand them in detail. Furthermore, the simultaneous occurrence of a multitude of processes makes the studied system so complex that the question arises whether it is really possible to discriminate the various processes and their influence on the simulated system response. Directly related to this is the need to improve the description of parameter estimation by inverse modeling. As noted by the authors (P540 L19.21) ’all parameters influence all processes’. I would therefore acknowledge if the issues


Introduction
Understanding evaporation is necessary in many fields of earth system sciences (Shuttleworth, 2007).In fact, soil evaporation is crucial in controlling the balance of soil-Correspondence to: M. Gran (meritxell.gran@idaea.csic.es)surface water and energy in arid and semiarid areas (Saito et al., 2006).The actual mechanisms controlling evaporation are intricate (Sakai et al., 2009).Soil evaporation may be controlled by the soil-atmosphere boundary layer when the soil is moist or by hydraulic conditions when it is dry (Schneider-Zapp et al., 2010).In the latter case, evaporation causes the soil to dry and heat up causing liquid, vapor and heat fluxes to interact.The presence of solutes increases the complexity of the system and exacerbates the consequences, leading to salinization.
A number of researchers have analyzed this problem from an experimental perspective (Wheeting, 1925;Scotter, 1974;Nassar and Horton, 1989;Scanlon, 1992).They conclude that water flux in dry and salinized soils is controlled by salinity and temperature gradients.Salinity causes water activity to drop, thus reducing vapor pressure in equilibrium with liquid water and driving vapor towards the saltier zone.Evaporation depends also on temperature and absorbs energy.Thereby, evaporation is affected by water flow and energy and solutes transport.The interaction of matric potential, temperature and salinity gradients under very dry conditions was studied by Gran et al. (2011), who observed salinity to decrease below the evaporation front, which they attributed to condensation of downward vapor flux.Unfortunately, experimental studies do not yield direct measurements of flow and phase change processes, which must be indirectly inferred from state variable measurements.This is not easy when the phenomena are complex and coupled.Therefore, quantitative understanding of the above processes requires mathematical modeling.dry, and that water flows in the liquid phase below the evaporation front.A good approximation to water table evaporation under isothermal conditions was obtained by Gowing et al. (2006), who divided the soil into liquid flow and vapor flow zones separated by an evaporation front.The width of this front is subject to debate: Gran et al. (2011) observed a sharp front, whereas Konucku et al. (2004) concluded that a sharp phase transformation could not be expected.Notwithstanding, these models do not consider the role of salinity.
The effect of high salinities was modeled by Nassar and Horton (1989), who simulated water transport in unsaturated nonisothermal salty soil on the basis of steady-state heat and mass transfer.Ironically, salinity effects in theirs and other models have commonly been analyzed assuming dilute solutions, which is not suitable since vapor pressure in equilibrium with a salty aqueous phase is very sensitive to salinity (Burns et al., 2006).This explains the difficulties encountered by Nassar et al. (1992) when modeling evaporation from salty solutions.
In addition to the effect of salinity or water activity, two other factors must be borne in mind when modeling evaporation from high salinity solutions.First, salts tend to precipitate in the pores (Nachshon et al., 2011) and/or to form a low permeability crust that should be modeled (Yakirevich et al., 1997).Crust formation was modeled but not compared with experimental data by Olivella et al. (1996a) and both modeled and compared with data by Fujimaki et al. (2006) though assuming isothermal conditions.Second, under hot and/or dry conditions, the residual saturation of water in soil can no longer be considered a lower bound for saturation (Milly and Eagleson, 1982;Rossi and Nimmo, 1994;Prunty, 2003).A modification of the retention curve must therefore be considered to represent water contents under oven dry conditions.
The emerging picture is complex.A nonisothermal multiphase flow model is necessary to distinguish advective and diffusive vapor fluxes.High concentrations near the surface require using non ideal solution chemistry (e.g.Pitzer) to simulate osmotic effects on vapor pressure and salt precipitation.Mass balances of water, air, heat and solutes are necessary and the effects of thermal, suction and salinity gradients must be simulated interacting simultaneously.Furthermore, traditional continuum mechanics descriptions (using Darcy and Ficks Laws, or concepts such as retention curve) may not suffice to include energy concepts.For instance, the role of film flows and the use of an enhancement factor for vapor diffusion are controversial.On one hand, molecular dynamics simulations (Odelius et al., 1997) and experimental observations (Hu et al., 1995;Miranda et al., 1997) suggest that adsorbed water may display a nearly crystalline structure (water would be virtually immobile).On the other hand, Rakhmatkariev (2006) argues that the mobility of water adsorbed onto muscovite is only slightly less than in the bulk liquid, and Shokri et al. ( 2010) describe a continuous water film feeding the evaporation front.There is also discussion in the use of an enhancement factor for vapor diffusion (Heit-man et al., 2008;Sakai et al., 2009) or the need for the Dusty Gas Model (Thorstenson and Pollock, 1989).Representing these concepts in the presence of high salinity and temperature gradients is currently under investigation.Therefore, the question is whether the traditional continuum mechanic model can represent the processes involved in soil salinization.To this end, the experiments of Gran et al. (2011) are ideal since they are relatively simple (independently characterized homogeneous sand and controlled boundary conditions) and all relevant state variables (water content, temperature and salinity) were measured.Yet, they appear to be sensitive to all of the above processes.
Therefore, the aim of this work is to test whether traditional models can be used to reproduce the experiments of Gran et al. (2011) in order to (a) evaluate the magnitude and direction of the water fluxes and gain a greater understanding of the downward vapor flow mechanism, (b) describe the evolution and location of condensation-evaporation, and (c) to assess the relevance of the matric potential, temperature and osmotic gradients in controlling the aforementioned water separation process.

Evaporation experiment and conceptual model
Laboratory experiments consisted of open sand columns initially saturated with a 14 g kg −1 epsomite (MgSO 4 × 7H 2 O) solution.Evaporation was forced by an infrared lamp so that radiation at the soil surface was similar to the summer radiation at mid-latitudes.The experiment continued until the overall saturation fell to 0.32, which was identified after preliminary tests to ensure the occurrence of a well developed crust and a deep evaporation front to allow the study of vapor fluxes.At this stage, the columns were dismounted to measure vertical profiles of temperature, volumetric water content and solute concentration.Two groups of four identical columns were dismounted at different times to obtain the time evolution of those profiles.These columns were dismantled sequentially after reaching saturations of 74 % (after 2 days of evaporation), 50 % (after 4 days), 40 % (5 days) and 32 % (12 days).
Both column experiments and results are described in Gran et al. (2011) and a diagram is shown in Fig. 1.The conceptual model emerging from these experiments requires tighly coupled processes.Capillarity brings about an upward liquid flux and soil drying.The liquid flux transports solutes by advection towards the top, where evaporation leads to a dramatic increase in concentration.A continuous salt crust is formed at the surface.This leads to changes in soil properties increasing thermal conductivity and reducing porosity and vapor diffusivity.Evaporation also reduces the water content and the unsaturated hydraulic conductivity causing the evaporation front to move downwards.this front is very narrow.A water separation process occurs at the front.On the one hand, concentrations are high above the front, where water flow is restricted to vapor phase.On the other hand, underneath the evaporation front, concentrations are diluted below the initial values.Dilution is probably a result of the vapor pressure increase at the front caused by evaporation.Below the front, the temperature and vapor pressure gradient lead to a downward vapor flux.That is, vapor flows not only upwards from the evaporation front but also downwards.Condensation of this downward vapor flux causes the dilution.

Processes and governing equations
The system is governed by thermohydraulic and geochemical processes.

Thermohydraulic processes
The thermohydraulic model focuses on the mass balance of water (liquid water and vapor) and air (dissolved in water and in the gas phase) in terms of pressure, and the energy balance in terms of temperature.The equations of water and air mass balance are: where subscripts l and g refer to liquid and gas and superscript w and a refer to water and air.ω is the mass fraction (kg kg −1 ) of a component in a phase, ρ is the density (kg m −3 ) of a phase, S is the hydraulic saturation (m 3 m −3 ), φ is the porosity (m 3 m −3 ), j (kg m −2 s −1 ) is the total flux (advective, diffusive and dispersive) and f is an external source/sink term (kg m −3 s −1 ) (i.e.top boundary condition and mineral hydratation-dehydratation).Note that the firsts two terms in the equations represent the change of mass of water (Eq. 1) or air (Eq.2) in the liquid and gas phase respectively and the third and fourth terms represent the fluxes of water (Eq. 1) or air (Eq.2) in liquid and gas phase respectively.
The energy mass balance is written as: where E is the specific internal energy, i c is the energy flux (J m −2 s −1 ) owing to conduction through the porous medium, the other fluxes (j El , j Eg ) are advective fluxes of energy (J m −2 s −1 ) caused by mass motions and f Q is an internal/external supply (J m −3 s −1 ) that accounts for boundary conditions at the top (i.e.heat entry from the lamp) and at the laterals (i.e.heat exit through column walls).
A state variable is associated with each mass balance: liquid pressure (P l ), gas pressure (P g ) and temperature (T ).Constitutive laws must be used to express the mass balance equations as a function of the state variables.The fluxes and constitutive laws that control these balances are shown in Table 1.Further details on mass balance equations and constitutive laws are described by Olivella et al. (1994Olivella et al. ( , 1996b)).Still, it is worth discussing some of the simplifications implicit in Table 1.First, gaseous components diffusion is simulated using Fick's Law.Although it is considered less accurate than the Dusty Gas Model (Thorstenson and Pollock, 1989), the difference is small when permeability or pressure is high.According to Webb and Pruess (2003) results, differences would be negligible in our case.Second, the effects

Oven dry conditions
As mentioned above, under oven dry conditions (that is, saturations well below residual saturation), the residual saturation can no longer be considered a lower bound for saturation.To reproduce the experimental data (Gran et al., 2011) a modification of the retention curve and relative permeability functions is necessary.To address this, Milly and Eagleson (1982) simply considered the residual saturation to be zero; Rossi and Nimmo (1994) proposed a different function to extend the capillary curve towards fully dry conditions; and Prunty (2003) used the zero value in standard retention curve models but modified the relative permeability function for the dry range.
The van Genuchten (1980) model is widely used under moist conditions but requires modification to represent oven dry.The assumption is that soil can reach full drying, i.e. if evaporation takes place in an oven at 105 • C or near the surface under a dry or hot atmosphere (Ross et al., 1991).The van Genuchten retention curve is: where P c is capillary pressure, P 0 is related to the capillary pressure required to desaturate the soil and λ is a shape parameter of the function.This equation permits to calculate the effective saturation (S e ) as a function of a minimum saturation S i and the actual saturation (S l ): In order to extend this curve for high suctions (i.e.conditions of drying by evaporation) and to represent the oven dry branch that goes from the residual saturation to the ovendryness, the minimum degree of saturation is expressed as follows: where there are three parameters: α, S 0 and P dry c .The latter can be identified with the capillary pressure for the Hydrol.Earth Syst.Sci., 15,[2077][2078][2079][2080][2081][2082][2083][2084][2085][2086][2087][2088][2089]2011 www.hydrol-earth-syst-sci.net/15/2077/2011/ dry material and can be considered equal to P ovendryness c = 1000 MPa.However, lower values may be considered if dryness is induced by atmospheric conditions that are less extreme than oven dryness.S 0 is the residual saturation, for which liquid water becomes discontinuous so that liquid permeability is zero.Parameter α scales the transition from the van Genuchten branch to the oven dry branch of the proposed retention curve.Its value may be chosen as α = 1/ln(P dry c /P 0 c ), where P 0 c is the capillary pressure for which the oven dry branch crosses S 0 (i.e., S i = S 0 ).It should not be chosen smaller than P 0 (or else, P c will remain low for S l much smaller than S 0 ) nor very large (or else, the oven dry branch will separate from the van Genuchten branch for S l much larger than S 0 ).We adopted here P 0 c = 0.03 MPa, for which α = 0.1, after some trial and error attempts.The remainig parameters (P 0 , λ and S 0 were obtained by fitting independent measurements of suction and saturation).This proposed retention curve is a continuous function with continuous derivatives.A similar form was already proposed by Fayer and Simmons (1995).
For the relative permeability function, we assume that the distribution of the liquid phase is unaffected by the above modification and only depends on saturation.Therefore, it is given by , S l > S 0 (7) where S ep is the effective saturation for permeability.Note that for saturations below S 0 the capillary pressure can be calculated from the retention curve, but the relative permeability is zero.This allows representing water isolated in the meniscus that can not flow as a liquid phase but can still evaporate.Figure 2 compares the proposed model and the original van Genuchten model, in terms of retention curve and shows the adopted relative permeability curve.The parameters used here (see Table 1) are obtained by manual calibration adjusting a numerical model using the data from an experimental retention curve.

Reactive transport
The mass balance used for reactive transport can be written as where vector c a (mol kg −1 ) is the concentration of aqueous species and L l is the linear operator for the advection, dispersion/diffusion, m l is the non-chemical sourcesink (mol m −3 s −1 ) and D l is the dispersion/diffusion tensor (m 2 s −1 ).The values used for the diffusion coefficient and dispersivity are 10 −9 m 2 s −1 (Lasaga, 1998) and 0.0015 m (between values used by Aggelopoulos andTsakiroglou, 2007 andZheng et al., 2008) respectively.R contains the rates of the kinetic reactions (r min ) for all the different mineral phases (r eps ,r hex ,r pent and r stark ), σ min is the mineral reactive surface and min is the ratio between the ion activity product and the equilibrium constant.The reaction rates enable us to calculate the liquid water provided by mineral hydratation-dehydratation.This amount of water is added to the water mass balance as a source/sink term (f w ) in Eq. ( 1).
where m w is the molecular weight of water.

Numerical model
The problem is considered one dimensional in a vertical direction.The grid is made up of 240 elements (for 24 cm column length).The medium is homogeneous except for the top 1.5 cm, where a reduced value of τ 0 (gas diffusion enhancement factor) was used to reproduce the increase in tortuosity caused by precipitated salts in the crust (recall that τ 0 includes the effects of tortuosity, constrictivity and diffusion enhancement).The gas diffusion enhancement factor (τ 0 ) (see Table 1) was calibrated to be 1.2 in the upper zone and 8 below.These values were required to fit the observed evaporation and water content profiles and, at the same time, the observed dilution.Boundary conditions (BC) for liquid, vapor and heat were chosen to reproduce the laboratory conditions (see Table 2).The top boundary is a mixed condition representing gas (air and vapor) and heat inflow-outflows.A radiative heat flux (from the lamp) was added at the top boundary.The lateral and bottom BC were of no-flow for water and solutes, but energy was permitted to dissipate across the insulating layer.
The column was initially saturated in water (initial gas and liquid pressures were P g = P l = 0.101325 MPa).The initial temperature and porosity were T 0 = 25 • C and φ =0.4,respectively.
Numerical simulations were carried out using the RETRASO-CODE BRIGHT (RCB) code, which couples the thermohydraulic model CODE BRIGHT (CB) of Olivella et al. (1996b) with the reactive transport model RETRASO of Saaltink et al. (2004), which incorporates the approach of Saaltink et al. (1998).Furthermore, geochemical calculations are performed with the object-oriented chemical module CHEPROO (Bea et al., 2009(Bea et al., , 2010)), which includes high salinity solutions using the equations of Pitzer (1973).The feedback of reactive transport in thermohydraulics is performed by a time lag approach.The code solves the thermohydraulic equations (Eqs. 1 to 3) for one time step.Results (Darcy fluxes, hydraulic saturation, etc) are used to calculate the reactive transport for the same time step (Eq.8).

Cap. Pre
Van Genuchten Proposed 1.E-6 Rel. P Energy flux j 0 e = 750 J s −1 j 0 e = 0 j 0 e = 0 j e = j 0 e + γ e (T 0 − T ) Where superscript 0 refers to prescribed boundary values, and subscripts l, g and s to liquid, gas and salinity, respectively.P g is gas pressure, ρ is density, ω w g is the mass fraction of vapor in gas (corresponding to the measured relative humidity above the column), γ g is leakage coefficient for gas advective flux, β g is the leakage coefficient for vapor non-advective flux, j 0 e is a prescribed heat flow from the lamp, γ e is the energy transfer coefficient for energy flux and E w g is the internal energy of water in gas phase per unit mass of water.
Subsequently, thermohydraulic properties such as porosity changes due to precipitation-dissolution (precipitated mass per unit volume is divided by the solid density of the salts and the volume variation obtained corresponds to the change in porosity) or water activity (computed from the osmotic coefficient (Felmy and Weare, 1986)) and the source/sink term (f w l ) are calculated using the reactive transport results.A new thermohydraulic time step is calculated using these new properties.

Results and discussion
Figure 3 displays the water saturation, temperature and salinity profiles computed for four different times along with the experimental results at the end of the experiment, after 12 days.Saturation profiles illustrate the progressive desaturation of the columns from the top.The water content drops over time to values near residual saturation at a depth that increases with time.Saturation at the top reaches oven dry conditions (water saturations close to zero).The bottom of this zone represents the location of the evaporation front.Below the front, water content continues to increase downwards, leading to a degree of saturation profile similar to that of the sand retention curve.The good match between model and experiment at the upper oven dry area (above 4 cm depth) supports the validity of the retention curve modification (note that without this modification reproducing the step in the saturation profile would not have been possible), which 3.3d 24 Fig. 3. Profiles of saturation, temperature and salinity measured at the end of the experiment (symbols) and computed (lines).The time evolution is shown for four different times (after 1.1 days, 3.3 days, 6.6 days and, at the end of the experiment, 12 days).Note that salinity is expressed as concentration in water to facilitate the analysis of mass transport processes.
improves the simulation of multiphase flow under very dry conditions.
Temperature rises during the experiment and displays a slope change at the evaporation front.The temperature gradient is larger above than below the front because the evaporation front acts as a heat sink.Another smaller temperature slope change can be detected at a depth of 1.5 cm for the temperature profiles at days 6.6 and 12.This is due to the increase in thermal conductivity caused by salt precipitate at the crust (note that conductivity depends on porosity, see Table 1).
The spatial distribution of concentration (expressed as salt mass per unit mass of water to facilitate the analysis of mass transfer processes) is noteworthy.Salinity is extremely high at the surface, where the water content is negligible, reaching salt solubility and producing precipitates.This high concentration zone grows with time, advancing in depth with the evaporation front.Immediately below, salinity drops sharply to values underneath the initial concentration.The minimum concentration is always located immediately below the evaporation front.Further down, salinity rises slightly with depth, but is still more diluted than the initial conditions.A difference between the experimental data and the numerical model is observed: the minimum in the simulated concentration is smaller than the measured one.The water content and temperature profiles do coincide with what might be expected (drier and warmer conditions at the surface than at depth) unlike the concentration profile.Most traditional models (e.g.Huinink et al., 2002) predict a maximum concentration at the evaporation front and a smooth monotonic reduction downwards toward initial concentration, caused by downwards diffusion.The radically different behaviour observed in our concentration profiles can be attributed to vapor fluxes.Actually, the time evolution of cumulative evaporation (Fig. 2 of Gran et al., 2011) evolves according to the traditional model (e.g.Sghaier et al., 2007), but the evolution of salt concentration profile does not.
The fact that the model reproduces qualitatively the observations suggests that it can be used to determine the role of water flow, heat and transport processes in the system.We discuss below the mechanism responsible for the dilution of the solution.
Figure 4 displays the profiles of water and heat fluxes for the same instants as in Fig. 3. Liquid water flows upwards because of capillarity throughout the experiment.An evaporation front, located where the liquid flux drops abruptly to zero, may be observed after 3.3 days.Water vapor flux profiles show that vapor flows both upwards and downwards from the evaporation front.This front advances deeper into the soil as the feeding liquid flux from the bottom diminishes over time.Above this front, water can no longer flow as a liquid.Nevertheless the water content continues to diminish towards the top of the column in response to some residual evaporation above the front.Note that the water content must be very low (very high suction) and/or the salinity very high (high osmotic effect) to ensure an upwards vapor flux despite the upwards increase in temperature.Both factors, high suction and salinity, contribute in this case to the upwards vapor flux above the evaporation front.None of these two factors contribute significantly below the evaporation front, so that vapor diffuses downwards according to the temperature gradient.Condensation of this downward vapor flux accounts for the decrease of the vapor flux at the bottom of the column and for the dilution of the solution below the evaporation front.Both upward and downward vapor fluxes are present throughout the experiment.The numerical model also enables us to study in detail heat fluxes.Figure 4 displays a conductive heat flux downwards throughout the column.This flux is larger in the upper zone, where the soil is dry, and decreases over time.Note the sudden fall in heat flux observed at all time steps at the evaporation front, which is due to the heat sink produced by evaporation.The advance of the evaporation front is observed very clearly and its location is controlled by the above heat flux, which decreases with depth and dryness.The total heat flux graph exposes that conductive and advective (latent) heat fluxes are similar in magnitude but in opposite directions above the front.Below, both fluxes are directed downwards although the advective flux is the dominant one.
Figure 5 displays the spatial distribution of the evaporation and condensation rates and the vapor mass fraction profile.The vapor mass fraction profile presents a marked change in the slope and provides evidence of two distinct gradients.Since evaporation occurs at this juncture, the increase in vapor and gas pressure generates vapor diffusion and advection both upwards and downwards.The graph on the right displays the evolution of the evaporation (negative values) and condensation (positive values).The evaporation rate is higher at the start and decreases as the evaporation front advances deeper into the soil.Note that condensation always occurs below the evaporation front and that its magnitude is substantially lower than that of the evaporation, causing the soil to dry.Condensation extends throughout the relatively dry portion of the column, from residual saturation to saturations around 0.6 (compare the inset of Fig. 5  the saturation profiles in Fig. 3).Vapor condensation occurs from the early stages and its maximum evolves decreasing as it advances into the soil because vapor diffusion is hindered by liquid water (recall Table 1, that gas diffusivity is multiplied by S 3 g ).This explains the decrease in concentration below its initial value.The fact that the model underestimates the reduction in concentration can be attributed to an excessive reduction of diffusivity or to further gas diffusion enhancement.

Sensitivity analysis
Further insight into the above processes and into the role of controlling parameters can be gained from a sensitivity analysis.Processes are strongly coupled, i.e. all the parameters affect all the processes.We focused on the sensitivity of the boundary heat dissipation (γ ) and the gas diffusion enhancement factor (τ 0 ), which proved to be more illustrative.Table 3 presents the values adopted for these parameters.
Figure 6 illustrates the effect of increasing heat dissipation across the walls.The effect of the heat dissipation through the bottom has also been studied but is not shown here (it presents a similar system response but less relevant).Increasing γ at the walls causes an increase in the rate of sensible heat dissipation, thus having less energy available for evaporation.Although Fig. 6 shows higher evaporation at the evaporation front, this is overcompensated by a higher condensation below the front.As a result the overall saturation increases and the evaporation front remains closer to the surface.It also leads to a lower overall temperature.Note that as more heat is dissipated the temperature gradient just beneath the evaporation front also increases.This causes an increase in the downwards vapor pressure gradient and in the downward vapor flux, whereas, the upward vapor flux diminishes.The shape of the condensation profile varies to give a bigger maximum concentrated just below the evaporation front.We can infer that the amount of heat dissipated through the walls controls the thermal gradient in the column.
The gas diffusion enhancement factor (τ 0 ) was homogenized increasing its value in the upper part of the column (Table 3).As a result, the soil dries faster and the overall degree of saturation diminishes.The evaporation front advances deeper into the soil and the area below is colder.The temperature profile changes from three to two different gradients: the temperature gradient near the surface increases and becomes uniform towards the evaporation front.As the lower temperature leads to less lateral heat dissipation, the temperature gradient below the evaporation front decreases.Accordingly, the downward vapor flux diminishes and hence, the condensation.By contrast, the upward vapor flux increases, showing that vapor can flow more easily towards the surface and explaining the downward displacement of the evaporation front and the column drying.As the evaporation front advances deeper into the soil, condensation takes place further down.Overall, the results applying a different value of τ 0 in the upper part of the column, fit to experimental data suggesting that variations in τ 0 must be included to improve our model.Table 3. Studied parameters for the sensitivity analysis: boundary heat dissipation (at the walls and at the bottom) and gas diffusion enhancement factor (τ 0 ).Compared to the base model (BM): the boundary heat dissipation, by means of γ value, has been doubled and increased by an order of magnitude alternatively and τ 0 value for the upper material (firsts 1.5cm) has been increased from 1.2 to 8 to equal the value for all the column.

Conclusions
The model reproduces quite accurately experimental observations of varied nature (temperature, water content and salt concentration).Most model parameters were either measured (retention curve) or derived from the literature (constitutive laws in Table 1).The only parameter that had to be varied significantly to obtain a good match was the vapor diffusion enhancement Reliability of model parameters and the good between observations outputs support the of the model, which prompted to analyze and quantify the computed processes.These confirm the initial about the highly coupled and nature of evaporation from a soil.In essence, evaporation is driven heat, but it can be limited by liquid and vapor flux processes and salinity.Some conclusions can inferred from model: Hydrol.Syst.Sci., 15, 2077-2089, 2011 M. Gran et al.: Modeling evaporation processes in a saline soil 2087 • Evaporation causes vapor pressure to increase at the evaporation front causing vapor to flow upwards and downwards.Both fluxes occur throughout the experiment, but the relative importance of the downwards flux increases over time.In our model the downward flux is half that of the upward flux at the end of the experiment.
• The evaporation front is very narrow, which contradicts the analysis of Konucku et al. (2004).Most evaporation is concentrated in less than 1cm.Some evaporation occurs above the front, but condensation starts immediately below.This finding may be due to the experimental conditions (loss of heat through the column walls).Without this heat loss, water vapor could have penetrated further into the soil, which is consistent with the findings of Scanlon and Milly (1994).
• Condensation of the downward vapor flux dilutes the solution beneath the evaporation front with the result that salinity drops below the initial value.This finding confirms the existence of a water separation process driven by evaporation, which contradicts the current models of soil salinization and formation of eflorescences.However, reproducing the concentration profiles required a large gas diffusion enhancement factor.
• Heat flows downwards from the surface to the bottom of the column mainly by conduction.Still, advection of latent heat is also relevant.Upwards advection almost compensates downwards conduction above the evaporation front, so that the total heat flux is not very sensitive to the location of the front.Immediately below the front, advection is comparable to conduction, but it diminishes further down as vapor condensate.In short, vapor diffusion leads to a highly relevant energy transport mechanism within the soil.
In summary, our model supports the traditional division of the soil into a virtually dry (and/or saline, as we have seen) region and a moist region, separated by an evaporation front.Vapor diffusion is the only relevant water flow mechanism in the upper region.This view may suffice for evaluating evaporation rates and water mass balances.However, assessing salt processes requires acknowledging that not only liquid water but also vapor diffusion occur below the evaporation front.
Further conclusions can be drawn from the sensitivity analysis: • Vapor diffusion is very sensitive to the heat boundary conditions.The amount of heat dissipated throughout the column walls controls the temperature gradient.Downward diffusion is enhanced by the lateral heat dissipation through the walls.
• It is the temperature gradient more than the temperature range, what governs the magnitude of vapor fluxes.Therefore, these processes can occur in salinized soils under temperatures lower than the ones observed here.
• The formation of a crust due to salt precipitation reduces porosity and increases tortuosity, which hinders evaporation.These have been simulated by increasing tortuosity at the crust (actually, reducing τ 0 from 8 in the column to 1.2 in the crust).This was necessary to reproduce evaporation rates together with observed temperatures and salinity profiles.
Finally, further research is warranted to resolve a number of issues.The dilution simulated by the numerical model is always lower than that measured in the experiments.The representation of vapor flux may not be accurate, diffusion enhancement factors are very large, which demands deeper analysis of this issue.

Fig. 1 .
Fig. 1.Diagram of the design of the evaporation column experiments and their conceptual model from day one (saturation) to the end of the experiment (oven dry conditions).The water fluxes are on the left, the salt fluxes are on the center and the energy ones on the right.The dashed and point lines show the location of the evaporation front and the minimum in salt concentration respectively.

MFig. 4 .
Fig. 4. Computed profiles of liquid, vapor and total water fluxes (above), and conductive and total (conductive plus advective) heat fluxes (below).Positive and negative values stand for upward and downward flows, respectively.The simulation results are shown for four different times (after 1.1 days, 3.3 days, 6.6 days and 12 days).The difference between the diffusive vapor flux and the total water flux, displayed in the top 4cm of the above graph on the right, is equal to the advective vapor flux.

Fig. 5 .
Figure5displays the spatial distribution of the evaporation and condensation rates and the vapor mass fraction profile.The vapor mass fraction profile presents a marked change in the slope and provides evidence of two distinct gradients.Since evaporation occurs at this juncture, the increase in vapor and gas pressure generates vapor diffusion and advection both upwards and downwards.The graph on the right displays the evolution of the evaporation (negative values) and condensation (positive values).The evaporation rate is higher at the start and decreases as the evaporation front advances deeper into the soil.Note that condensation always occurs below the evaporation front and that its magnitude is substantially lower than that of the evaporation, causing the soil to dry.Condensation extends throughout the relatively dry portion of the column, from residual saturation to saturations around 0.6 (compare the inset of Fig.5 and

Fig. 6 .
Fig. 6.Analysis of the effect of boundary heat dissipation.Computed profiles of saturation, temperature, concentration, water mass flux, evaporation and condensation after 12 days.

Table 1 .
Constitutive laws, parameters and values used in the numerical model.

Table 2 .
Retention and relative permeability curves (original van Genuchten and proposed) for the sand used in the column experiments.Liquid, vapor and heat boundary conditions and corresponding parameters.