Soil-atmosphere interaction in earth structures

Yujun Cui

Laboratoire Navier,Ecole des Ponts ParisTech,France

Keywords:Environmental chamber Soil water evaporation Suction-based evaporation model Embankment Field monitoring Numerical modeling

ABSTRACT The soil-atmosphere interaction was investigated through laboratory testing,field monitoring and numerical monitoring.In the laboratory,the soil water evaporation mechanisms were studied using an environmental chamber equipped with a large number of sensors for controlling both the air parameters and soil parameters.Both sand and clay were considered.In case of sand,a dry layer could be formed during evaporation in the near surface zone where the suction corresponded to the residual volumetric water content.The evaporative surface was situated at a depth where the soil temperature was the lowest.In case of clay,soil cracking occurred,changing the evaporative surface from one-dimensional to three-dimensional nature.The suctionbased evaporation model was adapted to take these phenomena into account by adopting a function of dry layer evolution in the case of sand and by adopting a surface crack ratio and a retative humidity ratio in the case of clay.In the field,the volumetric water content,and the suction as well as the runoff were monitored for an embankment constructed with lime/cement treated soils.It appeared that using precipitation data only did not allow a correct description of the variations of volumetric water content and suction inside the soils,the consideration of water evaporation being essential.It was possible to use a correlation between precipiration and runoff.The hydraulic conductivity was found to be a key parameter controlling the variations of volumetric water content and suction.For the numerical modelling,a fully coupled thermohydraulic model was developed,allowing analyzing the changes in temperature,volumetric water content and suction of soil,with the upper boundary conditions at the interface between soil and atmosphere determined using meteorological data.Comparison between simulations and measurements showed the performance of such numerical approach.

1.Introduction

According to IPCC(2014),many extreme weather and climate events have been observed since 1950,including extreme cold/hot temperatures,extreme high sea levels and an increase in heavy precipitation events.Most of these events are related to human activities.Indeed,IPCC(2018)indicated that human-induced warming reached approximately 1°C above pre-industrial levels in 2017.Note that the period 1850-1900 has been used as an approximation of pre-industrial temperatures.There is substantial evidence that human-induced global warming has led to an increase in the frequency,intensity and/or amount of heavy precipitation events at the global scale,as well as an increased risk of drought.It is worth noting that risks associated with increases in drought frequency and magnitude are projected to be substantially larger at 2°C than at 1.5°C.

According to CCR(Caisse Centrale de Réassurance,France)(https://catastrophes-naturelles.ccr.fr/-/002121_sech_2020), on the world scale,the first half of 2020 was the warmest with a mean temperature 1.12°C above the 1981-2010 average.The January-August 2020 period was the warmest ever observed since the beginning of meteorological measurements with 1.72°C above the normal in France.The year 2019 was marked by a winter that was generally drier than average,and by a very dry summer until mid-September.The year 2020 was characterized by a winter that was contrasted in terms of precipitation while being extremely mild.The summer,particularly in July and August,was characterized by signi ficant precipitation de ficits.For the sixth time in the last decade,the geotechnical drought event represents a major event in 2020.

These extreme weather and climate events can greatly affect buildings and geotechnical constructions.In France,the shallowly founded buildings of 216 communes in 17 departments were affected by the extensive drought from 1989 to 1990(Vandangeon,1992).According to the decree of November 1,2005,more than 870 communes were considered affected by the 2003 drought.Examination of the affected buildings showed that the differential settlement due to drying and the differential heave due to wetting were the main reasons for the damage of building structures.The similar phenomena were reported for other countries(Biddle,1983;Driscoll,1983;Holtz,1983;Ravina,1983;Williams and Pidgeon,1983;Gao,1995;Allman et al.,1998).Thereby,it is now widely recognized that the drought-related hazards have important economic,environmental and societal impacts,and deserve particular attention.

Basically,the effect of climate effect on geotechnical constructions represents a vast topic involving various disciplines such as geotechnical engineering,geology,hydrology,meteorology,etc.The research methodology can involve laboratory testing,field monitoring,theoretical development and numerical analysis.In this paper,a fundamental study on soil water evaporation using an environmental chamber was first presented.The analysis of the obtained experimental results allowed soil water evaporation models to be developed,for both sand and clay.Then,some field monitoring results in terms of soil thermo-hydraulic responses to varying meteorological conditions were presented,evidencing the signi ficant soil-atmosphere interaction in earth structures such as embankments.Finally,a numerical analysis was presented,showing the possibility of numerically assessing the thermohydraulic behavior of geotechnical constructions with consideration of meteorological data.

It is worth noting that most results presented in this paper were published by the author.However,all figures were re-plotted and most results were re-organized.Most importantly,the results were re-interpreted,allowing new conclusions to be drawn.

2.Investigation of soil water evaporation using environmental chamber

In order to analyze the soil-atmosphere interaction in geotechnical constructions,it is important to well understand the main mechanisms related to soil water evaporation.Previous studies showed that soil water evaporation is controlled by both atmospheric and soil conditions.The water evaporation rate was found to be greater at lower air relative humidity,higher air temperature,higher wind speed and stronger solar radiation(Cui and Zornberg,2008).The water evaporation rate was found to be dependent on soil water retention capacity,hydraulic conductivity,water content or suction,and dry density,as well as water table level(Wilson et al.,1994;Yanful and Choo,1997;An et al.,2018a).

For clayey soils,cracks can develop upon drying.Several studies focused on the factors affecting the cracking initiation,including soil layer thickness(Nahlawi and Kodikara,2006),temperature(Tang et al.,2010),wetting/drying cycles(Tang et al.,2011),fines content(Yesiller et al.,2000),etc.As the developed cracks modify the interface between atmosphere and soil in terms of morphology and water content distribution,the soil water evaporation rate is expected to be greatly affected.This was con firmed by Cui et al.(2013)through an evaporation experiment using an environmental chamber.

Therefore,to well investigate the soil water evaporation phenomenon,it appears necessary to well control the atmospheric conditions(relative humidity,temperature,wind speed,etc.)and the soil conditions(initial state,water table level,dry density,etc.).For clayey soils,in addition,it is necessary to monitor the soil cracking.Nevertheless,until now,few experiments were carried out with full control of atmospheric conditions and good monitoring of soil response(Wilson et al.,1997;Yanful et al.,2003).In this study,soil evaporation was investigated on both sand and clay using an environmental chamber which was initially developed by Tang et al.(2009),and then modi fied by Song et al.(2013).Analysis of the obtained results allowed the development of soil water models for sand and clay,respectively.

2.1.Environmental chamber

The environmental chamber used for evaporation experimentation is shown in Fig.1.It is an acrylic chamber of 1000 mm long,800 mmwide and 895 mm high.The chamber was equipped with a distributor for diffusing the air which was previously heated to a certain temperature,a collector for collecting the air after soil water evaporation,and a water supply for controlling the water table using a graduated tube.The soil sample was prepared by compaction to a desired dry density and the thickness of the sample was 250-300 mm.In the test,the air temperature and relative humidity were measured at the inlet of the chamber before the heated air was diffused into the chamber and at the outlet of the chamber in the air collector.

Different sensors were used to monitor the atmospheric and soil conditions,including high-capacity tensiometers for matric suction measurement,ThetaProbes for volumetric water content measurement,temperature sensors(PT1000)for soil temperature measurement,T3111 transmitters for air temperature and relative humidity measurements,and thermistors for air temperature measurement.Sensor Pyropen-D was installed at the chamber cover for measuring the soil surface temperature.A camera(Canon EOS400D)was used for monitoring the morphology of the soil surface(in the case of clay).A flowmeter(MAS-3120)was used for air flow rate measurement and an anemometer(Testo 435-2)was used for the wind speed measurement.More details about this environmental chamber can be found in Song et al.(2013).

2.2.Case of sand

Fig.1.Environmental chamber used for evaporation experiment(modi fied from Song et al.,2013).

Fontainebleau sand of 2.64 speci fic gravity was used.In order to investigate the effects of air temperature and air flow rate,two tests(tests 1 and 2)were carried out.The soil sample was prepared by compaction using a manual compactor to a dry unit mass of 1.7 Mg/m3.Prior to the evaporation test,the soil sample was saturated by allowing water in filtration from the bottom.An air temperature of 50°C in the heating tube and an air flow rate of 185 L/min were applied in test 1,while an air temperature of 200°C in the heating tube and an air flow rate of 130 L/min were applied in test 2.Note that the temperature in the heating tube was much higher than the temperature of the air diffused in the chamber.Test 1 lasted 11.5 d,while test 2 lasted 30 d.

Fig.2 shows the air-soil temperature pro files,with some representative selected times(0 d,4 d,8 d,and 11.5 d for test 1;and 0 d,4 d,8 d,12 d,and 30 d for test 2).The highest air temperature was obtained at the elevation of air distributor,which was normal because the air was heated in the heating tube before reaching the distributor.The air temperature was decreasing from the distributor to the chamber cover,indicating the room cooling effect.From the distributor to the soil surface,the air temperature was also decreasing.In the case of low heating tube temperature(test 1),the temperature at the soil surface was the lowest(Fig.2a),while in the case of high heating tube temperature(test 2),this phenomenon was not observed(Fig.2b).The soil temperature was observed to increase over depth under low heating tube temperature(Fig.2a),while it decreased slightly under high heating tube temperature(Fig.2b).

Fig.2.Temperature pro files in air and soil:(a)Test 1 and(b)test 2(data from Song et al.,2014).

As evaporation is an energy-consuming process,theoretically,the position of the lowest temperature corresponded to the evaporative surface below which there was liquid water and above which there was water vapor only.As in test 1,the soil surface temperature kept the lowest in the whole test duration(11.5 d),it could be deduced that the evaporative surface stayed close to the soil surface.On the contrary,in test 2,the lowest soil temperature was observed at the soil surface only at the beginning of the test;over time,the position of the lowest temperature was found to be deeper,sign of the development of deeper evaporative surface.This suggested the formation of a completely dry layer in the area close to the soil surface in the case of high heating tube temperature(test 2).

The air and soil temperatures were found to increase over time.This was because during soil water evaporation,the soil water content decreased and the soil suction increased,leading to a decrease in evaporation rate.As a result,the energy consumed by evaporation declined and more energy was available for heating soil and air.

The pro files of volumetric water content are presented in Fig.3,with some representative selected times(0 d,4 d,8 d,and 11.5 d for test 1;and 0 d,4 d,8 d,12 d,16 d,and 30 d for test 2).Due to evaporation,the volumetric water content was decreasing signi ficantly in both tests.Further examination showed that water loss occurred mainly in the near surface zone.At depths greater than 275 mm,the water loss was negligible,especially in test 1 with lower heating tube temperature and shorter test duration.Comparison of water loss during about 12 d between test 1 and test 2 showed that the water loss in test 2 was slightly larger,showing the effect of temperature on water evaporation.

The pro files of soil suction are shown in Fig.4,with some representative selected times(0.25 d,2 d,4 d,and 6 d for test 1;and 0 d,4 d,8 d,and 16 d for test 2).For test 1,the signi ficant variations occurred in the near surface zone(0-75 mm depth),while for test 2,a large variation in the near surface zone(0-25 mm depth)was followed by a smaller variation in the deeper zone.This suggested a larger in fluence zone in test 2,in agreement with the observation from Fig.3 for the volumetric water content variations.

It appeared that the suctions in the near surface zone were higher in test 1 than in test 2,which would suggest more water loss in this zone in test 1.This was contradictory to the results of volumetric water content shown in Fig.3.A tentative explanation would be related to the soil heterogeneity.Indeed,even though the soil was compacted at the same dry density in both tests,the water flow rates could be different in the two tests,and some fine sand particles would be migrated to the sample surface especially in test 1,leading to higher suctions by evaporation.Note that the sand particle migration is often reported as an internal erosion mechanism(see for instance Chang et al.,2020).This suggested that it should be cautious when further analyzing such data-it would be better to use the mean values of data from different tests.

The data of volumetric water content and suction allowed the water retention curve to be plotted(Fig.5).Note that in the determination,for each soil suction measurement,the corresponding volumetric water content was determined based on the volumetric water content pro files(see Fig.3).A large data scatter could be observed,con firming the sample heterogeneity mentioned previously.Fredlund and Xing(1994)’s model can be applied to describe the water retention property of the tested sand:

whereθwis the volumetric water content(%);θsis the volumetric water content in saturated state(θs=35.6%);θris the residual volumetric water content(θr=4%);ψis the matric suction(kPa);e is the base of natural logarithm(e=2.71828);anda,nandmare the soil parameters(a=40,n=1.4,andm=15).

Fig.3.Volumetric water content pro files:(a)Test 1 and(b)test 2(data from Song et al.,2014).

The air entry value could be estimated at 1-2 kPa.The suction corresponding to the residual volumetric water content was estimated at 40 kPa.Physically,the residual volumetric water content corresponded to a soil state with discontinuous water phase.In other words,the residual state was close to the completely dry state.Thus,the dry layer mentioned previously for test 2 referred to the soil with suction higher than 40 kPa.

The following equation was used to determine the evaporation rate(Aluwihare and Watanabe,2003):

whereEais the actual evaporation rate(mm/d),Ha-outletis the humidity at the outlet(Mg/m3),Ha-inletis the humidity at the inlet(Mg/m3),Qais the air flow rate through the chamber(L/s),ρlis the density of water(Mg/m3),andAis the area of the evaporative surface(m2).

The absolute humidity(Ha)is calculated as follows(Aluwihare and Watanabe,2003):

whereTais the air temperature(K);Ris the gas constant(287.04 J/(kg K));andeais the vapor pressure(Pa),which is calculated as

whereHris the air relative humidity(%);andesatis the saturated vapor pressure(Pa),which is calculated using the following equation(Brutsaert,1988):

Using Eq.(2),the actual evaporation rate was calculated for tests 1 and 2(Fig.6).It was observed that during the first 5 d,the evaporation rate in test 1 was between 1.75 mm/d and 2 mm/d,while the evaporation rate in test 2 was about 2 mm/d,indicating that the evaporation rate under low heating tube temperature(test 1)was lower than that under high heating tube temperature(test 2).Moreover,the evaporation rate in test 2 exhibited three-stage variation:(i)a first nearly constant rate,(ii)a second signi ficant declining rate,and(iii)a third stabilized rate.Owing to the short duration of test 1,only one stage was observed.

Fig.4.Pro files of suction:(a)Test 1 and(b)test 2(data from Song et al.,2014).

Fig.5.Water retention curve determined using data of volumetric water content and suction from tests 1 and 2(data from Song et al.,2014).

As mentioned previously through the observation of soil temperature variations(Fig.2),during the evaporation process,the evaporative surface was going deeper.From the water retention curve,it was found that this evaporative surface corresponded to a suction of about 40 kPa.Thereby,it is possible to determine the dry layer depthdfrom the suction pro file at a given time.In other words,Fig.7 can be considered as representative of soil state at a given time.

By considering the relative humidity values at the soil surfaceh_surfaceand at the evaporative surfaceh_drying,a functionf(d)can be de fined:

Thef(d)functions were determined for tests 1 and 2 using the experimental data and are shown in Table 1.A mean functionf(d)(test 1-2)was also determined considering both the results from tests 1 and 2.

Table 1Functions of f(d)determined from tests 1 and 2.

The suction-related model proposed by Campbell(1985)and Wilson et al.(1997)has been widely used in describing the soil water evaporation rate as follows:

wherehais the relative humidity of air(%);andEpis the potential evaporation(mm/d),which corresponds to the case with water level at the soil surface.

Substituting Eq.(7)into Eq.(8)leads to

The model proposed by Ta(2009)could be used to determineEp:

whereuis the wind speed(m/s)at a reference elevation(taken equal to 50 mm above the soil(or water)surface);ands1ands2are the two constants whose values were determined through free water evaporation experiments under different atmospheric conditions,ands1=0.022 ands2=0.031(see more details in Song et al.,2014).

Using Eqs.(9)and(10),the actual evaporation rate was calculated for tests 1 and 2.The calculated results were compared with the experimental ones in Fig.6.A good agreement was obtained for both tests,showing the performance of the proposed model and the importance of considering the in fluence of the dry layer throughf(d)function.Note that the meanf(d)function was used in the calculation.

2.3.Case of clay

To investigate the water evaporation mechanism for clayey soils,a test was performed on a clay taken from an embankment construction site in Héricourt,France.The clay was air-dried,ground and passed through 2 mm sieve.The soil sample was prepared by compaction using a manual compactor to reach a dry density of 1.4 Mg/m3,which was also the dry density of soil used for the embankment construction in Héricourt(Dong,2013).In order to ensure the soil homogeneity,the compaction was carried out in five layers of 50 mm thick each.The total thickness of the sample was thus 250 mm.More detail can be found in Song et al.(2016).

The actual evaporation rate determined is presented in Fig.8.As for test 2 on sand(Fig.6b),three zones were identi fied:a first zone with a constant value around 2.3 mm/d during the first 15 d,a second zone with a signi ficant decrease to 0.5 mm/d att=45 d,and a third zone with a rate close to zero.The similar phenomenon was observed by Yanful and Choo(1997)on a clayey soil.

Unlike for sand,for the tested clay,signi ficant cracking was observed.In order to quantitatively describe the cracking phenomenon,a surface crack ratio was adopted,which is the ratio of the surface of cracks to the total sample surface.The evolution of the surface crack ratio during the evaporation test is also shown in Fig.8.Three zones were also observed:(i)a first zone with a slight increase untilt=10 d,(ii)a second zone with a rapid increase fromt=10 d tot=25.5 d,and(iii)a third zone with a value around 25.3%.Comparison between the evolutions of actual evaporation rate and surface crack ratio showed that the initiation of cracking led to more water evaporation(the first zone).With further cracking,the evaporation rate was decreasing due to the less available water inside the soil(the second zone).Finally,when the soil water content became extremely low,the evolution of cracks reached a steady state and the evaporation rate approached zero.Further examination showed that the characteristic days separating different zones were not the same for the evaluations of surface crack ratio and actual evaporation,the evolution of actual evaporation curve lagging behind that of the surface crack ratio curve.This suggested that cracking was just one factor affecting the evaporation process,in addition to the factors related to soil conditions(suction and temperature)and atmospheric conditions(air relative humidity,air temperature,wind speed,etc.).To better appreciate the evolution of surface crack ratio,four photos of soil surface at different times are shown in Fig.8.

As cracking can greatly affect the water evaporation for clays,it is of paramount importance to take it into account while developing evaporation models.Compared to the case of sand,the problem is more complex because the interface between the atmosphere and the soil is no longer horizontal-the interface becomes three-dimensional(3D)with complex water content distribution:the water content is expected to be higher in deeper level of cracks with larger distance to the soil surface or soilatmosphere interface.Fig.9 schematizes such evaporation con figuration.

In Fig.9,the 3D cracked surface of clay was assimilated to a horizontal atmosphere-soil interface by introducing two parameters:the surface crack ratioRcand a relative humidity ratiok=h_crack/h_non-crackwhereh_crackis an equivalent relative humidity inside the desiccation cracks andh_non-crackis the relative humidity at the non-cracked surface.Normally,h_non-crackis much lower thanh_crackbecause of the higher water content in cracks.

Fig.6.Comparison of variations of evaporation rate between experiment and model:(a)Test 1 and(b)test 2(data from Song et al.,2018).

Fig.7.Evaporation with a dry layer(modi fied from Song et al.,2018).

An equivalent relative humidity at the soil surface can then be de fined as follows(Ta,2009):

or

Substituting Eq.(12)into Eq.(8)leads to

Eq.(13)allows water evaporation to be calculated with consideration of the effect of cracks.Fig.10 shows the comparison between the calculation and the measurement for the test conducted.A good agreement was obtained,showing the performance of the proposed model.

3.Field monitoring

Several studies were reported,aiming at investigating the soilatmosphere interaction in field conditions(Cui et al.,2008;Bittelli et al.,2012;Smerthurst et al.,2012;Bicalho et al.,2018).It was found that it is the near surface zone that is affected the most by the climate changes.In this zone,the changes in temperature,water content and suction are the most signi ficant,in agreement with the results obtained from the tests using environmental chamber(Figs.2-4).As a result,the thermo-hydro-mechanical(THM)behavior of soil also changes signi ficantly in this zone.In most cases,such climate effect has been accounted for in a simple way by considering the rainfall effect only(Cai and Ugai,2004;Regmi et al.,2017).Bittelli et al.(2012)reported that using rainfall data only is not suf ficient for assessing the landslide process,in particular for shallow clayey soils.It is thus necessary to consider the effects of full atmospheric conditions(solar radiation,air temperature,air relative humidity,wind speed,etc.)on the soil behavior(Hemmati et al.,2012;Cui et al.,2013;An et al.,2018b).In this regard,field monitoring with rich instrumentation of soil and rich measurements of air parameters is essential for further improving the analysis of the THMbehavior of earth structures.In this section,an example of such field monitoring is given,involving an embankment constructed with lime/cement treated soils.

3.1.Embankment and instrumentation

Within an ANR(Agence Nationale de Recherche)project-TerDouest(ANR-07-PCGU-006-10),an experimental embankment was constructed at Héricourt,in the northeast of France where a continental climate dominates with oceanic in fluences.The embankment was 107 m long by 4.9 m high with side slopes of 1:2(vertical:horizontal).The bottom was 25 m wide and the top was 5 m wide.The embankment was divided into two sections,constructed with a silty soil and a clayey soil,respectively.The silty soilwas classi fied as CL(Ip=18),an inorganic clay with low plasticity,while the clayey soil as CH(Ip=45),an inorganic clay with high plasticity.Both soils were treated with lime/cement in different dosages.The embankment consists of 17 layers with the two fill materials compacted at their standard optimum water contents.In total,5280 m3of silty soil,4710 m3of clayey soil,320 t of lime and 162 t of cement were used.

Fig.8.Evolutions of actual evaporation rate and surface crack ratio(data from Song et al.,2016).

The instrumentation layout was symmetrical for the two sections with the two treated fine-grained soils.Different sensors were installed in the embankment for monitoring the water content,suction,temperature and deformation in different positions.Four piezometers were installed to monitor the water table,two on the silty soil side and two on the clayey soil side.A system of runoff measurement was also installed to monitor the runoff from the side slope(An et al.,2017a).A meteorology station was installed on the top surface of embankment to record the meteorological data every 30 min,including solar radiation,precipitation,atmospheric pressure,wind speed and direction,air temperature and air relative humidity,at 0.5 m and 1.5 m above the top surface.In this section,only the soil suction and water content as well as the runoff are emphasized.

The Watermark soil suction sensors were used to monitor the matric suction.The working range of such sensor is 0-250 kPa.The reading frequency was once every 1-2 d.The time domain re flectometry(TDR)method was used to monitor the soil volumetric water content,together with the soil temperature.The reading frequency was every 3 h.More details can be found in Froumentin(2012).

3.2.Field monitoring data

Rich data were obtained from the field monitoring.However,for the purpose of illustrating the soil responses to the changes in meteorological conditions,only recordings of precipitation,volumetric water content,suction and runoff were presented here,as mentioned previously.

The variations of volumetric water content at points 1-4 are presented in Fig.11,together with the precipitation,for the monitoring period from July 2010 to April 2014.Points 1 and 4 were situated 0.25 m far from the side slope;point 2 was situated 0.75 m far from the side slope;and point 3 was situated 2.4 m far from the side slope.Note that the same instrumentation was adopted for the silty soil and the clayey parts.In other words,there were the same positions 1-4 for the volumetric water content measurement in both parts.As the soils were treated with different binders(lime or cement)and with different dosages,for the purpose of comparison,only the points with the same treatment have been selected for analysis,as in Boussa fir et al.(2018).

The case of silty soil treated with 2%lime is depicted in Fig.11a.In terms of rainfall,relatively frequent and high intensity events occurred in the period from July 2010 to November 2011 and in the period from May to October 2012.Other periods seemed quite dry.Note,however,the extreme event in late March 2014 with an intensity as high as 19 mm(in 30 min).The variations of volumetric water content showed clearly that the further the point from the side surface,the smaller the variation-the variation at point 3(2.4 m from the side slope)was much smaller than those at points 1 and 4.This indicated the signi ficant effect of hydraulic conductivity on the water transfer inside the soil.The low permeability of the treated compacted silty soil(1×10-9m/s)made the water transfer quite slow in the processes of in filtration during rainfalls and capillary movement during evaporation,limiting the variations of volumetric water content in the zone far from the side slope.Comparison between points 1 and 4 showed that the volumetric water content variation depended on elevation.Indeed,even though both points were situated at the same distance from the side slope(0.25 m),the variation at point 4 which had a lower elevation showed a higher mean value but lower fluctuation(or amplitude),suggesting a more intense in fluence of atmosphere in the higher elevation zone.This can be explained by the more contribution to the lower elevation zone by runoff during rainfalls and by capillary movement during evaporation.Two distinct decrease periods were observed at the two points,one from late March to early July 2011 and another from late May to early August 2013.Interestingly,in these two periods,there were no particularly less frequent and intensive rainfall events.This suggests a relatively signi ficant water evaporation effect.Thereby,while analyzing the variation of water content inside a soil,considering the rainfalls only is not suf ficient.It is of paramount importance to consider the full soil-atmosphere interaction.This was con firmed by the signi ficant decrease of volumetric water content after the extreme event in late March 2014.Compared to the initial value(23.4%),the mean values of volumetric water content at all points signi ficantly increased,suggesting that the monitoring period rather involved wet seasons.

The case of clayey soil treated with 4%lime is depicted in Fig.11b,at points 1 and 2.The similar observation could be made:(i)the larger variation amplitude in the zone closer to the side slope(point 1),showing the more signi ficant effect of atmospheric condition in the near surface zone;(ii)the higher mean volumetric water contents than the initial value(33.5%),indicating the wet monitoring period involved;and(iii)the lack of correspondence between the periods with signi ficant decrease in volumetric water content and the dry period without rainfall events,showing the importance of considering both water in filtration and evaporation while assessing the variations of water content in soil.

Comparison between the case of silty soil(Fig.11a)and the case of clayey soil(Fig.11b)showed that the variation amplitudes are larger in the case of silty soil.This could be explained by the higher hydraulic conductivity of the silty soil.The mean values of volumetric water content of the silty soil were lower than those of clayey soil,suggesting a higher water retention capacity and a lower water evaporation for the clayey soil.

Fig.9.Sketch of water evaporation from cracked clay(modi fied from Song and Cui,2020).

Fig.10.Comparison between the measured and predicted actual evaporation rates.

Fig.12 shows the variations of suction at three different distances from the side surface(0.25 m,0.5 m and 0.75 m),for the silty soil treated with 2%lime(Fig.12a)and the clayey soil treated with 4%lime(Fig.12b),in the period from July 2010 to October 2014(8 months longer than that for the volumetric water content recording).For the silty soil(Fig.12a),four distinct high suction periods were observed:(i)from April to July 2011,(ii)in late February 2012,(iii)from June to September 2013,and(iv)from June to July 2014.By referring to the variations of volumetric water content(Fig.11a),it appeared that the first and the third high suction periods coincided with the two periods with signi ficant volumetric water content decreases.However,the second short high suction period did not seem to correspond to the volumetric water content variations.Further examination showed that this short high suction period only involved the measurement at the position 0.5 m far from the side slope,the suctions at the other two positions(0.25 m and 0.75 m)being much lower.It could thus be inferred that the high suction measured in this period at position 0.5 m might be related to a technical problem(likely of electronical nature)and should be ignored in further analysis.Because the fourth high suction period occurred without volumetric water content recordings,no correspondence could be attempted.In addition to these four distinct high suction periods,a number of lower suction periods could be identi fied,such as the period at the beginning of the recording(early July 2010)and the period from late July to early September 2012,which indicated the soil state with relatively higher water content(the suction peaks were relatively less pronounced).

For the clayey soil(Fig.12b),six distinct high suction periods could be identi fied.The two low suction periods mentioned above for the silty soil became more distinct for the clayey soil,suggesting a higher sensitivity of clay to climate change than silt in terms of suction variation.The other four high suction periods corresponded to those identi fied on the silty soil.The suspected technical problem occurred also for the clayey soil,leading to the high suction period in late February 2012 with the recordings at position 0.25 m.

Comparison between the silty soil(Fig.12a)and the clayey soil(Fig.12b)showed that the variations of suction were more signi ficant in the silty soil.This con firmed that the hydraulic conductivity is an important factor in the response of soil to the changes in atmospheric conditions-the silty soil had a higher hydraulic conductivity and its suction changed more under the effects of rainfall/evaporation.

During the rainfall events,the changes in soil water content or suction are governed by the quantity of water in filtrating into the soil.It is thus essential to determine the water in filtration based on the measured precipitation.For this purpose,it is necessary to have the measurements of runoff.Basically,surface runoff occurs when the rainfall rate exceeds the in filtration capacity of the soil(Beven,2012).As direct measurement of runoff is dif ficult,indirect determination using correlations is often applied.One of the most widely used correlation is the Natural Resources Conservation Service(NRCS)runoff method(SCS,1985),which allows runoff to be estimated from rainfall.Following this method,the measured runoff on the embankment slope is correlated with rainfall using the hourly recordings(mm/h),as shown in Fig.13.It is interesting to note that the runoff remained quite low(＜0.1 mm/h)when the precipitation was lower than about 11 mm/h,and became signi ficant beyond 11 mm/h precipitation.Physically,the low runoff means that the water in filtration into the soil was dominating,while the signi ficant runoff means that the runoff became dominant facing the water in filtration into the soil.Interestingly,a linear correlation was established between runoff and precipitation.

It is worth noting that the turning point at 11 mm/h precipitation is valid only for the embankment considered.When the geometry of embankment changes,this value is expected to be different.

4.Numerical modeling of soil-atmosphere interaction

The results obtained from the tests in the environmental chamber and from the field measurements showed the necessity of considering the full soil-atmosphere interaction in the assessment of soil suction and water content changes.In this section,such an approach is presented through modeling the soil responses to the changes in atmospheric conditions for the Héricourt embankment presented in the previous section.As the embankment was constructed with treated soils,the vertical deformation was found quite small(＜0.3%).Thus,the volume change was not necessary to be accounted for,and only the coupled heat flow and water(liquid and vapor)flow needed to be taken into account.Note that when the volume change becomes a concern,a mechanical constitutive model needs to be incorporated(see for instance Hemmati et al.,2012;Cui et al.,2013).

Fig.11.Variations of volumetric water content and precipitation in the period from July 2010 to April 2014:(a)Silty soil with 2%lime and(b)clayey soil with 4%lime(data from Boussa fir et al.,2018).

The finite element method was adopted for the fully coupled thermo-hydraulic analysis.The boundary conditions at the soilatmosphere interface were first determined using the water evaporation model as those presented in Section 1,based on the meteorological data and the initially estimated soil surface temperature and suction.Then,the numerical tool with implementation of the fully coupled thermo-hydraulic model was used for the calculation of soil temperature,water content and suction.The calculated temperature and soil suction on the boundary at the soil-atmosphere interface were used to de fine the new boundary condition for the next time step.

4.1.Soil-atmosphere interaction model

Two processes are involved in the soil-atmosphere interaction:water transfer and heat transfer.To determine the water transfer,the mass balance at the soil surface was used(without vegetation,see An et al.,2017b):

wherePis the precipitation(m/s),Roffis the runoff(m/s),andInfis the in filtration(m/s).

Infwas used for de fining the boundary condition in terms of water transfer.It can be positive(in filtration)or negative(evaporation).PandRoffwere determined using direct field monitoring data.Eqs.(8)and(10)were used for the calculation ofEa.

To determine the heat transfer,the energy balance at the soil surface was used(An et al.,2017b):

whereRnis the net radiation flux(W/m2),Gis the soil heat flux(W/m2),LEis the latent heat flux(W/m2),andHis the sensible heat flux(W/m2).

Fig.12.Variations of suction in the period from July 2010 to October 2014:(a)Silty soil with 2%lime and(b)clayey soil with 4%lime(data from Boussa fir et al.,2018).

Gwas used for de fining the boundary condition in terms of heat transfer.It can be positive(soil heating)or negative(soil cooling).Rnwas calculated from the measured total radiation flux using the expression proposed by Evett et al.(2011):

whereRsiis the solar irradiance at the surface(W/m2);Rsois the solar radiation in the case of clear sky(W/m2);σis the Stefan-Boltzmann constant(σ=5.67×10-8W/(m2K4));acandbcare the cloud factors(ac=1.35 andbc=-0.35);a1andb1are the emissivity factors(a1=0.035 andb1=-0.14);αis the soil surface albedo(αcan be taken equal to 0.23 for the soils of Héricourt embankment);andedis the mean daily saturated vapor pressure,which was calculated from the mean daily dew point temperatureTd(°C)(Evett et al.,2011):

In the case of clear sky,the solar radiationRso(W/m2)is calculated as(Evett et al.,2011):

whereELmslis the elevation above mean sea level(m),andRsais the extraterrestrial solar radiation(W/m2)and was calculated by(Evett et al.,2011):

whereGscis the solar constant(0.08202 MJ/(m2min));dris the relative earth-sun distance;ϕis the latitude(m);δis the solar declination;w1andw2are the solar time angles(radian)at the beginning and end of the considered period,respectively,and were calculated by(Evett et al.,2011):

wherewis the solar time angle(radian)at the center of the period,andτis the length of the considered period(h).

The latent heat fluxLE(W/m2)and sensible heat fluxH(W/m2)were respectively calculated by(Blight,1997):

whereLvis the latent heat of water vaporization(J/kg),ρais the air density(kg/m3),cpis the speci fic heat of air(J/(kg K)),KHis the eddy diffusivity for heat through air(m2/s),andzis the elevation(m).

4.2.Coupled thermo-hydraulic model

The mass transfer of water is the sum of liquid flow and vapor flow:

Fig.13.Runoff measured on the side slope versus precipitation(data from An et al.,2017a).

where q is the water flow density(kg/(s m2)),qlis the liquid flow density(kg/(s m2)),and qvis the vapor flow density(kg/(s m2)).Darcy’s law was used to describe the non-isothermal liquid flow:

whereKis the unsaturated hydraulic conductivity(m/s),andφis the hydraulic head(m).

Fick’s law was used to describe the vapor flow by diffusion(Philip and De Vries,1957):

whereDatmis the molecular diffusivity of vapor in the air(m2/s),εis the tortuosity of soil,βis the cross-sectional area of soil,andρvis the density of vapor(kg/m3)calculated by Philip and De Vries(1957):

whereρ0is the density of saturated water vapor(kg/m3),gis the gravitational acceleration(m/s2),Mwis the molar mass of water molecule(kg/mol),andTis the absolute temperature(K).

Fourier’s law was used to describe the transferred heat flux through soil-atmosphere interface Q(W/m2):

whereλis the soil thermal conductivity(W/(m K)).

By substituting the fundamental equations(Eqs.(24)and(28))into the conservation equations for water and heat,respectively,the governing equations can be obtained.More details can be found in An et al.(2017b).

4.3.Simulating the thermo-hydraulic responses of Héricourt embankment

The proposed thermo-hydraulic model was implemented in the FreeFem++code(Hecht,2012)for analyzing the thermo-hydraulic behavior of Héricourt embankment.Examination of all data in the monitoring period showed that only the recordings during July 6-26,2011,are regular and complete.Thus,this short period was selected for the numerical analysis.The soil temperature and volumetric water content values on July 6,2011 at 14:42:52 were taken to de fine the initial conditions.The meteorological data were used to determine the in filtration(positive or negative)and the soil heat flux(positive or negative)using the soil-atmosphere interaction model.The obtained values were used to de fine the boundary conditions at the soil-atmosphere interface.A water table depth of 5 m below the ground surface was considered,where the soil suction was taken equal to zero.As only the temperature at the base of the embankment was available,this temperature was used to de fine the thermal bottom boundary condition by assuming that the temperature changed a little from the base of the embankment to the water table.

The silty soil part of the embankment was selected for the analysis.For simplicity,this part of embankment was considered as homogeneous.As proposed by De Vries(1963)and Cui et al.(2005),the soil thermal conductivityλ(W/(m K))was assumed to vary linearly with the volumetric water contentθw:

The soil water retention curve was determined in the laboratory using soil samples taken from the Héricourt site.van Genuchten(1980)’s model was used to fit the curve:

whereSeis the effective degree of saturation;θsis the saturated volumetric water content(θs=0.4);θris the residual volumetric water content(θr=0.004);anda,mandnare the fitting parameters(a=0.003 kPa-1,m=0.18,andn=1.8).

The hydraulic conductivity was also estimated from van Genuchten(1980)’s model:

whereKsis the soil saturated hydraulic conductivity(Ks=1×10-9m/s according to the laboratory measurement),andm1is a fitting parameter(m1=0.5).

Fig.15.Comparison of temperature changes between calculation and measurement in July 2011(data from An et al.,2017b).

Fig.16.Comparison of changes in volumetric water content between calculation and measurement in July 2011(data from An et al.,2017b).

The three monitoring points 1,3 and 4 presented in Fig.11a were considered.As a cover layer of about 0.2 mwas put on the side slope after the construction of the core part of embankment by compaction,this layer was considered in the numerical analysis(Fig.14).It was assumed that this layer was homogeneous.As it was not compacted,a relatively lower thermal conductivity of 0.25 W/(m K)and a relatively higher hydraulic conductivity of 1×10-8m/s were considered.

As regular recordings of every 3 h were performed for soil temperature and volumetric water content,while one value was recorded every 1-2 d for the soil suction,in order to well reveal the variations of soil thermo-hydraulic response during day and night times,only the soil temperature and volumetric water content were selected for the comparative analysis.

Fig.15 shows the comparisons of the variations of temperature between calculation and measurement at points 1,3 and 4.For clarity,a number of skip points equal to 5 were adopted for representing the measured results.Note that for some periods,the data were not available.It appeared that the proposed numerical approach allowed the soil temperature to be reasonably well predicted,even though a slightly higher temperature was calculated for points 1 and 4 in the first 2 d.Moreover,the measurement showed that the temperature at point 1 was slightly higher than that at point 4,and this was well reproduced by the numerical calculation.Note that the lower temperature at point 3 was due to its further distance from the side slope.The slightly higher temperature at point 1 than at point 4 was to be related to the slower water evaporation at point 1.Indeed,as the volumetric water content was lower at point 1(see Fig.11a),the corresponding suction must be higher.Thus,the water evaporation rate must be slower,as described in the first section.

Fig.16 shows the comparisons of the variations of volumetric water content between calculation and measurement.Again,for clarity,a number of skip points equal to 5 were adopted for representing the available measured results.For point 1,an overall good agreement was obtained except the beginning where a slightly lower value was predicted and the end where a slightly higher value was predicted.For point 3,the variation was small during the whole period and the calculation agreed well with the measurement.For point 4,slightly lower values were given by the calculation for the first 2 d and slightly higher values were predicted for the period from July 18 to the end.On the whole,all variations of volumetric water content were reasonably well predicted.

5.Conclusions

Facing the global warming,the effect of climate change on the earth structures has become an increasing concern in the community of geotechnical engineering.To address this concern,it appeared important to better understand the fundamental mechanisms involved in water evaporation process for different soils by performing laboratory tests,to upscale the laboratory condition to the field condition by performing field monitoring,and to develop numerical tools for analyzing and predicting the long-term behavior of earth structures by considering the soil-atmosphere interaction.This paper aims at illustrating such methodology through some representative studies.

For bare soils,the water evaporation process is governed by the atmospheric conditions and soil hydraulic conditions.Evaporation is enhanced by low air relative humidity,high air temperature,high wind speed,and low soil suction or high water content.This evaporation process can be well described by suction-related models such as the model proposed by Campbell(1985)and Wilson et al.(1997).However,depending on the soil involved,different physical phenomena should be taken into account.For sandy soils,as their water retention capacity is relatively low and a suction lower than 100 kPa may make them dry,the evaporative surface is not constantly at the ground surface-it becomes deeper and deeper upon drying.In that case,liquid water transfer by capillary effect in the soil below the evaporative surface and vapor transfer by diffusion in the dry layer above the evaporative surface should be accounted for while analyzing the evaporation process.This can be done by incorporating a function representing the difference of relative humidity between the evaporative surface and the ground surface into the suction-related models.For clayey soils,cracks occur in general upon drying.This changes the onedimensional(1D)evaporation problem to a 3D problem-the evaporative surface is no longer horizontal because cracks are involved with varying water content in them.To account for the effect of cracks in the evaporation process,two parameters can be introduced:a surface crack ratio which represents the ratio of crack surface to the total surface,and a relative humidity ratio which represents the ratio of equivalent relative humidity of cracks to the humidity of the non-cracked part.Comparisons between the models and measurements showed the relevance of such approaches.

The instrumentation for field monitoring should cover both the soil and the atmosphere,with emphasis put on the near surface zone of soil where the soil-atmosphere interaction is expected to be the most signi ficant.The common method of considering weather effect on soil behavior is to account for the precipitation effect.The field data showed that for a slope like the side slope of an embankment,this is far from being satisfactory.First,it is important to consider the runoff in the estimation of water in filtration when the rainfall is higher than a threshold value(about 11 mm/h for the Héricourt embankment).Second,there is no always a satisfactory correlation between the precipitation and soil water content change or soil suction change because of the effect of runoff and the water movement inside the soil under the effect of suction.In the zone far from the soil surface,the changes in water content or suction are small,indicating a limited effect of atmosphere.In the zone near the soil surface,the climate effect is signi ficant.On the side slope of embankment,the positions with higher elevations have a lower water content,higher suction and higher temperature.The lower water content and higher suction can be explained by the contributions of runoff and the water movement inside the soil by capillary effect,while the higher temperature can be explained by the lower evaporation due to the higher suction.Comparison between the silty soil and the clayey soil showed that the variations of water content and suction are smaller in the clayey soil,suggesting the importance of the soil hydraulic conductivity in soilatmosphere interaction-the higher the hydraulic conductivity,the more signi ficant the soil-atmosphere interaction.

By considering the water(liquid and vapor)flows and the heat flow in a coupled thermo-hydraulic model for the soil and a suction-related soil-atmosphere interaction model for the soilatmosphere interface,the soil thermo-hydraulic behavior can be well described.Indeed,comparison between calculations and measurements showed that such a numerical approach can well predict the variations of soil temperature,water content and suction using the meteorological data,in particular the smaller variations in the zone far from the soil-atmosphere interface,the higher temperature and lower water content at higher elevations on the side slope of embankment.It is worth noting that the developed numerical approach was successfully applied in analyzing another embankment constructed with lime-treated silty soil(An et al.,2018b)and in evaluating the long-term thermo-hydraulic behavior of a future site for the disposal of short-lived low and intermediate level nuclear waste(An et al.,2020).Nevertheless,it should be mentioned that in the analyzed embankment case,no retained water was considered.The presence of retained water can greatly change the water content and suction of soil,thus the soilatmosphere interaction.Moreover,as the analyzed embankment was not quite high(4.9 m),it was not necessary to distinguish the sunny and non-sunny slope sides.In case of much higher embankment,the solar radiation can be signi ficantly different for the two kinds of slope sides,and the consideration of such difference in the numerical analysis becomes necessary.

It is also worth noting that this paper did not address the effect of vegetation which is an important issue to deal with while investigating the climate effect on earth structures.The vegetation can not only greatly affect the soil water evaporation through interception and transpiration,but also modify the soil THM behavior through the development of root architecture.Furthermore,the volume change behavior was not accounted for in the numerical analysis because of the negligible vertical deformation recorded for the studied embankment with treated soils.In some cases,this mechanical aspect can be important,in particular for clayey soils.To incorporate such mechanical aspect,it is necessary to use a fully coupled THM model for soils(see for instance Hemmati et al.,2012),together with a soil-atmosphere interaction model for the soil-atmosphere interface.

Declarationofcompetinginterest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to in fluence the work reported in this paper.

Acknowledgments

The support from the French Research Agency(ANR)within the Project TerDOUEST ANR-07-PCGU-006-10 is greatly acknowledged.The author is also grateful to Dr Weikang Song,Dr Ni An and Dr Zi Ying for their help in improving the illustrations.

Notation

AArea of the evaporative surface(m2)

a,m,nFitting parameters for soil water retention curves

a1Emissivity factor(0.035)

acCloud factor(1.35)

b1Emissivity factor(-0.14)

bcCloud factor(-0.35)

cpSpeci fic heat of air(J/(kg K))

dDry layer depth(m)

DatmMolecular diffusivity of vapor in the air(m2/s)

eBase of natural logarithm

eaVapor pressure(Pa)

edMean daily saturated vapor pressure(Pa)

esatSaturated vapor pressure(Pa)

EaActual evaporation rate(mm/d or m/s)

EpPotential evaporation rate(mm/d)

ELmslElevation above mean sea level(m)

gGravitational acceleration(m/s2)

GSoil heat flux(W/m2)

HSensible heat flux(W/m2)

HaAbsolute humidity(Mg/m3)

HrAir relative humidity(%)

Ha-outletAbsolute humidity at the outlet of the chamber(Mg/m3)

Ha-inletAbsolute humidity at the inlet of the chamber(Mg/m3)

haRelative humidity value of air(%)

h_crackEquivalent relative humidity inside the desiccation cracks(%)

h_non-crackRelative humidity at the non-cracked surface(%)

h_surfaceRelative humidity value at the soil surface(%)

h_dryingRelative humidity value at the evaporative surface(%)

InfIn filtration(m/s)

kRelative humidity ratiohcrack／hnon-crack

KUnsaturated hydraulic conductivity(m/s)

KHEddy diffusivity for heat through air(m2/s)

KsSaturated hydraulic conductivity(m/s)

LELatent heat flux(W/m2)

LvLatent heat of water vaporization(J/kg)

MwMolar mass of water molecule(kg/mol)

m1Fitting parameter

PPrecipitation(m/s)

q Water flow desnity(kg/(s m2))

qlLiquid flow desnity(kg/(s m2))

qvVapor flow desnity(kg/(s m2))

Q Transferred heat flux through soil-atmosphere interface(W/m2)

QaAir flow rate through the chamber(L/s)

RGas constant(J/(kg K)or J/(mol K))

RcSurface crack ratio

RnNet radiation flux(W/m2)

RoffRunoff(m/s)

RsiSolar irradiance at the surface(W/m2)

RsaExtraterrestrial solar radiation(W/m2)

RsoSolar radiation in the case of clear sky(W/m2)

SeEffective degree of saturation

s1,s2Empirical constants

TAbsolute temperature(K)

TaAir temperature(K or°C)

TdMean daily dew point temperature(°C)

uWind speed at a reference elevation(m/s)

wSolar time angle at the center of the considered period(radian)

w1Solar time angle at the beginning of the considered period(radian)

w2Solar time angle at the end of the considered period(radian)

zElevation(m)

α Soil surface albedo(0.23)

β Cross-sectional area of soil

δ Solar declination

ε Tortuosity of soil

θrResidual volumetric water content(%)

θsVolumetric water content in saturated state(%)

θwVolumetric water content(%)

λ Soil thermal conductivity(W/(m K))

ρ0Density of saturated water vapor(kg/m3)

ρaDensity of air(kg/m3)

ρvDensity of vapor(kg/m3)

ρlDensity of liquid(Mg/m3or kg/m3)

σ Stefan-Boltzmann constant(5.67×10-8W/(m2K4))

τ Length of the considered period(h)

φ Hydraulic head(m)

ψ Matric suction(kPa)

ϕ Latitude(m)

GscSolar constant

drRelative earth-sun distance