M Le Hir

Y Wyart

G Georges

L Siozade

M Sergent

P Moulin

Ultrafiltration membrane fouling characterization by multi-scale methodology

Keywords: Characterization, multi-scale, fouling, nanoparticles, fluorescence

Introduction

Membrane characterization by tracer retention is widely used to determine the retention potential and integrity of membranes by quantification of tracers passed through this one. The use of fluorescent tracers makes this retention quantification quick and reliable. Moreover, additional information like the location of retained tracers [START_REF] Wu | Characterization of Ultrafiltration Membranes Fouled by Quantum Dots by Confocal Laser Scanning Microscopy[END_REF][START_REF] Marroquin | Location and Quantification of Biological Foulants in a Wet Membrane Structure by Cross-Sectional Confocal Laser Scanning Microscopy[END_REF]) and the composition of membrane fouling, with the use of complementary techniques, can be obtained, allowing to characterize membrane morphology and fouling. This study presents a multiscale characterization methodology to clearly determine the retention potential of ultrafiltration membrane regarding fluorescent nanoparticles (NPs) with calibrated size. Location of each NPs size in or/and on the membrane after filtration run allows to identify with precision the penetration profiles of NPs in membrane structure (skin and support) and to determine fouling mechanism(s) which can be correlated with macroscopic flux data observed during the filtration.

Material and Methods

Ultrafiltration (UF) membranes tested are PES multichannel hollow fiber with nominal pore size of 20 nm. Ideal suspensions of fluorescent calibrated NPs with size of 1,5, 10 and 100 nm were filtered. The size of the NPs used was chosen regarding the membrane pore size (smaller, close to and bigger than the membrane pore size). Zeta potential of NPs suspension and membrane were measured to ensure that no aggregation or adsorption occurs during filtration. The efficiency of the membrane was determined by NPs concentrations measurement of filtration fluxes by fluorimetry and Nanoparticle Tracking Analysis (NTA). The use of Nanosight NS300 (using NTA) allowed to obtain concentrations in number of particles per mL allowing to estimate the number of NPs remaining blocked in and/or on the membrane. Location of NPs fluorescent signal in and/or on the membrane was visualized by confocal laser scanning microscopy (CLSM) and penetration profiles of NPs have been determined. In the case of a NPs deposit on membrane surface, its thickness was measured with scanning electron microscopy (SEM) and CLSM. The use of NPs with different fluorescent characteristics allowed to distinguish each size of NPs in membrane material. Location of NPs fluorescent signals have been correlated with macroscopic permeate flux data to validate fouling mechanism(s) operating during the filtration experiment. Influence of Transmembrane Pressure (TMP) and Volume Concentration Factor (VCF) on NPs retention, fouling establishment and location was treated by experimental design.

Results and Discussion

Ultrafiltration of ideal monodisperse and polydisperse NPs suspensions was conducted. Membrane shows interesting retention rates (≥ 50% for all NPs size and for each condition tested) and a low recovery in the retentate reflecting a retention of NP by blocking at the membrane. Fouled membranes were imaged by CLSM on their section allowing to enjoy the lateral resolution of CSLM (about 500 nm), better than the axial one used in a previous study [START_REF] Wu | Characterization of Ultrafiltration Membranes Fouled by Quantum Dots by Confocal Laser Scanning Microscopy[END_REF]. Analysis in reflective mode allowed to visualize membrane and another in fluorescent mode allowed to identify NPs. MCBL visualization showed that 100 nm NPs were essentially retained on membrane surface, 10 nm NPs were stopped in the membrane skin (thickness ≈ 50 µm) and on membrane surface and 1.5 nm NPs are located in the entire material (Figure 1). These locations reflect respectively fouling models of cake filtration (100 nm), standard blocking followed by cake filtration (10 nm) and standard blocking (1.5 nm) evaluated with experimental data of permeate flux.

Conclusions

Multi-scale characterization techniques were used to determine UF membrane efficiency towards NPs, fouling establishment and its location in or/and on the membrane. The complementary techniques used allow to obtain concordant results of interest showing the reliability of characterization realized. It is the first time that the penetration profile of NP in UF membrane was determined with high precision (skin vs. support). Moreover, the number of NP in the permeate (i.e. in the drinking water) was evaluated against to answer at a societal point of view.

Figure 1 .

 1 Figure 1. Penetration profile of a) 1.5 nm NP, b) 10 nm NP and c) 100 nm NP in membrane material

Acknowledgment

"The project leading to this publication has received funding from Excellence Initiative of Aix-Marseille University -A*MIDEX, a French "Investissements d'Avenir" programme. It has been carried out in the framework of the Labex MEC." The authors are grateful to Aquasource for providing the new generation membranes used during this work.


I Schimmelpfennig

J M Schaefer

N A Kçar

S Ivy-O Chs

R C F Inkel

C Schlüchter

Holocene glacier culminations in the Western Alps and their hemispheric relevance

The natural variability of Holocene climate defines the baseline to assess ongoing climate change. Greenland ice-core records indicate warming superimposed by abrupt climate oscillations in the early Holocene, followed by a general cooling trend throughout the mid-and late Holocene that culminated during the Little Ice Age (LIA). Tropical precipitation changes correlate with these patterns throughout the Holocene.

Here we use mountain glaciers in the Alps to reconstruct the regional Holocene climate evolution and to test for a link between mid-latitude, North Atlantic and tropical climate. Our precise 10 Be chronology from Tsidjiore Nouve Glacier, western Swiss Alps, indicates a glacier culmination during the earliest Holocene ~11.4 k.y. ago, likely related to the Preboreal Oscillation. Based on our data, no Holocene glacier advance of similar amplitude occurred until ~3.8 k.y. ago, when the glacier reached LIA limits. 10 Be ages between 500 and 170 yr correspond to the LIA, while the youngest 10 Be ages overlap with the historically recorded post-

Page 2 of 15 LIA glacier positions. Integrating our data with existing records, we propose a hemispheric climate link between the Alps, North Atlantic temperature and tropical precipitation patterns for the Holocene, supporting the concept of a pervasive climate driver. These findings from northern mid-latitudes are consistent with the hypothesis formulated for the tropics that the earth's thermal equator, responding to North Atlantic temperature changes, might have migrated southward throughout the Holocene, reaching the southern turning point towards the end of the LIA.

I N T R O D U C T I O N

While there is increasing evidence that the climate during the current warm period, the Holocene, has been less stable than originally thought (e.g., [START_REF] Haug | Southward migration of the Intertropical Convergence Zone through the Holocene[END_REF], robust terrestrial data quantifying Holocene climate changes remain sparse mostly due to their moderate amplitude compared to the dramatic climate fluctuations during ice-ages. The Greenland temperature record indicates rapid warming after the Younger Dryas (YD, ~12,900 11,700 years ago), interrupted by early Holocene abrupt cold events such as the Preboreal Oscillation (PBO; [START_REF] Rasmussen | Early Holocene climate oscillations recorded in three Greenland ice cores[END_REF]. Warm early Holocene conditions in the North Atlantic area are followed by a long-term cooling [START_REF] Cuffey | Temperature, accumulation, and ice sheet elevation in central Greenland through the last deglacial transition[END_REF], culminating in the Little Ice Age (LIA, 1300 1850 common era, CE, e.g. [START_REF] Holzhauser | Glacier and lake-level variations in westcentral Europe over the last 3500 years[END_REF]. Hydrological changes in the tropics correlate with long-and short-term climate events in the North Atlantic [START_REF] Haug | Southward migration of the Intertropical Convergence Zone through the Holocene[END_REF][START_REF] Sachs | Southward movement of the Pacific intertropical convergence zone[END_REF], implying shifts in the Intertropical Convergence Zone (ITCZ) as a result of the thermal equator adjusting to northern temperature changes. If this scenario is correct, the mid-latitudes should have responded to polar cold and warm spells, but evidence remains controversial.

Here we use glacier fluctuations in the Western Alps to reconstruct regional Holocene climate patterns and to test the link between polar climate changes and mid-latitude mountain Publisher: GSA Journal: GEOL: Geology Article ID: G33169

Page 3 of 15 glaciers, which are highly sensitive to regional climate changes and in particular, to summer temperature variations [START_REF] Oerlemans | Extracting a Climate Signal from 169[END_REF]. We produced a detailed moraine map and a precise 10 Be chronology for one of the best resolved Holocene lateral moraine sequences in the Alps (Fig. 1); it was deposited by the Tsidjiore Nouve Glacier in the Valais, Switzerland, one of the most climate sensitive glaciers in the Swiss Alps (Röthlisberger, 1976). The overall goal of this study is to test for glacier fluctuations related to the early Holocene abrupt climate changes and to add constraints of climate's evolution during the mid-and late Holocene, including the LIA period. In the Alps, Holocene moraine chronologies are scarce [START_REF] Joerin | Multicentury glacier fluctuations in the Swiss Alps during the Holocene: The Holocene[END_REF], as the LIA advances wiped out the older moraines in most places (Röthlisberger and Schneebeli, 1979). A few pioneering surface exposure dating studies from the Swiss and Austrian Alps [START_REF] Ivy-Ochs | Latest Pleistocene and Holocene glacier variations in the European Alps[END_REF][START_REF] Schindelwig | Lateglacial and early Holocene dynamics of adjacent valley glaciers in the Western Swiss Alps[END_REF] report glacier advances exceeding the LIA position during the early Holocene, and correlate those advances to North Atlantic cold pulses.

G E O M O RP H I C SE T T I N G

The well-preserved sequence of Holocene lateral moraine ridges at Tsidjiore Nouve Glacier (Fig. 1) is composed of i) an outermost pre-LIA moraine: an up to 20 m high and 1 km long well-defined moraine ridge. At least five recessional sub-ridges, less than 1 m high, line its ice-proximal flank. Large boulders are embedded on all ridges; ii) a :

An up to 60 m high and 3 km long moraine wall, formed by at least five sub-ridges with few boulders. Such massive composite moraines fringing the non-vegetated glacier forefield are Röthlisberger and Schneebeli, 1979); iii) an outer and an inner post-LIA moraine embedded boulders. These two moraines are assigned to 1920CE and 1991CE, respectively (Abbühl et al., 2009;Fig. 1), based on direct measurements of the glacier length changes of numerous glaciers in this region (Glaciological reports, 1881 2009, Fig. DR1).

M E T H O DS

We collected 29 surface samples from large boulders protruding >0.80 m from the lateral moraine ridges of Tsidjiore Nouve Glacier (Figs. 1 andDR2). All samples were processed for 10 Be measurements at the Lamont-Doherty Earth Observatory Cosmogenic Nuclide Laboratory (method in [START_REF] Schaefer | High-frequency Holocene glacier fluctuations in New Zealand differ from the northern signature[END_REF]; http://www.ldeo.columbia.edu/cosmo/). 10 Be/ 9 Be ratios were measured at the Center for Accelerator Mass Spectrometry, Lawrence Livermore National Laboratory. To calculate the surface exposure ages, we used the 10 Be production rate calibrated in northeastern North America (Balco et al., 2009), whose value is confirmed by several recent calibration experiments elsewhere in the Northern Hemisphere (e.g., [START_REF] Fenton | Regional 10 Be production rate calibration for the past 12 ka deduced from the radiocarbon-dated Grøtlandsura and Russenes rock avalanches at 69[END_REF][START_REF] Björck | The Preboreal oscillation around the Nordic Seas: Terrestrial and lacustrine responses[END_REF]. 10 Be ages are reported in years before 2010 CE.

R ESU L TS

The 10 Be chronology at Tsidjiore Nouve is shown in Figure 1 and Table DR1. Ages from the pre-LIA moraine crest (n = 10) range from 10.9 to 11.9 k.y. with an arithmetic mean and standard deviation of 11.4 ± 0.4 k.y. (±0.7 k.y. including the production rate uncertainty for comparison with other climate proxies, shown in italic below; Fig. 2). Ages from the recessional ridges of the pre-LIA moraine (n = 7) are between 10.9 and 11.5 k.y., and the arithmetic mean is 11.2 ± 0.2 k.y. (±0.6 k.y.; Fig. 2). The error-weighted mean ages, 11.44 ± 0.06 k.y. for the crest and 11.25 ± 0.08 k.y. for the recessional ridges are consistent with the stratigraphy.

ages in the central part of the outmost crest are 3790 ± 200 (ARO-12), 3320 ± 170 (ARO-13) and ridge, and the two younger boulders from the ice-proximal side (Fig. DR3), yielding an arithmetic mean of 3260 ± 180 years. (ii) We obtain ages of 480 ± 30 (ARO-10) and 180 ± 10 years (ARO-9) for boulders on an inner crest; ARO-9 is from a more ice-proximal area than boulder ARO-10 (Fig. DR3). On the lowest-elevation segment of the composite moraine, in a stratigraphically similar position as ARO-10, sample ARO-23 yields an age of 580 ± 30 years.

The five samples from the outer post-LIA moraine yield ages of 170 ± 10 years (ARO-66), 150 ± 10 years (ARO-67), 150 ± 10 (ARO-63), 110 ± 10 years (ARO-65) and 100 ± 10 years (ARO-64). One sample from the inner post-LIA moraine yields an age of 120 ± 10 years (ARO-62).

D ISC USSI O N

The glacier chronology at Tsidjiore Nouve yields constraints for five periods of Holocene glacier fluctuations (Fig. 3A), which we compare to other terrestrial records in the northern midlatitudes, and to North Atlantic temperature and tropical precipitation records (Fig. 3):

(i) Earliest Holocene: The mean age of the pre-LIA moraine crest of ~11.4 k.y. indicates that the Tsidjiore Nouve Glacier was close to its LIA-limits at or shortly after the YD/earliest Holocene transition. The mean age of the recessional ridges is slightly younger (~11.2 k.y.)

indicating that the glacier might have retreated slowly from the crest over a few centuries.

The age range for this moraine overlaps with both the late YD and the PBO (Fig. 3A), a brief cold-spell identified in the Greenland ice cores, peaking 11.40 k.y. ago and terminating 11.27 k.y. ago [START_REF] Rasmussen | Early Holocene climate oscillations recorded in three Greenland ice cores[END_REF][START_REF] Kobashi | 4±1.5°C abrupt warming 11,270 yr ago identified from trapped air in Greenland ice[END_REF] 3F). Further evidence for the PBO in European and North American lacustrine and glacial records (e.g., [START_REF] Björck | The Preboreal oscillation around the Nordic Seas: Terrestrial and lacustrine responses[END_REF][START_REF] Hu | Abrupt climatic events during the last glacial-interglacial transition in Alaska[END_REF] suggest a wide geographic footprint of this short cold spell.

(ii) Early and mid-Holocene: It is characteristic for the Alps that no moraines are preserved during the early and mid-Holocene, suggesting that glaciers were smaller during this period than during the LIA [START_REF] Holzhauser | Glacier and lake-level variations in westcentral Europe over the last 3500 years[END_REF]. In Greenland, a warming of 4 ± 1.5 °C within ~14 years is recorded at the end of the PBO [START_REF] Kobashi | 4±1.5°C abrupt warming 11,270 yr ago identified from trapped air in Greenland ice[END_REF]. The manifestation of this abrupt warming might have driven glaciers in the Alps far inside their LIA configuration, where they remained during the entire mid-Holocene. Warm early to mid-Holocene climate is also supported by recent studies from the European Alps: Mont Miné Glacier, a neighboring glacier of Tsidjiore Nouve, was smaller than today 9 k.y. ago, before a shortlived, minor advance 8.2 k.y. ago killed many trees growing upterminus [START_REF] Nicolussi | The 8.2 ka event -calendar dated glacier response in the Alps[END_REF]; Fig. 3C); Swiss glaciers were smaller than today Publisher: GSA Journal: GEOL: Geology Article ID: G33169

Page 7 of 15 during the early and mid-Holocene, based on sub-fossil wood and peat records (Joerin et al., 2006, Fig. 3D); a speleothem record implies that the Upper Grindelwald Glacier (Swiss Alps)

was in a quasi-continuously retracted position until 5.8 k.y., followed by a period of larger glaciers (Luetscher et al., 2011, Fig. 3E). Re-advances of glaciers in the Alps during the second half of the Holocene agree with the general cooling trend in Greenland [START_REF] Cuffey | Temperature, accumulation, and ice sheet elevation in central Greenland through the last deglacial transition[END_REF]; Fig. 3H) and correlate with a precipitation decrease in the tropics [START_REF] Haug | Southward migration of the Intertropical Convergence Zone through the Holocene[END_REF]; Fig. 3I), both beginning after the thermal maximum ~7.5 k.y. ago and culminating in the LIA.

moraine provide less robust age constraints, but indicate that the glacier was similar to LIA size 3.8 3.2 k.y. ago. We cannot entirely exclude the possibility of 10 Be accumulation prior to final boulder deposition. However, several arguments support the view that these three ages represent the formation of this moraine segment. First, the chronology agrees with the stratigraphic order. Second, the 14 C-dating of a larch log found in the basal till of the proglacial streambed yielded an age of 3440-2770 cal yr B.P. (Table DR2), consistent with the two younger 10 Be ages, implying an advance of the glacier at that time (Röthlisberger, 1976). This scenario is further supported by 14 C dates on two fossil soils that formed on this moraine 1510 1080 cal yr B.P. (lower horizon) and 1180 800 cal yr B.P. (upper horizon;

Röthlisberger, 1976, Table DR2). The 14 C dates assign a conservative minimum age to this moraine, supporting the pre-LIA 10 Be ages. Interestingly, both soil-dates fall within the Medieval Warm Period (MWP; Fig. 3B). Other studies in the European Alps support the onset of substantial glacier re-advances to at Upper Grindelwald Glacier implies glacier advances close to its Holocene peak extent from 3.8 k.y. ago onward [START_REF] Luetscher | Holocene glacier history from alpine speleothems[END_REF], Fig. 3E). [START_REF] Patzelt | Die postglazialen Gletscher-und Klimaschwankungen in der Venedigergruppe (Hohe Tauern, Ostalpen)[END_REF] 3 ago [START_REF] Wiles | E: Behavior of the Upper Grindelwald Glacier inferred from petrographic and stable isotope changes in speleothems (Luetscher et al[END_REF]. In the tropics, the period between ~3.8 and 2.8 k.y. ago is marked by centennial-scale, large-amplitude variations in precipitation, including several precipitation minima similar to those during the LIA [START_REF] Haug | Southward migration of the Intertropical Convergence Zone through the Holocene[END_REF]; Fig. 3I).

(iv) Little Ice Age: Lasting from about 1300-1850 CE, the LIA in the Swiss Alps showed three maxima, c. 1350, 1650 and 1850 CE [START_REF] Holzhauser | Glacier and lake-level variations in westcentral Europe over the last 3500 years[END_REF]. The three younger ages on the inner sub-ridge of the composite moraine correspond to the years 1430 ± 32, 1534 ± 28, 1829

± 11 CE, all within the LIA period. A drawing of Tsidjiore Nouve Glacier from 1836 CE (Fig. DR4) testifies the advanced position of the glacier terminus to be more than 1 km downvalley of its modern position, and is in general agreement with the boulder age corresponding to 1829 ± 11 CE (ARO-9).

(v) Little Ice Age Termination: Our six post-LIA boulder ages correspond to the years between 1840 ± 13 and 1906 ± 11 CE, broadly overlapping with the historic record (Glaciological reports, 1881 2009; Fig. DR1). Offsets of ~15 to 100 years between the 10 Be ages and the historically recorded deposition years of the post-LIA moraines (Fig. 1) are probably due to a small 10 Be signal (~300 -2000 atoms 10 Be g -1 ) inherited from previous exposure. This result underlines the potential to resolve (sub-)centennial climate events during the last millennium with 10 Be dating of moraines, previously only reported from New Zealand [START_REF] Schaefer | High-frequency Holocene glacier fluctuations in New Zealand differ from the northern signature[END_REF] and Peru [START_REF] Licciardi | Holocene glacier fluctuations in the Peruvian Andes indicate northern climate linkages[END_REF].

C O N C L USI O N
We make the case that Holocene glacier fluctuations in the western Alps (i.e., midlatitude temperatures) are related to North Atlantic temperature changes and the tropical hydrological cycle, confirming a hemispheric climate link. Our observations thus reinforce the thermal equator and the connected westerly winds have shifted south over several millennia during the Holocene in response to North Atlantic temperature changes, reaching its southernmost position during the LIA [START_REF] Sachs | Southward movement of the Pacific intertropical convergence zone[END_REF].

A C K N O W L E D G M E N TS
We thank R. Schwartz, J. Hanley and J. Frisch for help with sample preparation, A.

Putnam for assistance during figure preparation, and the staff of LLNL-CAMS for the excellent measurements. We acknowledge support by the CRONUS-Earth project (NSF Grant EAR-0345835), by the International Balzan Foundation and the German Academic Exchange Service (DAAD). This is Lamont-Doherty Earth Observatory publication 7540.

Comments from Jason Briner and two anonymous reviewers greatly improved this manuscript.

R E F E R E N C ES C I T E D Abbühl, L.,

Geology

Table DR1: Sample details, analytical data and surface exposure ages. Sample 10 Be/ 9 Be ratios, measured at the Center for Accelerator Mass Spectrometry of the Lawrence Livermore National Laboratory, were normalized to one of the indicated standards, KNSTD (until mid-2007; 10 Be/ 9 Be = 3.15 x 10 -12 ) or KNSTD07 (after the year 2007; 10 Be/ 9 Be = 2.85 x 10 -12 ). Measurements normalized to KNSTD were corrected by a factor of 0.9048 during age calculation. The ages are calculated using the 10 Be production rate with a value of 3.85±0.19 atoms (g yr) -1 (Balco et al. 2009, normalized to standard 07KNSTD) and the scaling method 'Lm' (time-dependent version of [START_REF] Lal | Cosmic ray labeling of erosion surfaces: in situ nuclide production rates and erosion models[END_REF] according to [START_REF] Balco | A complete and easily accessible means of calculating surface exposure ages or erosion rates from 10Be and 26Al measurements[END_REF]. They are reported in calendar years before 2010 CE. 1! analytical uncertainties range between 1.5% and 3%, except for samples younger than 500 years, which have uncertainties between 3% and 10%. Table DR2: Radiocarbon dates on fossil soils within superimposed deposits of the 'LIA' composite moraine (Figs. 1 and DR5a,b) and detrital logs found in basal till of the proglacial streambed at Tsidjiore Nouve Glacier (Fig. DR5c) and glacierclimatic significance (Röthlisberger, 1976, Röthlisberger andSchneebeli, 1979). All calibrated ages below are given as 2! intervals, referenced to the year 1950 CE (=BP, before present). They are calibrated with OxCal 4.1 [START_REF] Bronk Ramsey | Bayesian analysis of radiocarbon dates[END_REF][START_REF] Bronk Ramsey | OxCal Program 4.1[END_REF] relative to the IntCal09 calibration data set [START_REF] Reimer | IntCal09 and Marine09 radiocarbon age calibration curves, 0-50,000 years cal BP[END_REF] (Röthlisberger, 1976;Röthlisberger and Schneebeli, 1979). A: Photograph of the glacier and the moraine sequence in the year 1972 with white arrow pointing to the location of the excavated fossil soils. B: Schematic illustration of the two fossil soils found in superimposed moraine deposits in the 'LIA' composite moraine and corresponding 14 C ages. According to the authors, the younger fossil soil was covered by solifluction during a cold period. C: Radiocarbon dated larch log embedded in basal till and washed out by the glacial stream (3440-2770 cal years BP). D: Left lateral moraine sequence with pilot age estimates, which are based on 4 radiocarbon dates at Tsidjiore Nouve Glacier (summary in Table DR2) in combination with radiocarbon dates from the nearby glaciers Findelen, Mont Miné, and Ferpècle (Röthlisberger, 1976;Röthlisberger and Schneebeli, 1979). The deposition ages of 2500, 1500 and 900 years are presumably mostly inferred from 14 C ages at the three nearby glaciers. The white arrow in the back points to terminal moraine of the year 1817 CE.

!

  Fig.3G). Chronologically, two scenarios are possible: Either the pre-LIA moraine corresponds to a late YD position of the glacier, followed by a slow oscillatory retreat during the PBO marked by the recessional ridges; or the pre-LIA moraine was deposited during the PBO, interrupting the glacier retreat extensive late glacial position. Traditionally, the YD has been related to theIvy-Ochs et al., 2009, and references herein) with Egesen snowlines in the region of Tsidjiore Nouve Glacier (Valais) being ~200 m below the LIA maximum[START_REF] Maisch | Zur Gletschergeschichte des alpinen Spätglazials: Analyse und Interpretation von Schneegrenzdaten[END_REF], far downstream of the pre-LIA moraine. Adopting this scenario here, it appears more likely that the pre-LIA moraine is related to the PBO. Independent evidence in the Alps further supports the PBO signature (Fig.3):[START_REF] Burga | Gletscher-und Vegetationsgeschichte der südrätischen Alpen seit der Späteiszeit (Puschlav, Livigno, Bormiese[END_REF] reports an increase in non-arboreal pollen during the earliest Holocene in a peat bog at Palü Glacier 18 O changes in Swiss lake sediments showing a strong similarity with the Greenland 18

  .0 cal kyr ago based on 14 C-dates of the advance of two glaciers and climatic deterioration in a pollen record in Austria. Glacier culminations at that time are also documented in other regions of the Northern Hemisphere. According to Menounos et al. (2009) -14 C-dates from proglacial detrital and in situ wood. In southeastern Alaska, 14 C-dates of a tree-ring series indicate a glacial expansion ~3.4 3.0 k.y.

  Abbühl, L., Akçar, N., Strasky, S., Graf, A., Ivy-Ochs, S., and Schlüchter, C., 2009, A zero-

Figure 2 .

 2 Figure 2. Summed probability curves of 10 Be ages from the crest and the recessional ridges of the

Figure 3 .

 3 Figure 3. Comparison of Holocene climate records in the European Alps (upper panel), in

  Fig. DR1: Post-LIA length measurements of several glaciers in the Valais, Switzerland. The red and dark blue graph shows length changes of Arolla Glacier (direct neighboring glacier of Tsidjiore Nouve; Glaciological reports, 1881-2009). The dashed thick line of this graph indicates lack of continuous measurements between the years 1856 and 1886. Thin superposed lines are length measurements of seven other glaciers in the Valais starting around the year 1880, with the thick black line representing their average: 1 Fee Glacier, 2 Schwarzberg Glacier, 3 Lang Glacier, 4 Rossbode Glacier, 5 Trient Glacier, 6 Allalin Glacier, 7 Findelen Glacier. All glaciers show similar fluctuation patterns, i.e. during the general retreat from the LIA maximum a two-fold re-advance occurred between 1890 and 1930 CE, followed by another re-advance between 1970 CE and 1990 CE. The thick light blue curve is the summed probability plot of the individual 10 Be ages (thin light blue curves) from the outer post-LIA moraine.

Fig. DR2 :

 DR2 Fig. DR2: Examples of sampled boulders. A and B: ARO-1 (view west) and ARO-59 (view east), respectively, both embedded on the crest of the pre-LIA moraine. C: ARO-16, protruding from a recessional ridge of the pre-LIA moraine, view east.

Fig. DR4 :

 DR4 Fig. DR4: Upper panel: Drawing by Bühlmann from the year 1836 CE (Skizzenbücher Bd. 10, 256, Graphische Sammlung ETH, Zürich) showing the Tsidjiore Nouve Glacier on the right, the Arolla Glacier on the left and the village Arolla in the front. Glacier terminus positions are close to their Little Ice Age peak extent. Lower panel: Google Earth picture of approximately the same view in the year 2009 CE.

Fig. DR5 :

 DR5 Fig. DR5: Figures from the earlier work on the Holocene moraine sequence at Tsidjiore Nouve Glacier(Röthlisberger, 1976; Röthlisberger and Schneebeli, 1979). A: Photograph of the glacier and the moraine sequence in the year 1972 with white arrow pointing to the location of the excavated fossil soils. B: Schematic illustration of the two fossil soils found in superimposed moraine deposits in the 'LIA' composite moraine and corresponding 14 C ages. According to the authors, the younger fossil soil was covered by solifluction during a cold period. C: Radiocarbon dated larch log embedded in basal till and washed out by the glacial stream (3440-2770 cal years BP). D: Left lateral moraine sequence with pilot age estimates, which are based on 4 radiocarbon dates at Tsidjiore Nouve Glacier (summary in TableDR2) in combination with radiocarbon dates from the nearby glaciers Findelen, Mont Miné, and Ferpècle(Röthlisberger, 1976; Röthlisberger and Schneebeli, 1979). The deposition ages of 2500, 1500 and 900 years are presumably mostly inferred from 14 C ages at the three nearby glaciers. The white arrow in the back points to terminal moraine of the year 1817 CE.

Table DR3 :

 DR3 . Statistics of the 10 Be ages from the crest and the recessional ridges of the pre-LIA moraine including different mean ages with respective uncertainties and the reduced chi-squares. For the 'weighted mean', both mean age and uncertainty were weighted by the inverse variances (e.g.[START_REF] Taylor | An Introduction to Error Analysis[END_REF]. The relative high chi-square value is a result of the exceptionally small 1! analytical errors of the individual ages (~1.5%). Hence, for both crest and recessional ridges the boulder age distributions indicate that significant errors due to complex geological conditions (inheritance from prior exposure, erosion, or boulder instability) or snow cover are unlikely.

	Dated material	Uncal. age ( 14 C-yr)	Cal. age (yr BP) Glacier-climatic significance
	Upper fossil soil in moraine Lower fossil soil in moraine	1075 ± 80 1380 ± 85	1180-800 1510-1080	Timing of 'favorable' (warm) climate, minimum age for moraine deposition
	Fossil log	2940 ± 150	3440-2770	Timing of glacier advance,
	Fossil log	8400 ± 200	10120-8780	following a period of warmer climate
		10 B ages on crest (n=10)	10 Be ages on recessional ridges (n=7)
	Arithmetic mean ± standard deviation 11,440±380 years		11,230±210 years
	-incl. production rate uncertainty	11,440±680 years		11,230±590 years
	Weighted mean	11,440±60 years		11,250±80 years
	Peak age	11,470 years		11,190 years
	Median age ± interquartile range	11,430±730 years		11,230±290 years
	Reduced X 2	3.9		0.9


François-Xavier Coudert

Fabien Cailliez

Rodolphe Vuilleumier

Alain H Fuchs

Anne Boutin
email: anne.boutin@lcp.u-psud.fr

Water nanodroplets confined in zeolite pores

In a zeolite framework, the tetrahedrally bonded Si and Al atoms are referred to as the T atoms.

May 1, 2019

We provide a comprehensive depiction of the behaviour of a nanodroplet of 20 water molecules confined in the pores of a series of 3D-connected isostructural zeolites with varying acidity, by means of molecular simulations. Both Grand Canonical Monte Carlo simulations using classical interatomic forcefields and first-principles Car-Parrinello molecular dynamics were used in order to characterise the behaviour of confined water by computing a range of properties, from thermodynamic quantities to electronic properties such as dipole moment, including structural and dynamical information. From the thermodynamic point of view, we have identified the all-silica zeolite as hydrophobic, and the cationic zeolites as hydrophilic; the condensation transition in the first case was demonstrated to be of first order. Furthermore, in-depth analysis of the dynamical and electronic properties of water showed that water in the hydrophobic zeolite behaves as a nanodroplet trying to close its hydrogen bond network onto itself, with a few short-lived dangling OH groups, while water in hydrophilic zeolites "opens up" to form weak hydrogen bonds with the zeolite oxygen atoms. Finally, the dipole moment of confined water is studied and the contributions of water self-polarisation and the zeolite electric field are discussed.

Introduction

Water in zeolites represents a model system for a wide range of experimental as well as theoretical investigations, aimed at understanding the effect of confinement on the structure, [1] dynamics [2,3] and thermodynamics [4,5,6,7,8] of molecular fluids. Understanding the behaviour of nano-confined water is also relevant to such important problems as protein stability, [9,10,11,12,13] transport in ion channels, [START_REF] Beckstein | Proc. Natl. Acad. Sci[END_REF][START_REF] Hummer | [END_REF]16] nano-fluidic devices [17] and the thermodynamics of colloidal assemblies. [18,19] From a practical point of view, water plays a key role in many applications of zeolites such as ion-exchange and liquid phase separation.

Several computer simulation studies of water in zeolites have been reported in the past few years, using a variety of techniques including energy minimisation (with effective shell model-ionic potentials), [20,21,22,23,24] Car-Parrinello first-principles Molecular Dynamics (MD), [25,26,27,28,29] classical MD or Monte Carlo (MC) simulations using conventional two-body or polarisable forcefields. [7,8,30,31,32,33,34,35,36,37,38,39] Here we report, for the first time to our knowledge, a comprehensive description of the behaviour of 20 water molecules confined in the pores of a series of 3D-connected isostructural zeolites with varying acidity, by means of both classical Grand Canonical MC simulations and first-principles Car-Parrinello MD simulations. We characterise the behaviour of confined water by a range of properties, from thermodynamic quantities to electronic properties such as dipole moment, including structural and dynamical information.

2 Methods and systems

Systems studied

To understand the behaviour of water confined in the nanopores of a zeolite and the influence of the cationic content (or acidity) of the zeolite on the properties of water, we have studied three zeolites with the same structure and different Si:Al ratios. The chosen system was the LTA structure type, which displays supercages (also called α cages) of diameter 13 Å, connected to one another in a cubic symmetry by 8-ring windows of diameter 6 Å, and smaller sodalite cages of diameter 7 Å (see Figure 1). As opposed to previous simulations of 1D water chains confined in zeolitic materials such as Li-ABW, bikitaite, [25,26] and VPI-5, [27] the medium-size, 3D-connected pores of LTA make it a good candidate for studying the effect of confinement on liquid water condensation in open framework nanoporous materials. Furthermore, the small size of the LTA unit cell (approximately 12 Å) makes it a reasonable system for first-principles molecular dynamics simulations.

Synthesis of all-silica zeolites remains an experimental challenge. Although a synthesis of the reportedly all-silica LTA zeolite was published in the 1970's, [40] the proposed procedure (thermal dealumination) is known to result in materials exhibiting nonnegligible concentrations of defects [41] and no direct structural characterisation was reported in this study. To our knowledge, this synthesis was never reproduced. More recently, the all-silica LTA zeolite was synthesised and characterised by the groups of Corma [42] and Patarin. [43,44] Because the work reported in the present communication was started before these experimental results where published, the structure we used is an hypothetical structure taken from Ref. [START_REF] Baerlocher | Atlas of zeolite framework types[END_REF]. Its unit cell structure (Si 24 O 48 ) is cubic (P m 3m symmetry and a = 11.9 Å) and contains one supercage and one sodalite cage. It is very similar in both unit cell parameter and atomic positions to the experimental structures reported by Patarin and Corma. When it comes to cationic zeolites with the LTA structure, a large number of materials have been synthesised. Most natural and synthetic zeolites feature sodium as an extraframework cation. We have thus studied two different cationic LTA zeolites, containing 1 and 9 sodium cations respectively. The unit cell formulas are NaSi 23 O 48 Al and Na 9 Si 15 O 48 Al 9 , corresponding to Si:Al ratios of 23 and 1.67 respectively. To make our discussion easier to read, these two materials will hereafter be named Na 1 -LTA and Na 9 -LTA. Na 1 -LTA corresponds to a single aluminium atom and a single Na + cation per unit cell. As all T atoms * in the LTA structure are symmetrically equivalent, there is only one possible framework structure for this material, which is represented in Figure 2. As there has been no synthesis of a LTA zeolite with this particular Si:Al ratio, we have used the same structure and unit cell parameter as the all-silica LTA zeolite mentioned above.

The choice of Na 9 -LTA for the present study was dictated by our intention to model a zeolite with a high cationic content (low Si:Al ratio) and the constraints of the first-principles approach. This latter arises from the Loewenstein rule, [START_REF] Loewenstein | [END_REF] which states that no two aluminium atoms of a zeolite structure should be linked directly by an oxygen atom. This results in highly cationic zeolites having unit cells multiple times larger than the unit cell of the LTA structure, which are then too costly to perform first-principles molecular dynamics simulations. For example, the unit cell parameter of the well-known NaA zeolite (Na 12 Si 12 Al 12 O 48 , Si:Al = 1) is twice that of the pure-silica LTA structure. We thus performed a systematic exploration of all possible distributions of n aluminium atoms among the 24 T atoms of a single LTA unit cell that conform to the Loewenstein rule. These calculations proved that the only values of n compatible with both constraints are in the range 1 ≤ n ≤ 9. They furthermore showed that for n = 9, corresponding to the lowest Si:Al ratio, there are only two possible distributions. These two structures, though not symmetrically equivalent, are very similar (see Figure 1S). We have thus taken the first one as a starting point for our simulations, along with a unit cell parameter taken from the experimental structure of zeolite NaA (a = 11.26 Å). [47] 2.2 Monte Carlo simulations

Models and forcefields.

The development of models and forcefields to describe the interaction of zeolites with various guest molecules is a very active field of research. [48] In this study, the focus is put on the properties of confined water and their evolution with the nature of the zeolite. No chemical change occurs during water condensation in any of the LTA zeolites studied and Monte Carlo simulations were thus performed in the classical limit. We used existing forcefields, previously designed in the group of Fuchs and known to yield good transferability for both all-silica [7,49] and cationic zeolites. [50] For Na 1 -LTA and Na 9 -LTA systems, an average tetrahedrally bonded (T) atom was used, instead of distinguishing aluminium and silicium atoms. [50] It has been shown that this yields a reasonable description of sodium cation distributions in faujasites. [50,51] The zeolite framework was considered fixed. Its interaction with extraframework cations and guest molecules included a Coulombic term and a repulsion-dispersion term; the latter was modeled by a Buckingham (or exp-6) potential for the cations [52] and a Lennard-Jones potential [50] for water molecules. A summary of the forcefield parameters (partial charges, Lennard-Jones and Buckingham parameters) can be found in Tables 1 and2). The rigid, nonpolarisable TIP4P model [53] was used for water. Lorentz-Berthelot rules have been used to obtain cross potential parameters.

Grand Canonical Monte Carlo simulations.

Adsorption isotherms were obtained using Grand Canonical Monte Carlo simulations [START_REF] Nicholson | Computer Simulation and the Statistical Mechanics of Adsorption[END_REF][START_REF] Frenkel | Understanding molecular simulation: from algorithms to applications[END_REF] to compute the average number of adsorbed water molecules for several values of the chemical potential. The chemical potential of water was related to the thermodynamic pressure by using the ideal gas law in the vapour phase, and by computing the µ(P ) relation in liquid water as described in Ref. [7]. Preinsertion, [START_REF] Mackie | [END_REF] orientation [57,58] and jump [50] bias moves were used to accelerate the convergence of the Monte Carlo simulations.

All Monte Carlo runs were performed at 300 K on a 2×2×2 supercell to minimise the influence of periodic boundary conditions. Typical runs included 10 to 200 million steps to equilibrate the system followed by a production run of 30 to 100 million steps. Longrange electrostatic interactions were taken into account using the Ewald summation technique. In order to improve the efficiency of the calculations, the electrostatic and repulsion-dispersion interaction energies between the rigid zeolite framework and both water molecules and extraframework cations were precomputed on a grid (with a grid mesh of 0.2 Å) and stored for use during the simulation.

Grand Potential Computation.

The question of the nature of the transition in adsorption/desorption experiments is often subject to debate. Is this transition a real phase transition, between a gas-like and a liquid-like fluid? Likewise, when hysteresis occurs, its origin also raises several questions. Is it due to a bad sampling of the system or can it be related to the existence of multiple metastable states? In order to bring an answer to these issues, knowledge of the thermodynamic potential of the system is necessary. Such a quantity is not straightforward to obtain and requires the use of sophisticated techniques. Recently, Wang and Landau developed such a technique that is based on the computation of the density of states (DOS) on-the-fly, using biased Monte Carlo simulations. [59,60] (Once the DOS is known, the calculation of the thermodynamic potential is straightforward.) This method has been used successfully to study a wide range of systems, including Ising models, [59] Lennard-Jones fluids [61] and clusters, [62] or proteins. [62,63] It has also been revisited in the context of expanded ensembles. [64,65,66] In expanded ensembles, substates of the system are defined by the introduction of a new variable often called reaction coordinate ξ. The aim of those EXEDOS (Expanded Ensemble Density Of States) simulations is not to determine the density of states as a function of the energy, like in the original Wang-Landau algorithm, but as a function of the reaction coordinate, in order to compute variations of thermodynamic quantities along ξ. When dealing with adsorption, the thermodynamic potential is the Grand Potential Ω, since simulations are conducted in the Grand Canonical ensemble. We have applied the EXEDOS methodology to compute Ω as a function of the number of adsorbed water molecules in the zeolite framework. To this end, we have computed the density of states of the system as a function of N , noted Q(N ).

The probability to generate randomly a configuration with N water molecules is proportional to Q(N ). A flat distribution is thus obtained if each state is visited with probability proportional to 1/Q(N ). The acceptance probabilities of insertion and destruction moves have been modified accordingly (rotations, translations, and jump moves are still accepted according to the conventional Metropolis criteria):

p(N → N + 1) = min 1, V (N + 1)Λ 3 Q(N ) Q(N + 1) e -β∆E (1)
p(N → N -1) = min 1, N Λ 3 V Q(N ) Q(N -1) e -β∆E (2)
Λ is the de Broglie wavelength and ∆E is the energy change between the new and the old configuration. Q(N ) is determined iteratively during the simulation. Initially, Q is uniformly set to 1 for the entire range of values accepted for N . Q is then updated at each insertion or destruction trial move between two configurations containing respectively N 1 and N 2 water molecules. If the move is accepted, Q(N 2 ) is multiplied by a factor f (with f > 1) ; if the move is rejected, Q(N 1 ) is multiplied by f . At the same time, a N -histogram H(N ) is accumulated. The simulation is stopped when each bin has been visited at least 100/ √ ln f times, following Zhou and Batt. [67] At this time, the multiplicative factor is decreased, so that f → √ f . The histogram is reset to zero and a new cycle is started, until the histogram is flat again. This process is repeated until f reaches a small enough value. f is set to e 1 for the initial cycle, and 20 cycles are made, so that the last value of f is approximately equal to e 10 -6 . For small values of f (typically f < e 10 -4 ), the "minimum visit per bin" criterion may be too restrictive and the simulation is stopped if the total number of MC moves is superior to 200000/ √ ln f , as suggested by Poulain et al. [62] Once Q(N ) has been determined, one is able to calculate Ω for the system, using:

Ω(µ, V, T ) = -k B T ln N Q(N )e βµN (3)
For a given value of the chemical potential µ, the Landau free energy between two states containing respectively N 1 and N 2 molecules is:

Ω L (µ, N 2 ) -Ω L (µ, N 1 ) = -k B T ln Q(N 2 ) Q(N 1 ) e βµ(N2-N1) (4)
This allows to determine Ω L to within an additive constant. However, in our particular case, we can determine this constant, since we know that N = 0 only corresponds to one configuration for any value of µ, so that Ω L (µ, 0) = 0 kJ.mol -1 . Using the µ(P ) relationship for bulk water, we can also present variations of Ω (or Ω L ) as a function of the external pressure of the reservoir.

First-principles simulations

Molecular dynamics.

The molecular dynamics of anhydrous and hydrated zeolite were performed with the Car-Parrinello approach [68] as implemented in the CPMD software package, [69] using the Kohn-Sham formulation of DFT. We employed the gradient-corrected BLYP functional. [70,71] A plane-wave basis set with an energy cutoff of 70 Ry was used and only the valence electrons were explicitly considered, with semilocal norm-conserving Trouillier-Martins pseudopotentials taking into account the interactions between core and valence electrons. All molecular simulations were performed using periodic boundary conditions on a single cubic unit cell of the zeolite. Hydrogen atoms were deuterated in order to improve decoupling of the dynamics of the ionic and electronic subsystems and allow for a larger time step in the integration of the equations of motion. A timestep of 7 a.u. (0.17 fs) and a fictitious electron mass of 1000 a.u. were used. Molecular dynamics simulations of the anhydrous zeolite were started from the hypothetical structure for the purely siliceous LTA zeolite available in the literature. [START_REF] Baerlocher | Atlas of zeolite framework types[END_REF] For the hydrated zeolites, starting configurations were extracted from the previously performed Grand Canonical Monte Carlo simulations. A first equilibration run of at least 1.2 ps was performed on each system, during which the temperature was controlled to 300 K by a Nosé-Hoover thermostat [72,73] with a frequency of 3800 cm -1 . Production runs of 5 to 6 ps in the microcanonical NVE ensemble were then performed, from which were extracted all the quantities reported in this article.

Analysis of electronic properties.

The total polarisability M (t) of the periodic unit cell was defined using the Berry phase approach of Resta. [74,75] The IR absorption coefficient of the system, α(ω), was then computed by calculating the autocorrelation function of the polarisability, within the framework of linear response theory:

α(ω) ∝ ω tanh (βhω/2) n(ω)V ∞ -∞ dt e -iωt ij µ i (t)µ j (0) (5)
where n(ω) is the refractive index. A harmonic quantum-classical correction was then applied to the classical correlation function to account for the quantum aspects of the hydrogen motion, as was shown to be adequate in previous works. [76,77,78] The resulting infrared spectrum was then smoothed by a linear filter as in Ref. [77,78]. Dipole moments of individual water molecules were also computed: in order to do so, the electronic wavefunctions are replaced by the maximally localised Wannier functions [79] and the Wannier centers are assigned to each water molecule. [80,81,77] The calculation of both the macroscopic polarisability and the Wannier functions was performed every 1.7 fs of the dynamics, for each system.

To characterise the presence of hydrogen bonds and their strength from an electronic point of view, we propose to study atomic polar tensors on hydrogen atoms. These tensors and particularly their strong anisotropy can be shown to be a signature of hydrogen bonding in liquid water. [START_REF] Vuilleumier | Liquid water dynamical partial charges and infrared intensities[END_REF] The atomic polar tensor of a given atom n is defined as α

(n) ij = ∂M i ∂R n,i
, where i and j can be x, y or z

where M i is the component of the total polarisation of the system along axis i and R n is the position of atom n.

In this study, the atomic polar tensors were calculated in the Density Functional Linear Response Theory, [START_REF] Gonze | [END_REF] using Putrino et al. [84] implementation in the CPMD code [69] of the linear response to an external electric field. The atomic polar tensors of all atoms in a given configuration of the system were computed by three calculations of the linear response of the system to an external electric field, according to the following formula:

α (n) ij = ∂F (n) j ∂E i (7)
where F (n) is the force acting on atom n and E is the external electric field. For each system studied, this analysis was performed on 200 configurations extracted at regular interval from the NVE dynamics.

From the APT tensor we derived an atomic charge, namely the Born charge, as the trace of the APT tensor. To better characterise the APT which is a nonsymmetric 3 × 3 tensor we also determined the singular values of the APT from singular value decomposition. The three singular values (ξ 1 , ξ 2 and ξ 3 ) have the dimension of a charge and their differences characterise the anisotropy of the atomic polar tensor. For Na 9 -LTA, the isotherm is a typical type I isotherm, according to the IUPAC classification scheme, [START_REF] Rouquerol | Adsorption by Powders and Porous Solids[END_REF] commonly found in the case of hydrophilic aluminosilicate zeolites. It is characterised by a strong water uptake at low pressure (high value of Henry's constant). In contrast, the isotherm for the all-silica zeolite is typically that of water intrusion into an hydrophobic nanoporous solid: [49,[START_REF] Eroshenko | [END_REF] a negligible water uptake is observed up to a pressure of 100 MPa (i.e. much higher than the saturated vapour pressure of 3500 Pa), at which an abrupt condensation transition takes place. Moreover, in that case, a hysteresis is observed: water extrusion takes place at 45 MPa.

Results and discussion

In the case of Na 1 -LTA, three different stages can be distinguished during water adsorption. At low pressure, 2 to 5 water molecules per unit cell are involved in small clusters around sodium cations, equally distributed in each of the eight unit cells of the simulation box (see Figure 4a). Around 250 Pa, a first transition occurs that corresponds to the filling of half the unit cells of the simulation box (Figure 4b). This is accompanied by a redistribution of the cations, so that 4 unit cells are filled with 2 cations and water, whereas the other 4 are almost empty. This is the stable situation until 2500-3000 Pa where a second transition occurs leading to the total filling of all the supercages. This transition is accompanied once again by the displacement of sodium cations so that, in the end, the system is homogeneous with each unit cell filled with water and one cation (Figure 4c). It is interesting to point out the importance of simulating eight unit cells instead of one, because it would otherwise not have been possible to observe such cationic redistribution and inhomogeneous pore filling.

The change in shape of the adsorption isotherms from all-silica to Na 9 -LTA can be regarded as an evolution from a purely hydrophobic nanoporous solid to a hydrophilic one, the intermediate Na 1 -LTA case being that of a hydrophobic medium containing strong hydrophilic defects. If we ignore the two-step filling mechanism (which is a consequence of the mobility of sodium cations), the adsorption isotherm in Na 1 -LTA is of type IV: adsorption at low pressures followed by a plateau, and a condensation transition occurring below the saturated vapour pressure P 0 . [START_REF] Rouquerol | Adsorption by Powders and Porous Solids[END_REF] To our knowledge, this is the first time that such a "pure" type-IV isotherm is predicted in a nanoporous solid. In previous works, [39,87] we observed that the inclusion of a silanol-type defect could result in type-IV isotherm, but the condensation transition occurred above P 0 .

We now discuss the effect of the cations on water density. Values for the three systems are shown in Table 3. Porous volumes have been computed using an approach similar to that of the calculation of Connoly surfaces, and is described in detail elsewhere. [7,[START_REF] Cailliez | Proceedings of 4th International FEZA Conference[END_REF] In the case of all-silica and Na 1 -LTA systems, a negligible amount of water molecules is found in sodalite cages. Thus, the porous volume only takes into account supercages. In the case of Na 1 -LTA and Na 9 -LTA, we also compute a porous volume and a density that takes into account the presence of sodium cations (considered as a sphere with a radius being half the σ value of sodium, i.e. 1.29 Å). Water in all-silica LTA is strongly depleted with respect to bulk water, as observed earlier in other siliceous zeolites. [7,[START_REF] Cailliez | Proceedings of 4th International FEZA Conference[END_REF] The presence of only one sodium cation also lowers water density. This is because the presence of one cation is not sufficient to modify strongly water-water interactions (see Section 3.2.1), and because of the volume occupied by the sodium cation. On the other hand, the presence of the nine cations results in a large increase of water density, up to 0.99 g.cm -3 , which is similar to bulk water density. Cations modify water structure in the LTA framework, so that a larger amount of water can be introduced. However, as will be shown in Section 3.2, water structure in Na 9 -LTA is different from the bulk water.

On the nature of the condensation transition.

In order to determine the nature of the condensation transition in all-silica LTA, we have computed the evolution of the grand potential Ω with pressure. Figure 5a shows the profile of Ω L as a function of the water content for various pressures. Below roughly 20 MPa, only one stable state exists, that corresponds to an empty LTA zeolite. When pressure is increased, a filled metastable state appears, that is however less stable than the empty state. Around 70 MPa, the filled configuration becomes the most stable state. The existence of two stable states is characteristic of a firstorder transition between gas-like and liquid-like phases. The conversion between those two states is possible on a wide range of pressure (20 MPa → 190 MPa), as seen in Figure 5, but the system has to overcome a macroscopic energy barrier. Monte Carlo simulations proceed via microscopic moves. In a finite time MC run, the most stable state may not be reached, and this is why hysteresis can be observed. In the present case, although the equilibrium pressure is approximately 70 MPa, intrusion and extrusion occur respectively at 85 MPa and 40 MPa in the simulations. This is to our knowledge the first time that the Wang-Landau algorithm is used in the context of fluid intrusion in porous solids, and that the existence of a firstorder transition is really demonstrated. We believe that this represents an important progress for the study of adsorption (intrusion) processes, since we are now able to derive thermodynamic quantities that were hardly accessible until now.

On the hydrogen bonding of water and the zeolite

Structure and orientation dynamics of water.

To study the interaction of water with the zeolite and add a microscopic insight to the thermodynamic vision of water developed in the previous section, we have first studied the structure and orientation of water confined in the all-silica, Na 1 and Na 9 LTA zeolites. The arrangement of water molecules in the α cages of the zeolites is irregular and reveals none of the features exhibited by gas phase water clusters. [START_REF] Cambridge Cluster Database | [END_REF] In order to probe the nature of the hydrogen bond network of confined water, we have studied the orientation of the molecules around neighbouring water molecules and zeolite oxygen atoms. The key quantities for that task are the OO distances and HOO angles (respectively noted d and α in the following, see Figure 6), where O and H are atoms from a same water molecule, and O is a neighbouring oxygen atom. Following the common practice, we used a geometric criterion [START_REF] Luzar | [END_REF] to characterise hydrogen bonds and study quantitatively the evolution of water structure between the different zeolites: we consider the hydroxyl group OH of a water molecule hydrogen bonded to an acceptor atom O if d ≤3.5 Å and α ≤ 30°. The limits of the criterion are based on the position of the minima of the 2D distribution function of liquid water orientation. This criterion is widely used in the literature, as well as a few stricter variants imposing an additional restriction on the OH distance. Using this criterion, we have studied the evolution of the number of both water-water and water-zeolite hydrogen bonds when the Si:Al ratio of the zeolite is lowered (see Table 4a). First, it can can be seen that the number of water-zeolite hydrogen bonds in the all-silica zeolite is fairly small and that the number of water-water hydrogen bonds is still large, while smaller than for same-density bulk water (around 1.6 per water molecule). The confined water in the hydrophobic, all-silica zeolite can thus be seen as a nanodroplet of liquid water with its hydrogen bonds network largely closed onto itself, and a few OH groups oriented towards the zeolite. When we move to cationic zeolites, the number of water-zeolite bonds increases while water-water hydrogen bonds diminish; the net result in the total number of hydrogen bonds per water molecule is a slight increase with acidity. As the acidity of the zeolite framework and the number of cations grow, the droplet of confined water "opens up" and forms more hydrogen bonds with the zeolite, at the expense of the water-water hydrogen bonds framework.

To further understand the interaction between the zeolite and the confined water, the classical Monte Carlo simulations were complemented with first-principles molecular dynamics simulations, in order to provide valuable information on the intermolecular and vibrational dynamics of the system as well as to probe its electronic properties (the details of these Car-Parrinello DFT simulations can be found in Section 2.3). Table 4b shows the number of water-water and water-zeolite hydrogen bonds counted by the same geometric criterion, for all three zeolites studied. Comparison to Table 4a shows a similar evolutions for both classical models and DFT-based simulations. The absolute number of hydrogen bonds counted is, however, larger in the latter case. This is in keeping with observations that Kohn-Sham DFT liquid water described with the BLYP exchange-correlation function is overstructured at ambient conditions. But the use of molecular dynamics simulations also allows us to get a picture of the dynamics of these hydrogen bonds, in the form of the autocorrelation function of the existence of hydrogen bonds. Figure 7 shows this autocorrelation function for water-zeolite hydrogen bonds in all three zeolites. It can be seen that the OHO alignments, counted as hydrogen bonds by the geometric criterion, decay very rapidly in the case of all-silica zeolite and Na 1 -LTA, while the alignments in Na 9 -LTA subside longer. The hydrogen bond lifetime for each zeolite was obtained by an exponential fit of the 0.05 to 1.5 ps region of the autocorrelation function, where they display monoexponential behaviour. The lifetimes of water-zeolite hydrogen bonds, reported on Figure 7, range from 125 fs to 520 fs and are smaller than the lifetime of water-water hydrogen bonds both in bulk water (1.5 ps, following the same procedure described above) and in confined wa-ter (between 1.3 ps and 1.6 ps depending on the number of cations). The very short lifetime of the "geometric H bonds" in both the all-silica zeolite and Na 1 -LTA is incompatible with the chemist's view of a hydrogen bond. It is thus plausible that what the geometric criterion counts as hydrogen bonds in these cases are occasional, short-lived OH• • •O zeolite alignments, corresponding to dangling OH groups of the confined water. On the other hand, the lifetime of 520 fs observed in Na 9 -LTA is compatible with weak water-zeolite hydrogen bonds. To characterise these two different situations further, we gain insight from the first-principles dynamics and turn to a electronic characterisation of hydrogen bonds.

Electronic characterisation of hydrogen bonds: the atomic polar tensor.

By using first-principles molecular dynamics, we can clearly depict the electronic properties of the confined water and thus help understanding its interaction with the confining zeolites. This is needed to characterise hydrogen bonds in a way that will go past the limits of simple geometric criteria. To probe the electronic nature of the hydrogen bonds, we studied the anisotropy of the atomic polar tensor of the hydrogen atoms. In the presence of a hydrogen bond, the hydrogen atom involved exhibits a fairly anisotropic atomic polar tensor. This anisotropy can be characterised by the three singular values of the atomic polar tensor (see Section 2.3). For example, in the case of the water dimer shown in Figure 8, it can be seen that the ellipsoid defined from the singular value decomposition of the APT tensor of the hydrogen atom involved in the hydrogen bonding exhibits a clearly oblong shape, with one singular value much larger than the others and the corresponding singular vector directed along the axis of the hydrogen bond. Physically, the ellipsoid shown in Figure 8 shows the magnitude of dipole moment change upon a small displacement of the hydrogen atoms in the direction of the principal axis. For an hydrogen involved in an H bond, a displacement parallel to the H bond can be seen to lead to a much larger change of dipole moment than for a displacement perpendicular to it. This arises from noticeable charge transfer across the H bond when the hydrogen is moved along the H bond.

We have thus applied this method to study the strength of both water-water and water-zeolite bonds identified by the geometric criterion; Figure 9 shows the average eigenvalues in each case, as well as for bulk water. It can be seen that the strength for water-water hydrogen bonds in the confined water is fairly similar to that of bulk water: the largest singular value, ξ 1 , is close to 1 e. Singular values ξ 2 and ξ 3 , however, are smaller than in the bulk.

The water-zeolite hydrogen bonds, on the other hand, can be seen to be significantly weaker than the water-water bonds. Moreover, we see a gradual increase of ξ 1 when we go from the all-silica LTA ( ξ 1 0.36) to Na 1 -LTA ( ξ 1 0.40) and then Na 9 -LTA ( ξ 1 0, 77). Both findings are consistent with the lifetimes reported above, but quantify the nature and strength of these hydrogen bonds in a different perspective: the difference in lifetime, for example, between water-zeolite bonds in Na 9 -LTA and the bulk water-water hydrogen bond lifetime is markedly more important than the difference in the anisotropy of the atomic polar tensor. The picture they draw for the evolution of water-zeolite interaction with increasing cationic content of the zeolite, however, is the same: while the water confined in the all-silica zeolite can be considered as a droplet forming no hydrogen bonds with the zeolite, but exhibiting a few dangling OH groups, water confined in the same-structure zeolite with a Si:Al ratio of 1.67 forms relatively weak hydrogen bonds with the oxygen atom of the zeolite framework. This difference in the water-zeolite interaction strengthens the explanation of the thermodynamic properties of water condensation in both all-silica and cationic zeolites, as presented in Section 3.1: the hydrophobicity of the internal surface of allsilica zeolites, due to the extreme weakness of the water-zeolite attraction compared to the water-water dispersion energy, leads to intrusion in the liquid phase, while the larger water-zeolite interaction of acidic zeolites leads to adsorption in the gas phase.

The Na 1 -LTA, on account of the thermodynamic, dynamic and electronic results shown here, can be seen as a situation midway between the two cases just discussed. Because of the strong electric field created by the presence of an aluminium atom and a cation in the zeolite, water molecules are strongly attracted in the vicinity of the cation, leading to a first water uptake in the gas phase (less than 500 Pa). However, in spite of this local attraction, the internal surface is still only weakly attractive for water, as confirmed by the study of the nature of the water-zeolite interaction in Na 1 -LTA. This leads to the condensation of water at a pressure close to the saturation vapour pressure.

Infrared spectrum and vibrational spectroscopy of the all-silica zeolite.

To link our simulations to experimentally observable properties, we have also studied the vibrational properties of the all-silica zeolite and the water confined inside its pores. No experimental data is available today on this system, but vibrational spectroscopies have been widely used to investigate the behaviour of water at interfaces. [91,92] We have calculated the infrared spectrum of both the anhydrous zeolite and the hydrated system. Through the use of the Berry phase approach and maximally localised Wannier functions described in Section 2.3, it was also possible in that latter case to separate the contributions of water and zeolite to the total spectrum. We have first shown that the calculated infrared spectrum of the anhydrous all-silica LTA (Figure 10, lower panel) is in reasonable agreement with the experimental data [93] (Figure 10, upper panel) of the material synthesised by Fyfe et al. [40] The general shape of the spectrum is reproduced, the positions of the absorption bands present relatively large shifts (up to 200 cm -1 ). These are attributed to the fact that the simulation was performed at constant volume, with a unit cell parameter that was not optimised, leading to constraints on the vibrations modes of the zeolite, as well as the dealumination procedure used in Ref. [40], which is known to yield materials with nonnegligible concentrations of defects. [41] For the hydrated zeolite, we have shown that the presence of water has no significant effect on the infrared spectrum of the zeolite. Moreover, the vibration modes of the zeolite and the confined water are not correlated and the total spectrum of the system matches the sum of the two separate spectra for water and the zeolite, at frequencies higher than 300 cm -1 . [29] We have then focused our study on the OH stretch region of the confined water IR spectrum, shown in Figure 11. It can be seen that, while the ν OH peak of bulk water is large and structureless, that of water confined in allsilica LTA exhibits a narrower high-frequency peak on top of a larger band. This peak is attributed to dangling OH groups, which are not part of the water hydrogen bonds network and are oriented towards the zeolite. Such a peak has been reported for water near a variety of hydrophobic surfaces [94,95] and can thus be used as a signature of hydrophobicity. As an example of systems where it has been observed, the lower panel of Figure 11 shows the IR SFG (Sum Frequency Generation) spectrum of water near an H 2 O/CCl 4 interface, from Scatena et al. [95] In the present case of water confined in an all-silica zeolite, we have further confirmed the attribution of this high-frequency ν OH peak to dangling OH groups by studying the correlation between the average stretching frequency of water OH groups and the percentage of simulation time during which they are involved in water-water hydrogen bonds (according to the geometric criterion). Figure 12 shows that there is a clear correlation between these two quantities: OH groups with ν OH ≥ 2250 cm -1 , which compose the high-frequency peak, are markedly less involved in the water-water hydrogen bonds network than the other groups.

Dipole of the confined water

The dipole moment of the water molecules is an important quantity in order to understand the behaviour of confined water and to help assess the quality of the simulation models used. The histogram of the dipole moments of water molecules confined in the α cage of the three zeolites under study are reported in Figure 13, as well as the distribution obtained for a bulk water simulation including 32 water molecules in a periodic unit cell. It can be seen that the distributions for confined water peak around 2.9 D, which is close to the peak value of bulk water at 3.1 D. The self-polarisation of confined water is thus already quite large in these relatively small nanopores; from the point of view of the dipole moment, the confined water is liquid-like. Moreover, Figure 13 and Table 5 reveal that the dependence of the dipole moment on the acidity of the zeolite is very small. The histograms for the all-silica, Na 1 -LTA and Na 9 -LTA zeolites have the same peak value, and the only statistically significant difference is that zeolites with more cations present a less symmetric distribution, with a more important high-dipole wing, attributed to water molecules in close contact with the cations.

The relatively small influence of the cationic content of the zeolite is somewhat counterintuitive, as the presence of ionic species in close proximity to the water molecules might be expected to polarise the water molecules significantly. Indeed, calculations performed on a system comprised of a single water molecule close the Na + cation in Na 1 -LTA have yielded a water dipole moment of 2.4 D, which represents a large increase (0.6 D) from the gas phase value. These facts lead to the logical conclusion that the dipole moment of the confined water is almost entirely due to the self-polarisation of water, while the influence of both the zeolite and the extraframework cations is faint. Figure 14 presents the histogram of the dipole moments of 15 and 20 water molecules confined in all-silica LTA, which have their maximum at 2.7 D and 2.9 D respectively. This shows that the dependence of the dipole on the water content of the zeolite is more important that the dependence on the number of cations, and that the self-polarisation of water is the main factor determining the dipole in the systems we study here. This conclusion was further strengthened by the analysis of the quadrupole moments of the confined water molecules, presented in Table 5. The quadrupoles are only weakly affected by the changes in the cationic content of the zeolite while they have a more pronounced dependence on the number of water molecules.

Partial charges from population analysis and their use as forcefield parameters

There are many different ways to calculate charges from electronic structure. We have here investigated the atomic charges of zeolite atoms from four different definitions.

One of our goal was to compare the obtained values from these definitions with charges used in the empirical potentials that were shown to reproduce correct thermodynamics. Two very common methods are the Mulliken population analysis and Hirshfeld population analysis. However, of possible definitions of atomics charges, ESP charges fitted on the electrostatic potential seem the best suited for extracting charges for building an empirical potential. On top of these common atomic charge definitions we have also computed the Born charges obtained as the trace of the APT tensor discussed above. Table 6 shows the resulting charges for atoms of the all-silica zeolite. All these charges give a large spectrum of values. None seem however to agree with the range of empirical charges (-0.6 e to -0.8 e for O) that can correctly describe the hydrophobicity of silicalite-1. [7] We have also checked the dependance of Mulliken charges on the electronic method used. To this end we selected a few snapshots from the first principle molecular dynamics of empty all-silica zeolite and computed Mulliken charges with a gaussian basis set, using the CP2K program [96] † instead of a plane wave basis set. The resulting charges are -0.6 e for oxygen atoms and 1.2 e for silicon atoms, yet again different from other approaches.

Nevertheless, if Mulliken charges are computed for the hydrated all silica zeolite as well, no variation is observed with respect to the empty zeolite. This is consistent with observations made previously in this paper (from the vibrational dynamics analysis for instance) that the coupling between water and the all-silica zeolite is very small. This is also in line with the charges on zeolite atoms independent of filling commonly assumed in empirical forcefields.

Turning now to water, the distribution of Mulliken and Hirshfeld charges are shown in Figure 15. For these two definitions, a decrease of atomic charges is found from an isolated water molecule to bulk liquid water or confined water. Both Mulliken and Hirshfeld population analysis thus give rise to atomic charges that is in contradiction with the increase of water dipole moment. This may be understood by noting that H bond leads to an increase of electron density in the intermolecular region between O and H. Instead of an increase of molecular polarisation, the Mulliken and Hirsfeld population analysis attribute some of that density to the partner hydrogen of a donating oxygen, leading to a decrease of absolute values of charges arising from this apparent intermolecular charge transfer.

On the contrary, the Born charges (see Figure 16) show an increase of atomic charges upon going from an isolated molecule to a confined water molecule and to a bulk water molecule. This increase of Born atomic charges is furthermore in near quantitative agreement with the dipole moment variation: assigning these charges to hydrogen and oxygen atoms, we find dipole moments of 1.4 D, 2.77 D and 3.2 D for an isolated molecule, a confined molecule and a bulk water molecule, respectively. These charges also appear to be in the range of partial charges found in empirical forcefields for liquid water.

From considerations of infrared intensities, [START_REF] Vuilleumier | Liquid water dynamical partial charges and infrared intensities[END_REF] it is possible to define an alternative to the Born charge from the APT tensor by

ξ 2 = ξ 2 1 + ξ 2 2 + ξ 2 3 3 , (8)
where ξ 1 , ξ 2 and ξ 3 are the singular values of the APT tensor as discussed above and ξ has thus the dimension of a charge. These charges are found to be slightly lower that the Born charges (ξ H = 0.39 e for confined water and 0.54 e for bulk liquid) but show the same trend, with an increase of molecular polarisation from confined to bulk water. APT derived charges and in particular Born charges appear as interesting candidates to define meaningful atomic partial charges, at least for molecular systems. As mentioned before the Born charges of the zeolite atoms appear however too large with respect to the range of charges used in classical forcefields.

Concluding remarks

The use of both classical GCMC and first-principles CPMD simulations has allowed to describe in some details the behaviour of a nanodroplet of water confined in the supercage of LTA zeolites. Each method has its (known) strengths and drawbacks, and provides information in a physical chemistry domain where the other does not. Whenever the same kind of information can be obtained by both techniques, the agreement was found satisfactory. This was the case for instance for the number of hydrogen bonds between water molecules and the framework atoms as well as between water molecules themselves. Work is in progress to get a better insight into the hydrogen bond reorientation dynamics of water in such confined systems. [START_REF] Laage | [END_REF] The all-silica version of LTA as well as two model systems with 1 and 9 sodium cations were chosen for this study. From the thermodynamic point of view, we have identified the all-silica LTA zeolite as hydrophobic, and the cationic Na 9 -LTA zeolites as hydrophilic. The condensation transition in the first case was demonstrated to be of first order. The adsorption isotherm in the case Na 1 -LTA and Na 9 -LTA was of Type IV and I respectively. Furthermore, the in-depth analysis of the dynamic and electronic properties of water showed that water in the hydrophobic zeolite behaves as a nanodroplet trying to close its hydrogen bonds onto itself, with a few shortlived dangling OH groups, while water in the hydrophilic zeolites opens up to form weak hydrogen bonds with the oxygen atoms of the framework. We conclude that a thoughtful combination of methods, in which each technique is used for what it was designed for, can provide a nice insight into the physics and chemistry of such a complex interfacial system as water confined in zeolite nanopores.

Molecule

Charges (e) TIP4P water q O = -1.04 q H = +0.52 all-silica LTA framework q O = -0.7 q Si = +1.4 Na 1 -LTA framework q O = -0.7125 q T = +1.3833 Na 9 -LTA framework q O = -0.8125 q T = +1.25 Na + q Na = +1 Figure 16: Distribution of partial charges of the water oxygen and hydrogen atoms according to the Born population analysis methods. The distributions were computed for water confined in an all-silica LTA (red), for bulk water (black) and for a single gas phase water molecule (green).

3. 1

 1 Thermodynamics of water in LTA zeolites 3.1.1 Water adsorption isotherms. Water adsorption isotherms for the three different zeolite models are shown in Figure 3.

Figure 1 :

 1 Figure 1: Schematic view of the LTA structure showing an α cage (center) with its eight surrounding sodaliate cages. α cages are linked to one another by octogonal windows.

Figure 2 :Figure 3 :Figure 4 :Figure 5 :Figure 6 :

 23456 Figure 2: Left: Structure of the Na 1 -LTA zeolite, containing one aluminium atom (in green) and one extraframework Na + cation (in blue). The cation is represented in its most stable cristallographic site, known as Na(3).[47] Right: Structure of Na 9 -LTA zeolite with nine aluminium atoms and 9 Na + cations. Eight cations are in sites Na(1), and one is in a site Na(3).

Figure 7 :

 7 Figure7: Autocorrelation function C(t) of the existence of water-zeolite hydrogen bonds, as defined per the geometric criterion described in the text, for water confined in different zeolites of identical framework type: all-silica (red), one cation (blue) and nine cations (red). The lifetime of the bond, indicated on each curve, is obtained by an exponential fit of the 0.05-1.5 ps region of the autocorrelation function.

Figure 8 :

 8 Figure 8: Graphical representation of the Atomic Polar Tensor (APT) of atoms in a water dimer. The ellipsoids are drawn from the singular value decompositions of the APT tensors. They describe the magnitude of dipole moment changes upon displacement of the atoms in the principal directions of the ellipsoids. The lengths of the principal axis of these ellipsoids are thus the singular values of the APT tensors. It can be noted that the hydrogen involved in the hydrogen bond exhibits a significant anisotropy of the APT tensor.

Figure 9 :Figure 10 :

 910 Figure 9: Average singular values of the symmetric of the atomic polar tensor (see Section 2.3) for hydrogen atoms of water confined in the three zeolites (all-silica, one cation and nine cations): ξ 1 in black, ξ 2 in red and ξ 3 is blue. The unit of the singular values is the atomic unit of charge, e. The upper panel features the singular values for hydrogen atoms involved in water-water hydrogen bonds; the values for bulk water are indicated for comparison. The lower panel features singular values for hydrogen atoms oriented towards the zeolite.

Figure 13 :

 13 Figure13: Distribution of the dipole moment of water molecules in confined in the all-silica LTA (red), Na 1 -LTA (blue) and Na 9 -LTA (green). The distribution of dipole moments in bulk water is drawn in black.

Figure 14 :

 14 Figure14: Distribution of the dipole moment of water molecules in bulk water (black, full line), of 20 water molecules confined in an all-silica LTA zeolite (red, dashed line) and of 15 water molecules inside the all-silica zeolite (red, dotted line). The dipole moment of a water molecule in gas phase is also presented (green line).

Figure 15 : 29 - 1

 15291 Figure15: Distribution of partial charges of the water oxygen and hydrogen atoms according to the Mulliken and Hirshfeld population analysis methods. The distributions were computed for water confined in an all-silica LTA (red), for bulk water (black) and for a single gas phase water molecule (green).

Table 1 :

 1 Forcefield parameters used for modelling repulsion and dispersion interactions in Monte Carlo simulations.

	Interaction	Parameter Value
	Na-Na	σ	2.584 Å
			50.27 K
	O water -O water	σ	3.154 Å
			78.03 K
	O zeolite -O zeolite σ	3.000 Å
			93.53 K
	Na-O zeolite	α	61.1 × 10 6 K
		β	4.05 Å-1
		γ	76.52 × 10 4 K. Å6

Table 2 :

 2 Partial charges borne by atoms in the Monte Carlo simulations.

Table 3 :

 3 Filling of LTA at saturation, porous volume, and computed density of water. The porous volumes and densities shown in brackets take into account the excluded volume due to sodium cations.

	System	N/u.c. V P /u.c. ( Å3 ) ρ conf (g.cm -3 )
	LTA all-silica 21	754	0.83
	Na 1 -LTA	20	754 (745)	0.79 (0.80)
	Na 9 -LTA	30	983 (901)	0.91 (0.99)

Table 4 :

 4 Number of hydrogen bonds among the water molecules confined inside the α cage of LTA zeolites with different numbers of cations, using a geometric criterion (see text), from Monte Carlo simulations (upper table) and first-principles molecular dynamics (lower table). For each zeolite, the number of water-water hydrogen bonds and water-zeolite hydrogen bonds is indicated at well as the ratios per water molecule.

	(a) From Monte Carlo simulations
	zeolite	water-water H bonds	water-zeolite H bonds	all H bonds
		total per H 2 O	total per H 2 O	per H 2 O
	all-silica	25.6 1.35	3.9	0.20	1.55
	one cation	23.2 1.22	6.0	0.32	1.54
	nine cations 17.6 0.93	13.0 0.68	1.61
	(b) From first-principles molecular dynamics
	zeolite	water-water H bonds	water-zeolite H bonds	all H bonds
		total per H 2 O	total per H 2 O	per H 2 O
	all-silica	30.8 1.54	4.1	0.21	1.75
	one cation	27.2 1.43	7.1	0.37	1.80
	nine cations 22.9 1.21	14.0 0.74	1.95

Table 5 :

 5 Dipole moment and quadrupole eigenvalues of water confined in different LTA zeolites. The values for bulk water, ice and an isolated water molecule are also shown for comparison purposes.

	system dipole (D) quadrupole eigenvalues (esu.cm 2 )
	15 H 2 O in all-silica LTA 2.79	-3.04 -0.17 3.21
	20 H 2 O in all-silica LTA 2.83	-3.18 -0.14 3.32
	20 H 2 O in Na 1 -LTA 2.83	-3.14 -0.12 3.27
	20 H 2 O in Na 9 -LTA 2.88	-3.06 -0.12 3.18
	1 H 2 O in Na 1 -LTA 2.38	-2.57 -0.26 2.82
	isolated water molecule 1.84	-2.47 -0.10 2.57
	bulk liquid water 3.0	-3.19 -0.19 3.38
	ice 3.44	-3.69 -0.12 3.81

Table 6 :

 6 Average partial charges on the oxygen and silicon atoms of the all-silica zeolite, depending on the analysis method used: Mulliken, Hirshfeld, ESP and Born (see Section 3.4).

	method oxygen silicon
	Mulliken -1.007 2.175
	Hirshfeld -0.293 0.586
	ESP -1.901 3.801
	Born -1.557 3.114

)

  Figure 11: OH stretch band of the infrared spectrum of 20 water molecules confined in all-silica LTA zeolite (upper panel) and 32 bulk water molecules (middle panel). The lower panel shows the Vis-IR SFG spectrum of water near an H 2 O/CCl 4 interface, from Scatena et al.[95], exhibiting a high-frequency peak considered as a signature of hydrophobicity. The simulated systems (upper and middle panel) are deuterated, hence the difference in vibration frequency.Figure12: Correlation between the average stretching frequency (horizontal axis) of a water OH group confined in an all-silica LTA and the percentage of simulation time during which this OH group is involved in water-water hydrogen bonds (vertical axis, red circles) and is oriented towards the zeolite (vertical axis, blue triangles).

	Fraction of time involved in H bonds	Confined water 2100 2200 OH stretch frequency (cm 2300 -1 ) 0% 20% 40% 60% 80% 100% Hydrogen bond with water OH group oriented towards zeolite	2400
		Bulk water
		3000	3200 IR Energy (cm 3400 -1 3600	3800

† DFT calculations were performed with the BLYP functional and Goedecker-Tetter-Hutter pseudopotentials in the GPW approach, with DZVP basis set and a cutoff of 300 Ry for the description of the electronic density in plane waves.

Acknowledgments

This work was supported by the French Agence Nationale de la Recherche, under contracts # JC05 51932 and BLANC06-3 144027. We acknowledge the Institut du Développement et des Ressources en Informatique Scientifique (IDRIS) and ENSCP for computer time allocation. We thank Isabelle Demachy, Florent Calvo and Nicolas Desbiens for fruitful discussions.


Microalgae, which constitute third generation biomass, allow the production of biofuels through various conversion routes. Because of the high amount of water in this biomass, hydrothermal liquefaction is particularly appropriate for the conversion of microalgae into liquid fuels. Microalgae mainly consist of proteins, lipids and carbohydrates. Lipids are easily recoverable and are often extracted by different methods. Thus, the valorization described here can follow a first step of lipid extraction.

Hydrothermal liquefaction and partial oxidation of microalgae : influence of lipid content

Complementary to classic hydrothermal liquefaction of microalgae, an additional process is proposed by using partial oxidative conditions. This process may change the reaction pathways and increase the formation of some recoverable products such as organic acids.

Experiments are carried out in a batch reactor in the range 10-30 MPa, 473-623 K. Aqueous microalgae are introduced in the reactor at a known concentration. Then, the reactor is heated at the operating temperature and air may be added to adjust the oxygen content in the gas phase. The amount is calculated to represent from 0% to 30% of the total organic carbon. After the reaction, remaining phases are recovered and analyzed. In the aim to study the influence of the algae composition, experiments are performed with whole microalgae, same microalgae after lipid extraction, pure algae proteins and pure algae carbohydrates.

Presented results compare the yields in recoverable products during hydrothermal liquefaction with and without oxygen of biomass according to operating conditions. A particular attention is paid to the influence of microalgae composition. Moreover, the analysis of the products is used to understand the reaction mechanisms in order to model and simulate the process in further studies.

02116839 Hydrothermal liquefaction and partial oxidation of microalgae

	Tsilla Bensabath, Olivier Boutin, Jean-Henry Ferrasse

Tsilla Bensabath, Olivier Boutin, Jean-Henry Ferrasse To cite this version: Tsilla Bensabath, Olivier Boutin, Jean-Henry Ferrasse. Hydrothermal liquefaction and partial oxidation of microalgae : influence of lipid content. EUBCE 2018, May 2018, Copenhague, Denmark. hal-Aix-Marseille Université, CNRS, ECM, M2P2 UMR 7340, FR-13451 Marseille, France