Bennici Simona
email: simona.bennici@uha.fr

Polimann Téo

Ondarts Michel

Evelyne Gonze

Vaulot Cyril

Le Pierrès Nolwenn

Long-term impact of air pollutants on thermochemical heat storage materials

Keywords: thermochemical storage, pollution, energy recovery, solar energy, zeolites, material characterisation H, relative humidity, FID, flame ionisation detector, FTIR, Fourier-transform infrared spectroscopy, WDXRF, wavelength dispersive x-ray fluorescence, XPS, x-ray photoelectron spectroscopy, BE, binding energy, BET, Brunauer Emmett and Teller, SSA, specific surface area, MAS, magic angle spinning, NMR, nuclear magnetic resonance, TG-DSC (thermogravimetry-differential scanning calorimetry), IR, infrared, S BET , BET specific surface

Heating of buildings is a highly energy-demanding task. Therefore, improving energy management is crucial for more environmentally friendly development of future constructions, and the use of thermochemical heat storage technologies in solid materials is one of the most promising solutions. Nevertheless, certain drawbacks such as the impact of air pollutants on the long term durability of thermochemical heat storage materials (ageing up to 30 years) need to be studied to implement installations. In certain thermochemical heat storage systems, air passes through porous materials to carry water and heat, and the air's pollutants can interact with the solid material, thus decreasing its storage capacity.

In the present work, adsorption tests were performed (under dynamic conditions) on zeolitebase materials, using different model molecules (i.e. toluene, styrene and hexaldehyde) representative of air pollutants. Strong competition between the adsorption of pollutants and water on the storage material was demonstrated. Water molecules were able to displace the molecules of pollutants, previously adsorbed on the material surface, thus delaying the saturation of the material and extending its life. Nevertheless, a lowering of the water adsorption and heat storage capacity was measured for a high number of cycles. The cycled materials have been thoroughly characterised in their physicochemical properties to identify the modifications of the materials (i.e., de-alumination) and correlate them with the heat storage capacity. The impact of pollutants was demonstrated to be stronger on the pure zeolite, while the salt present on the MgSO 4 /13X composite seems to protect the material surface from the pollutant effect.

Highlights

-The zeolite-base materials are highly stable in hydration/dehydration cycles -The heat storage capacity decreases in presence of pollutants -Pollutant molecules are displaced by water molecules during hydration -The microstructure of the material is modified during cycles

To decrease the dependence on fossil energy sources (expensive, limited and polluting), the solution is to orient development towards the use of renewable energies, which can be implemented in more environmentally friendly solutions [START_REF] Aydin | The latest advancements on thermochemical heat storage systems[END_REF][START_REF] Yan | A review of promising candidate reactions for chemical heats storage[END_REF]. Solar energy is the most abundant renewable energy source. Unfortunately, solar energy is intermittent, and, consequently, the need to bridge the gap between energy supply and energy demand remains a challenge. For this reason, identifying, developing, and optimising the energy storage technologies is the entryway to a renewable energy-based economy [START_REF] Yan | A review of promising candidate reactions for chemical heats storage[END_REF][START_REF] Kuznik | A review on recent developments in physisorption thermal energy storage for building applications[END_REF]. Solar power peak availability occurs during the day and summer, periods during which the heating demand in buildings (for sanitary hot water and building heating) is at its lowest. Storage systems are therefore necessary to store the excess heat produced during sunny periods, to be used during the periods of highest demand [START_REF] Shah | Seasonal thermal energy storage system for cold climate zones: A review of recent developments[END_REF]. Thermal energy storage can be applied by means of different specific technologies: water, phase-change materials, and thermochemical reactions [START_REF] Yan | A review of promising candidate reactions for chemical heats storage[END_REF][START_REF] N'tsoukpoe | A review on long-term sorption solar energy storage[END_REF]. These technologies can be based on different physicochemical phenomena. The most widely used systems are based on the exploitation of sensible heat, linked to the heat exchanged due to changes in the system's temperature [START_REF] Hadorn | Solar and Heat Pump Systems -Status of the IEA SHC Task 44 & HPP Annex 38[END_REF]. Currently, water is widely used for sensible heat storage; however, for building heating applications and long-term storage, excessively large volumes are required, and its application is unrealistic for a small family home, for example [START_REF] Van Essen | Characterization of MgSO 4 hydrate for thermochemical seasonal heat storage[END_REF]. The second solution is based on latent heat storage; in this case the heat stored/released by the material stems from a phase change. The third way, which is considered in the present study, is thermochemical heat storage, based on reversible chemical reactions (for example hydration/dehydration of a material) [START_REF] Yan | A review of promising candidate reactions for chemical heats storage[END_REF][START_REF] Kuznik | A review on recent developments in physisorption thermal energy storage for building applications[END_REF]. Systems based on this principle have the advantage of presenting negligible heat loss during storage and having a higher heat density than the other storage systems mentioned above. The main advantage of using this type of system is the possibility of storing heat during the period of energy availability (e.g. solar radiation) to perform the charging step (i.e. dehydration of the solid material) and release it (i.e. hydration of the solid) for utilisation in the desired application (house heating), as shown in Fig. 1.

Fig. 1. Principle of the adsorption heat storage system functioning (a) charging step: endothermal dehydration of the adsorbent using air heated by a solar collector; (b) discharging step: heating of the ventilation air by exothermal hydration of the adsorbent

The solid storage material is at the heart of such thermochemical heat storage systems, and it must possess certain specific properties [START_REF] Tatsidjodoung | A review on long-term sorption solar energy storage[END_REF], such as high energy density, high affinity for the sorbate, water (to ensure the exothermic chemical reaction), high thermal conductivity (to facilitate heat transfer and recovery), a high heat transfer coefficient with the heat transfer fluid, and a desorption (dehydration) temperature as low as possible. Moreover, the storage materials must be environmentally friendly, non-toxic, non-corrosive, inexpensive, and present good thermal and chemical stability under operating conditions (temperature, pressure, etc.). Various authors have explained that the main disadvantage of such materials is the limitation in heat and mass transfer [START_REF] Tatsidjodoung | A review on long-term sorption solar energy storage[END_REF]. The main limitation in solid adsorbents is the low thermal conductivity of sorbents and the possibility of crust formation (in case salts are used), which limits water diffusion. Despite these drawbacks, several sorption couples have been identified and studied for the implementation of the process. Three main families of porous materials (mesoporous silicates, zeolites and metal aluminophosphates) have been widely reported in the literature [START_REF] Yan | A review of promising candidate reactions for chemical heats storage[END_REF][START_REF] Kuznik | A review on recent developments in physisorption thermal energy storage for building applications[END_REF]. Among them, zeolites are the best choice for large-scale and long-term applications, due to their high heat storage capacity and durability [START_REF] Jänchen | Shaping adsorption properties of nano-porous molecular sieves for solar thermal energy storage and heat pump applications[END_REF][START_REF] Chen | Catalytic decomposition performance for O 3 and NO 2 in humid indoor air on a MnO x /Al 2 O 3 catalyst modified by a costeffective chemical grafting method[END_REF].

In certain types of thermochemical heat storage installations (such as that considered in the present study), air passes through the porous material to carry water (the sorbate) and heat (Fig. 1). It is therefore necessary to study the impact of air pollutants (indoor and outdoor air) [START_REF] Kirchner | Etat de la qualité de l'air dans les logements français[END_REF][START_REF] Sundell | On the history of indoor air quality and health[END_REF][START_REF] Edwards | VOC source identification from personal and residential indoor, outdoor and workplace microenvironment samples in EXPOLIS-Helsinki, Finland[END_REF] on the ageing of these materials. In the long term (functioning over 30 years), these pollutants can interact with the solid and eventually decrease the material's storage capacity.

Outdoor and indoor air pollutants are mainly composed of solid particles, inorganic (ozone, nitrogen oxides) and organic pollutants (volatile organic compounds and carbonyl compounds). Their concentrations vary from a few µg.Nm -3 up to a few hundred µg.Nm -3 [START_REF] Sundell | On the history of indoor air quality and health[END_REF][START_REF] Edwards | VOC source identification from personal and residential indoor, outdoor and workplace microenvironment samples in EXPOLIS-Helsinki, Finland[END_REF]. The present study focuses only on chemical pollutants. Moreover, zeolites are known to be good adsorbents for atmospheric pollutants [START_REF] Brosillon | Mass Transfer in VOC Adsorption on Zeolite: Experimental and Theoretical Breakthrough Curves[END_REF][START_REF] Lee | Adsorption and thermal regeneration of acetone and toluene vapors in dealuminated Y-zeolite bed[END_REF], and their affinity for these molecules might be a drawback for their application in real thermochemical heat storage systems.

In this study, the impact of air pollutants on the durability of thermochemical heat storage materials has been examined, because, to our knowledge, no research concerning the impact of pollution on these systems and their heat storage performance has been reported to date in the literature. Critical pollutants were chosen and ageing tests were performed to model the storage material implementation in a long-term application and in successive hydration/dehydration cycles. Adsorption tests (performed under static and dynamic conditions) were carried out on zeolite-base materials (13X pure zeolite and a MgSO 4 /13X composite), using different molecules representative of the pollutants (toluene, styrene and hexaldehyde), chosen as model molecules.

Then the cycled materials were deeply characterised to shed light on the changes in their properties and correlate them with the loss in heat storage capacity.

Material and methods

Zeolite, composite and pollutants

Commercial 13X zeolite (supplied by Chemiewerk Bad Köstritz GmbH, brand Köstrolith 13XBF K distribution), shaped in spheres 1.2 to 2.0 mm in diameter, was chosen as the sorbent solid material, the sample was labelled 13X. This type of zeolite has been widely studied in thermochemical heat storage processes [START_REF] Helaly | Theoretical and experimental investigation of the performance of adsorption heat storage system[END_REF][START_REF] Henninger | New materials for adsorption heat transformation and storage[END_REF].

A second sample, a composite, was prepared by depositing 15 wt% MgSO 4 (99.5% purity, purchased from Sigma-Aldrich) on the 13X zeolite, used as the support, using the impregnation method under vacuum, with the same procedure and the same vacuum apparatus as reported in Ref. [START_REF] Chen | Catalytic decomposition performance for O 3 and NO 2 in humid indoor air on a MnO x /Al 2 O 3 catalyst modified by a costeffective chemical grafting method[END_REF], and labelled MgSO 4 /13X. A large comparison of the possible thermochemical heat storage materials (mainly salts hydrate) was carried out at ECN (Energy research Center of the Netherlands and reported in 2004 [START_REF] Visscher | Compacte Chemische Seizoenopslag Van Zonnewarmte; Eindrapportage[END_REF] whereby MgSO 4 was identified as a potential material for long-term heat storage thanks to its theoretical high storage density of 780 kWh.m -3 . During hydration/dehydration, MgSO 4 can be present in different hydrated forms [START_REF] Van Essen | Characterization of MgSO 4 hydrate for thermochemical seasonal heat storage[END_REF][START_REF] Whiting | Heats of water sorption studies on zeolite-MgSO4 composites as potential thermochemical heat storage materials[END_REF] and under practical conditions it cannot release the total theoretical amount of heat stored. Nevertheless, even if not completely hydrated/dehydrated, magnesium sulphate represents a good choice for heat storage, when its dispersion on a support allows to overcome the diffusion water issues present in the bulk salt [START_REF] Visscher | Compacte Chemische Seizoenopslag Van Zonnewarmte; Eindrapportage[END_REF][START_REF] Whiting | Heats of water sorption studies on zeolite-MgSO4 composites as potential thermochemical heat storage materials[END_REF][START_REF] Hongois | Development and characterisation of a new MgSO 4 -zeolite composite for longterm thermal energy storage[END_REF]. Moreover, MgSO 4 is also non-toxic, relatively cheap and non-corrosive, which are very important properties for a material that might be used at large scale in households' applications.

The possibility to efficiently deposit MgSO 4 on zeolites was already demonstrate in previous researches [START_REF] Whiting | Heats of water sorption studies on zeolite-MgSO4 composites as potential thermochemical heat storage materials[END_REF]. In the same research, loading 15wt% of MgSO 4 on the zeolites was shown to be the best compromise to obtain the highest increasing of the heat storage capacity and preserve the zeolite porosity.

Three representative air pollutants were selected on the basis of Ref. [START_REF] Kirchner | Etat de la qualité de l'air dans les logements français[END_REF]: toluene (99.9% purity, Fisher Scientific, Leicestershire, UK), styrene 99.5% purity, Acrôs Organics) and hexaldehyde (96% purity, Acrôs Organics). The three pollutants have been selected among the others present in indoor air due to their different physico-chemical properties and chemical reactivity. Toluene is a synthetic intermediate for the manufacture of many products: benzene, xylenes, phenol, etc. It is also largely present in solvents for paints, varnishes, printing inks, glues, waxes, etc. Toluene is therefore present in large quantities in indoor air. Moreover, it is almost not soluble in water. Lee et al. [START_REF] Lee | Adsorption and thermal regeneration of acetone and toluene vapors in dealuminated Y-zeolite bed[END_REF] enlightened the possible adsorption of toluene on zeolite. In this study, industrial effluents loaded with toluene and acetone were treated with a Y-type zeolite. Toluene desorption was obtained in the 150-180°C interval, indicating the formation of chemical bonds between toluene and zeolite. Knowing that in the solar storage process the expected temperature in the charge phase is 165°C, the study of toluene will allow to verify if the charging temperature is higher enough to avoid the accumulation of this compounds on the storage material. Styrene is present in small quantities in most indoor environments due to its desorption from plastic materials, glues, paints, solvents, etc. It is not miscible with water as toluene, but present the particularity of being higly reactive: it polymerises and oxidizes very easily. The polymerization reaction is accelerated by the action of light and heat. Barson et al. [START_REF] Barson | The polymerisation of styrene on 13X zeolite[END_REF] showed that starting from 30°C styrene polymerizes on zeolite 13X. For this reason, styrene has been chosen as model pollutant to investigate the possible blocking of the zeolite pores due to the formation of polymerisation products and sub-products. Finally, hexanal is an organic compound of the aldehydes family. It is present in most perfumes and aromas, because it is used in the extraction process. Hexanal was selected for this study because is a persistent pollutant present in relatively high quantities in indoor air, and because, being very miscible with water, it might reach in depth the micropores of the storage material during the adsorption phase.

Pollutant adsorptionbreakthrough curves and cycling tests

The acquisition of the breakthrough curves of the three pollutants on the zeolite and composite, as well as hydration/dehydration cycles, were performed in continuous flow mode.

These breakthrough curves are widely applied methods to evaluate the performance of adsorption storage systems, as described in [START_REF] Hongois | Development and characterisation of a new MgSO 4 -zeolite composite for longterm thermal energy storage[END_REF] or [START_REF] Kuznik | Numerical modelling and investigations on a full-scale zeolite 13X open heat storage for buildings[END_REF]. In our case, dry (relative humidity H < 1%) and humidified air (H ≈ 60%) flows) were used to represent the charging and discharging steps of the real system. The experimental set-up used is shown in Fig. 2.

Fig. 2. Scheme of the experimental set-up.

The carrier gas (atmospheric compressed air) was purified and dried using an oil separator, a desiccant cartridge, an activated carbon filter and a particulate filter. In the experiments conducted in wet conditions, as well as for the hydration steps during cycling, a fraction of the air flow was saturated passing through a temperature-controlled water bubble column (Variocool VC 1200 LAUDA). The dry air and humidified air were then mixed to obtain the targeted humidity value. The dry and humidified air flows were controlled by mass flowmeters (from Brooks, working in the 0-30 NL.min -1 range). The humidity value was measured continuously downstream by means of a D2 General Eastern Chilled Mirror Hygrometer. The parameters were followed on-line on a monitor control implemented by means of a Siemens Sipart DR 20 S/K controller. The temperature was measured upstream and downstream of the reactor by means of Pt100 sensors and K-type thermocouples.

The pollutants were introduced in the multi-reactor device by mixing the main air flow with a nitrogen flow saturated with the pollutants' molecules, to reach the targeted pollutant concentration of 100 mg.Nm -3 in a total flow of 1 L.min -1 . The nitrogen flow was saturated with each pollutant in a temperature-controlled bubbled column (at 25°C). The flow rate in each column was controlled by a GFC thermal mass flow controller (Aalborg). The concentration of the pollutant at the inlet and outlet of the reactors was measured by means of a flame ionisation detector (FID9000 Model, EPA) with a 1-min acquisition interval.

Twenty-four reactors were fed at the same time with the air/pollutant mixture. Each reactor was filled with 1 ± 0.02 g of sorbent (zeolite 13X or MgSO 4 /13X composite). The use of 24 parallel and identical reactors made it possible to sample the material contained in one of the reactors while continuing the ageing (cycling) test on the others.

The following temperature program was applied during the cycling tests:

1. Hydration: isothermal step at 30°C for 90 min; 2. Dehydration: increasing temperature step (at 4°C/min rate) up to 165°C followed by an isothermal step at 165°C for 30 min; 3. Cooling/pre-adsorption step, decreasing the temperature step to 30°C in about 90 min. These temperatures were chosen to represent as well as possible the real operating conditions of the system in a building: 165°C it is a typical temperature that can be reached by evacuated tube solar air collectors, as explained in [START_REF] Kuznik | Numerical modelling and investigations on a full-scale zeolite 13X open heat storage for buildings[END_REF][START_REF] Tatsidjodoung | Experimental and numerical investigations of a zeolite 13X/water reactor for solar heat storage in buildings[END_REF], and the 30°C temperature is close to the building inside air temperature during the discharging phase, and was chosen as it can be stabilised precisely from an experimental point of view. The 4°C/min rate was chosen for the temperature increase as it is of the same order of magnitude as the value that could be observed in a realistic system, as shown in [START_REF] Tatsidjodoung | Procédé de stockage d'énergie solaire thermique par adsorption pour le chauffage des bâtiments : Modélisation et Simulation numérique[END_REF].

Thirty-three cycles were performed on each sample with humid air and humid air/pollutant flows, to verify the impact of cycling and of the presence of pollutants on the storage materials. The number of cycles corresponded to the amount of pollutant that would pass through the reactor during different utilisation periods in a real application. Assuming an air flow estimated for a real installation for heating a single house (≈ 500 000 m 3 /year) and the average air pollutant concentration reported in Table 1, one cycle in the multi-reactor device corresponds to 19 months use in a real installation, five cycles to 8 years, 10 cycles to 16 years and 33 cycles to 52 years use.

Heat and water storage capacity measurements of the fresh and aged samples

The heat and water adsorption/desorption capacities of the fresh and cycled samples (in presence or absence of pollutants) were recorded using a Sensys TG-DSC (Thermogravimetry coupled to differential scanning calorimetry) apparatus, equipped with a Wetsys flow humidity generator, both from Setaram. Prior to the hydration process, the samples were dehydrated at the same temperature (165°C) as that applied during the cycles performed in the multi-reactor device described in section 2.2, by increasing the temperature from 30 to 165°C at 5°C min -1 under dry air flow (30 mL min -1 ), followed by an isothermal step of 60 min at 165°C. The sample was then cooled at 30°C and hydrated under humid air (H = 60%). The hydration heat (kJ g dry sample -1 ) of the various samples was obtained by integration of the surface of the exothermic calorimetric peak obtained during hydration. The calorimetric and thermogravimetric signals were treated by mean of AKTS Advanced Thermokinetics software.

Physicochemical characterisations

Fourier-transform infrared spectroscopy (FTIR) was performed on the powders of the 13X and MgSO 4 /13X samples by means of an FTIR iN10 MX (Thermo Scientific) analytical system.

Wavelength dispersive x-ray fluorescence (WDXRF) was performed in a Philips MagiX (2.4-kW) spectrometer to evaluate the Si/Al ratio in the fresh and cycled samples. Prior to analysis, the samples were shaped in pellets 13 mm in diameter and 1 mm thick by pressing the powders five 5 min, under 7-tonne pressure.

X-ray photoelectron spectroscopy (XPS) measurements were taken on a VG Scienta SES 2002 spectrometer equipped with a monochromatic Al Kα X-ray source (Al Kα = 1486.6 eV)

and a hemispherical analyser. The XPS spectra were recorded with a pass energy of 100 eV and at 10 -9 Pa pressure, in the analysis chamber; 24 mm 2 was the surface of the spot analysed and 9 nm the penetration depth of the x-ray photoelectrons. The binding energies (BEs) were calibrated on the basis of the C 1s peak (284.7 eV) of the atmospheric contaminant carbon.

The peaks were fit using Gaussian-Lorentzian functions by means of the XPS-CASA software.

The surface area and porosity of the fresh and cycled samples were studied with a Micrometrics ASAP 2020 instrument by adsorbing N 2 at -196°C. CO 2 adsorption at 0 °C was also performed to finely investigate the microporosity. Prior to analysis, the samples were first pretreated at 130°C under vacuum for 12 h on the degas port and then again at 130°C for 2 h under vacuum on the analysis port. The free volumes were measured separately after analysis to avoid the pollution of the microporosity by He. The specific surface area (SSA) was calculated applying the BET (Brunauer, Emmett and Teller) theory [START_REF] Brunauer | On a theory of the van der Waals adsorption of gases[END_REF] (within a range determined by the Rouquerol curve [START_REF] Thommes | Physisorption of gases, with special reference to the evaluation of surface area and pore size distribution (IUPAC Technical Report)[END_REF]), while the micropore volume (V micro ) was calculated by applying the t-plot method [START_REF] Lowell | Ch.9 Microporosity in "powder surface area and porosity book[END_REF]. Finally, the pore size distribution was determined by applying the DFT model and considering the N 2 adsorption on an oxide surface and according to a cylindrical pore geometry. 27 Al magic angle spinning (MAS) nuclear magnetic resonance (NMR) spectroscopy spectra were recorded on a Bruker AVANCE II 400 MHz spectrometer.

The 13X samples (fresh and cycled 33 times in presence of toluene) were analysed in a 2.5mm-diameter cylindrical rotor at a 25-kHz spinning frequency.

Results and discussion

The aim of this study was to determine the impact of pollutant adsorption on the performance of the storage system considered during the building's lifetime (at least 30 years of use), and to understand the cause of this impact by identifying the structural and morphological modifications that can occur on the solid storage material.

At first the adsorption capacity of 13X zeolite was determined for the three pollutants selected (toluene, styrene and hexanal) by acquiring the breakthrough curves, under dry and wet conditions (Table 2). These preliminary experiments were carried out only on the bare zeolite to verify the pollutants adsorption behaviour in presence of water vapour. In dry conditions, the zeolites presented the highest adsorption capacity for hexanal, when compared to styrene and toluene; in wet conditions the total adsorption of all pollutants was drastically decreased by displacement phenomena: zeolite 13X shows a higher affinity for water that can then replace the pollutants' molecules, thus reducing their adsorption on the storage material [START_REF] Tao | Effect of moisture on the adsorption of volatile organic compounds by zeolite 13X[END_REF][START_REF] Azambre | Adsorption and Desorption of a Model Hydrocarbon Mixture Over HY Zeolite Under Dry and Wet Conditions[END_REF]. Figure 3 shows this behaviour for styrene (reported here as an example). At first, when all adsorption sites are free, the pollutant and the water molecules adsorb on the adsorbent (13X) at the same time. Then, when the number of free sites decreases, competition between styrene and water takes over, and water molecules displace the previously adsorbed styrene molecules, replacing them. This phenomenon is visualised by the light blue peak in Fig. 3, generated by a higher concentration of styrene on the out-stream flow (C out ) due to the release of the molecules of the pollutant displaced. The same behaviour was observed for the molecules of the three pollutants tested. This displacement behaviour has been previously reported in the literature [START_REF] Azambre | Adsorption and Desorption of a Model Hydrocarbon Mixture Over HY Zeolite Under Dry and Wet Conditions[END_REF][START_REF] Gong | A qualitative analysis of the effects of water vapor on multicomponent vapor-phase carbon adsorption[END_REF] and is called the "roll-over" phenomenon.

Fig. 4. Hydration calorimetric curve for the 13X sample at 30°C and H = 60%. The surface in yellow represents the total hydration heat (in kJ•kg -1 dry sample ). The blue curve represents the average hydration energy (in kJ•mol H2O ) at the different time intervals.

Even if the adsorption of the pollutants on the zeolite-based materials is strongly decreased due to the presence of water vapour during the hydration/dehydration cycles applied to the thermochemical storage system, the long-term impact (over several years' utilisation) needs to be investigated. Cycling tests were therefore performed under the experimental conditions detailed in Section 2.2 by means of the multi-reactor device shown in Fig. 2. Then the heat and water adsorption capacities of the fresh and cycled samples were measured in the TG-DSC apparatus as described in Section 2.3. The typical hydration curve is reported in Fig. 4 (for the 13X fresh sample, shown as an example). The integration of the surface under the curve (in yellow) represents the heat released (in kJ•kg -1 dry sample ). The results are compared in Fig. 5a,b, respectively, for the 13X and MgSO 4 /13X samples cycled 1 and 33 times. In the case of the bare zeolite 13X, the heat of hydration of the sample cycled in presence of pollutants was always lower than that measured for the fresh 13X sample (see Fig. 5a).

Moreover, by increasing the number of cycles (from 1 to 33), the heat storage capacity (expressed as kJ•kg -1 dry sample ) decreased further. A different behaviour was observed for the 13X sample cycled in presence of only water vapour (with no pollutant); in this case the zeolite presented higher hydration heat values. In contrast, even if the heat storage capacity of the MgSO 4 /13X fresh sample, measured by calorimetry during hydration, was higher than that of the cycled samples, the material cycled 33 times showed a higher hydration heat than that cycled only once in presence of pollutants. The pollutant molecules seem to first strongly adsorb on MgSO 4 /13X and then be displaced by the water molecules during the successive cycles, thus "reactivating" the sample that regains its initial heat storage capacity.

As a general trend, the negative effect of the pollutants on the heat storage capacity was in the styrene > hexanal > toluene order, on both samples. The decrease in heat storage capacity is directly correlated with the decrease in water sorption (Fig. 6a,b). The water molecules are adsorbed with similar energy on the sample whatever the coverage extent is, as illustrated in The samples were then characterised in detail to provide an understanding of the modification in the water sorption and heat capacity of the two samples after cycling in presence or absence of pollutants.

At first, the possible presence of pollutants accumulating on the samples' surface was checked by IR spectroscopy and XPS analysis. No presence of the pollutants on the samples was detected. The FTIR spectra acquired for the 13X sample after cycling did not significantly differ from that of the fresh sample. In particular, neither stretching nor bending vibration peaks of the -CH 2 and -CO groups [START_REF] Ozturk | Infrared and SEM analyses of polyethyleneglycol-400 adsorbed on zeolites NaA, CaA, NaX and NaY[END_REF] of the adsorbed pollutants was observed, and the infrared spectra of the cycled 13X zeolites reported only the vibrations of the bare zeolite. The absence of the accumulation of pollutants on the samples' surface was also confirmed by XPS. Once again, the surface atomic composition of the fresh and cycled (in presence of styrene) 13X samples, did not differ. In particular, no accumulation of carbon-containing species was measured.

The decrease of the heat capacity and water up-take upon cycling is apparently not related to the deposition of pollutants on the heat storage material surface. The microstructure of the fresh and cycled 13X and MgSO 4 /13X samples were then investigated by analysis of the N 2adsorption isotherms acquired at -196°C. The main results are reported in Table 3 in terms of BET specific surface area (S BET ), micropore specific surface area (S micro ), total pore volume (V tot ) and micropore volume (V micro ), for all samples analysed. The N 2 -adsorption isotherms are reported in Fig. 7 for the 13X samples cycled from 1 to 33 times in presence of pure water vapour (a), styrene (b) and for the MgSO 4 /13X composite cycled in presence of styrene (c).

Figure 8 represents the associated micropore distribution, and Fig. 9 the meso/macropore distribution.

The microstructure of zeolite 13X was clearly preserved during cycling with clean water vapour, as shown by the very low variation in the surface area and porous volume values. The associated isotherms (Fig. 7a) are all type I [START_REF] Sing | Reporting physisorption data for gas/solid systems with special reference to the determination of surface area and porosity[END_REF], typical of microporous samples with small external surface areas, for which the limiting nitrogen uptake is driven by the micropores' accessibility. The micropores' distribution (Fig. 8a) did not vary by cycling the sample in presence of pure water vapour, and no macroporosity was detected. The behaviour upon cycling differed in presence of pollutants in the air flow: the specific surface area and porous volume of the samples decreased with the number of cycles (Table 3). This behaviour was not related to modification of the micropores distribution (determined by analysis of CO 2 adsorption isotherms) that did not vary for the samples cycled several times. In all cases, the type-I N 2 -adsorption isotherm characterising the 13X zeolites changed towards type-IV isotherms (Fig. 7b for cycling in presence of styrene), with pronounced hysteresis, which can be associated with nitrogen capillary condensation into mesopores, typical of many industrial mesoporous adsorbents. These observations are also confirmed by the pore size distribution determined by analysis of the N 2 -adsorption isotherm hysteresis (reported as an example in Fig. 8b and 9b for the cycles performed in presence of styrene), where a clear shift towards larger pores is observed. In particular, the appearance of meso/macroporosity is shown in Fig. 9b. The specific surface area of the MgSO 4 /13X composite was much lower than that of the parent zeolite, due to the salt deposition that, blocking the micropores' access, diminished the access and adsorption of nitrogen molecules, as already reported in the literature [START_REF] Kirchner | Etat de la qualité de l'air dans les logements français[END_REF]. After the first cycle, the salt distribution and arrangement probably was modified due to the hydration/dehydration phenomena letting supposing that the new salt space configuration free the pore access; then nitrogen could reach the pores and the porous volume, and the S BET measurement significantly increased. The MgSO 4 /13X sample was less affected by the presence of pollutants (Table 3) and all N 2 -adsorption isotherms of the cycled samples were type I; no creation of macro/mesoporosity was observed (Fig. 9). Only the presence of styrene impacted the impregnated zeolite at a higher level, for which the specific surface and micropore volume decreased from 489 to 258 m 2 •g -1 , respectively, and from 0.17 to 0.09 cm 3 •g -1 between cycle 1 and cycle 33 (Table 3). The modification in the samples' microstructure (in particular for the 13X sample upon cycling) could not be assigned to pollutant accumulation on the sample surface, as shown by XPS and FTIR analysis. Then the bulk sample composition was analysed using WDXRF for the 13X fresh and cycled samples. This investigation was inconclusive and a constant value of 1.25 for the Si/Al ratio was measured for all samples. Indeed, even if in-situ de-alumination takes over during cycling, no mass transfer is involved in the process and the extracted aluminium might remain in the sample matrix, thus resulting in an unmodified Si/Al ratio.

To finally verify the true chemical environment of aluminium atoms in the 13X zeolite, 27 Al NMR spectra were acquired for the fresh sample and the sample cycled 33 times in presence of toluene (Fig. 10a, b, respectively). A main peak at an approximately 60-ppm chemical shift was visible on both samples and assigned to tetrahedrally coordinated aluminium atoms (Al(OSi) 4 ) [START_REF] Muller | Quantitative multinuclear MAS NMR studies of zeolites[END_REF], connected to the main zeolite structure. In case of de-alumination, octahedrally coordinated aluminium is present in extra-framework positions in the zeolite matrix. In the 13X zeolite cycled 33 times in presence of toluene (13X-Tol sample), a weak peak at around 3 ppm confirms the presence of de-alumination, even if to a small extent.

Slow de-alumination undertaken by cycling in presence of water vapour and pollutants therefore seems to be the cause of the modification of the microstructure of the 13X sample and to the related loss in heat storage capacity and water uptake. The presence of MgSO 4 appears to protect the zeolite from de-alumination, and hydration/dehydration cycling on the MgSO 4 /13X sample, even in presence of pollutants, does not have a strong effect on the sample's microstructure, probably because the pollutants cannot reach the internal pores, due to blockage of the pores' access by the salt deposited.

Conclusions

Zeolite-based heat storage materials are confirmed to be a very good option for open sorption thermochemical heat storage systems, even when the carrier gas (air) is polluted by organic molecules such as toluene, styrene and hexanal. The strong affinity of zeolite 13X and MgSO 4 /13X composite for water allows displacement processes of the pollutant molecules that help to avoid the fast poisoning of the solid material, thus enhancing their durability (in terms of maintaining the heat capacity). Nevertheless, during several years' utilisation, hydration/dehydration cycling in presence of pollutants can modify the microstructure in the zeolite framework, due to slow and long-term de-alumination phenomena, and the consequent decrease in the heat storage capacity and water sorption. The impact of pollutants has been demonstrated to be stronger on the pure zeolite, while the deposited salt (MgSO 4 ) seems to protect the material surface from the pollutant effect, and no dramatic microstructural modifications could be detected on the salt-containing composite.

The choice of pure zeolites or composites (prepared by impregnating the zeolites with various salts, such as MgSO 4 as in the present research) for implementation in a thermochemical heat storage device must then be made not only by comparing the intrinsic heat storage capacity of the materials measured in ideal conditions (absence of pollutants), but taking into account the average concentration of organic molecules in the air stream that will be used in the actual installation. In ideal conditions, bare zeolite remains the best choice in terms of total heat storage, but addition of an optimised amount of hydrated salts can improve the durability of the storage material and the long-term efficiency of the thermochemical heat storage systems.

Funding: -Région Auvergne Rhône Alpes for the PhD thesis grant of T. Polimann -Région Grand Est for providing funding for the acquisition of the TG-DSC equipment within the "STOCKFATAL project.

-Carnot MICA and Energies du futur Institutes, for funding the Projet InterCArnot "StockEner".

Fig. 4 Fig. 5 .

 45 Fig. 4 for the 13X sample. The energy of water molecule adsorption on the 13X sample at

Fig. 6 .

 6 Fig. 6. Water sorption capacity for the 13X (a) and MgSO 4 /13X (b) samples cycled 1 and 33

Fig. 7 .

 7 Fig. 7. N 2 -adsorption isotherm for the 13X samples cycled from 1 to 33 times in presence of

Fig. 8 .

 8 Fig. 8. Micropore for the 13X samples cycled from 1 to 33 times in presence of pure water

Fig. 9 .

 9 Fig. 9. Meso/macropore for the 13X samples cycled from 1 to 33 times in presence of pure

Fig. 10 .

 10 Fig. 10. 27 Al MAS NMR spectra of (a) 13X and (b) 13X cycled 33 times in presence of toluene

Table 1 .

 1 Examples of organic and inorganic indoor/outdoor air pollutants and typical concentrations expressed in µg of pollutant per m 3 of air

	Indoor pollutants	Outdoor pollutants
	Pollutant	Concentration	Ref. Pollutant	Concentration	Ref.
		(µg.m -3 )		(µg.m -3 )
	Formaldehyde	19.4	[11] Benzene	0.3	[28]
	Hexaldehyde	17.0	[26] Ethylbenzene	1.9	[29]
	Nitrogen dioxide	32.5	[26] Nitrogen dioxide	15.0	[28]
	Styrene	2.1	[27] Sulfur dioxide	1.0	[28]
	Tetrachlorethylene	8.54	[27] Toluene	12.9	[26]
	Toluene	12.2	[26] Xylene	0.9	[29]

Table 2 .

 2 Pollutant adsorption capacity on 13X zeolite in dry and wet conditions

	Pollutant	Concentration	Adsorption capacity
		(mg•Nm -3 )	(mg•g -1 )
			H < 1%	H = 60%
	Hexanal	100	575	1.0
	Styrene	100	208	0.1
	Toluene	100	113	1.1

Fig. 3. Styrene/water breakthrough curves in wet conditions at 30°C. C out and C in , concentration at the exit and entrance of the reactor, respectively.

Table 3 .

 3 Microstructural properties of the fresh and cycled 13X and MgSO 4 /13X samples determined from N 2 -adsorption isotherms

	Sample
					13X						MgSO 4 /13X
		Cycle	S BET	a	S micro	b	V tot	V micro	Cycle	S BET	a	S micro	V tot	V micro	b
	Pollutant	n.	(m 2 •g -1 )	(m 2 •g -1 )	(cm 3 •g -1 )	(cm 3 •g -1 )	n.	(m 2 •g -1 )	(m 2 •g -1 )	(cm 3 •g -1 )	(cm 3 •g -1 )
	-	-	662	612		0.28	0.23	-	367	330	0.22	0.17
	Pure water	1	647	602		0.28	0.23
	vapour	4	670	629		0.28	0.24
		33	668	620		0.28	0.23
		1	424	415		0.18	0.18	1	517	481	0.23	0.18
	Toluene	5	318	272		0.22	0.22	5	466	433	0.21	0.16
		10	244	201		0.23	0.23
		33	179	139		0.17	0.17	33	499	462	0.23	0.17
		1	341	339		0.13	0.13	1	489	457	0.22	0.17
	Styrene	5	331	298		0.24	0.24	5	430	399	0.21	0.15
		10	360	313		0.24	0.24
		33	238	196		0.18	0.18	33	258	228	0.14	0.09
									1	492	457	0.22	0.17
	Hexanal								4	477	443	0.22	0.17
		5	401	346		0.23	0.23
		10	211		182		0.16	0.16
		33	179	139		0.17	0.17	33	439	405	0.19	0.15
	a Calculated by applying the BET equation [32]

b Calculated using t-plot method

[START_REF] Lowell | Ch.9 Microporosity in "powder surface area and porosity book[END_REF]

Acknowledgements

The authors would like to thank the Carnot Institutes MICA and Energies du futur (France) for supporting a part of this study within the STOCKENER project; Region Grand Est (France) for providing funding for the acquisition of the TG-DSC equipment in the STOCKFATAL project; and Region Auvergne Rhône-Alpes for the PhD thesis grant of T.

Polimann.

All physicochemical characterisations were performed on the IS2M technical platforms. The authors are very grateful to S. Hajjar (XPS) and L. Michelin (XRF) for their contribution.


Alessio Elia
email: alessio.elia@medaustron.at

Andreas Franz Resch

Antonio Carlino

Till Tobias Böhlen

Hermann Fuchs

Hugo Palmans

Virgile Letellier

Dreindl Ralf

Osorio Jhonnatan

Markus Stock

David Sarrut

Loïc Grevillot

A GATE/Geant4 beam model for the MedAustron non-isocentric proton treatment plans quality assurance

Keywords: GATE/Geant4, proton scanned beam delivery, Monte Carlo, beam modelling, non-isocentric treatment

Purpose: to present a reference Monte Carlo (MC) beam model developed in GATE/Geant4 for the MedAustron fixed beam line. The proposed model includes an absolute dose calibration in Dose-Area-Product (DAP) and it has been validated within clinical tolerances for non-isocentric treatments as routinely performed at MedAustron. Material and Methods: the proton beam model was parametrized at the nozzle entrance considering optic and energy properties of the pencil beam. The calibration in terms of absorbed dose to water was performed exploiting the relationship between number of particles and DAP by mean of a recent formalism. Typical longitudinal dose distribution parameters (range, distal penumbra and modulation) and transverse dose distribution parameters (spot sizes, field sizes and lateral penumbra) were evaluated. The model was validated in water, considering regular-shaped dose distribution as well as clinical plans delivered in non-isocentric conditions.

Results: simulated parameters agree with measurements within the clinical requirements at different air gaps. The agreement of distal and longitudinal dose distribution parameters is mostly better than 1 mm. The dose difference in reference conditions and for 3D dose delivery in water is within 0.5% and 1.2%, respectively. Clinical plans were reproduced within 3%.

Conclusion:

a full nozzle beam model for active scanning proton pencil beam is described using GATE/Geant4. Absolute dose calibration based on DAP formalism was implemented. The beam model is fully validated in water over a wide range of clinical scenarios and will be inserted as a reference tool for research and for independent dose calculation in the clinical routine.

Introduction

For proton beam delivery, the air gap between the treatment head window (nozzle exit window) and the patient can enlarge the lateral penumbra, as well as lead to lower accuracy of the pencil beam algorithm, especially if the Range Shifter (RaShi) is considered [START_REF] Soukup | A pencil beam algorithm for intensity modulated proton therapy derived from monte carlo simulations[END_REF][START_REF] Saini | Dosimetric evaluation of a commercial proton spot scanning Monte-Carlo dose algorithm: comparisons against measurements and simulations[END_REF]. One possible solution to reduce the air gap in clinical practice is to move the patient towards the nozzle exit.

Monte Carlo (MC) simulation tools are often used for independent dose calculation as they can provide useful insights that cannot be obtained otherwise (Tourovosky et al 2005[START_REF] Paganetti | Clinical implementation of full monte carlo dose calculation in proton beam therapy[END_REF][START_REF] Perl | Topas: An innovative proton monte carlo platform for research and clinical applications[END_REF][START_REF] Grevillot | GATE as a GEANT4-based Monte Carlo platform for the evaluation of proton pencil beam scanning treatment plans[END_REF], Magro et al 2015). However, any independent dose calculation engine should be carefully characterized and validated before its implementation in clinical applications. Typical characterization of a proton beam can be performed based on measured depth dose curves, spot maps measured at different air gaps and dose measurements in reference condition [START_REF] Parodi | Monte carlo simulations to support start-up and treatment planning of scanned proton and carbon ion therapy at a synchrotron-based facility[END_REF][START_REF] Clasie | Golden beam data for proton pencil-beam scanning[END_REF][START_REF] Grassberger | Characterizing a proton beam scanning system for Monte Carlo dose calculation in patients[END_REF]. In the literature, different methods for proton pencil beam modeling have been proposed depending on the types of the MC code and on the characteristics of the beam line [START_REF] Paganetti | Clinical implementation of full monte carlo dose calculation in proton beam therapy[END_REF][START_REF] Grassberger | Characterizing a proton beam scanning system for Monte Carlo dose calculation in patients[END_REF], Grevillot et al 2011[START_REF] Testa | Experimental validation of the topas monte carlo system for passive scattering proton therapy[END_REF][START_REF] Fracchiolla | Characterization and validation of a monte carlo code for independent dose calculation in proton therapy treatments with pencil beam scanning[END_REF].

For active scanning systems, accurate simulations of a proton pencil beam require a deep understanding of the pencil beam optics properties as the correlation between proton position and angular spread (beam emittance) cannot be neglected a priori [START_REF] Paganetti | Clinical implementation of full monte carlo dose calculation in proton beam therapy[END_REF]. Omitting the simulation of the nozzle geometry allows an empirical modeling of the pencil beam with an additional reward in computational time [START_REF] Grassberger | Characterizing a proton beam scanning system for Monte Carlo dose calculation in patients[END_REF], Grevillot et al 2011[START_REF] Fracchiolla | Characterization and validation of a monte carlo code for independent dose calculation in proton therapy treatments with pencil beam scanning[END_REF]. Nevertheless, such approach may require additional corrections [START_REF] Soukup | A pencil beam algorithm for intensity modulated proton therapy derived from monte carlo simulations[END_REF][START_REF] Schwaab | Experimental characterization of lateral profiles of scanned proton and carbon ion pencil beams for improved beam models in ion therapy treatment planning[END_REF] to account for large-angle single scattered primary and secondary particles produced in the nozzle components [START_REF] Grassberger | Characterizing a proton beam scanning system for Monte Carlo dose calculation in patients[END_REF], Pedroni et al 2005[START_REF] Gabriel O Sawakuchi | Monte carlo investigation of the low-dose envelope from scanned proton pencil beams[END_REF].

In this work, we provide a beam model of the MedAustron fixed beam lines using a full nozzle modeling approach that will be used as a reference for future independent dose calculation. Particular attention is given for the first time to the modeling and to the validation of the pencil beam properties in non-isocentric conditions as routinely exploited for patients treatments at MedAustron since December 2016. In addition, we exploit the Dose-Area-Product (DAP) formalism for absolute dose calibration in terms of delivered protons, differently from the traditional Monitor Unit (MU) calibration [START_REF] Fracchiolla | Characterization and validation of a monte carlo code for independent dose calculation in proton therapy treatments with pencil beam scanning[END_REF][START_REF] Paganetti | Monte carlo calculations for absolute dosimetry to determine machine outputs for proton therapy fields[END_REF], Mirandola et al 2015). The three-dimensional dose distribution simulated in water is carefully validated in terms of ranges, distal penumbra, modulation, field sizes and lateral penumbra, evaluated for different cubic targets in water, including non-isocentric conditions and the use of a range shifter (RaShi). Clinical plans are recalculated in GATE/Geant4 and compared to patient specific quality assurance data.

Materials and Methods

2.A. Beam Delivery System

The synchrotron-based Beam Delivery System (BDS) at MedAustron is, in clinical mode, configured to deliver protons at 255 different energies in a range between 62.4 and 252.7 MeV (range in water from 3 to 38 cm in 1 mm steps till 18.8 cm and 2 mm steps otherwise) [START_REF] Stock | The technological basis for adaptive ion beam therapy at medaustron: status and outlook[END_REF]. Using the slow extraction method [START_REF] Bryant | Proton-Ion Medical Machine Study (PIMMS)[END_REF], the particles are extracted from the synchrotron spill by spill and are guided by the High Energy Beam Transport (HEBT) line towards the Irradiation Room (IR) (spill length from 1 s to 5 s). The maximum number of particles extracted per spill is 2E10. A quasidiscrete spot scanning technique has been implemented for active delivery [START_REF] Haberer | Magnetic scanning system for heavy ion therapy[END_REF]. The nozzle geometry is shown in figure 1. The Water Equivalent Thickness (WET) of the full nozzle is 2.4 mm. Material composition and dimensions of each nozzle element were adopted from the manufacturer or measured wherever possible.

2.B. Experimental Data

Experimental data used for beam modeling were acquired during medical physics commissioning of the horizontal beam line (HBL) in the irradiation room 3 (IR3). The Full Width at Half Maximum (FWHM) was extracted from two dimensional lateral dose profiles measured with a scintillating screen detector (Lynx, IBA-dosimetry, Germany) oriented perpendicular to the beam axis.). Spot maps of 9 or 25 regularly spaced spots were delivered over a 20x20 cm2 field size. FWHM were averaged over the entire spot map and were determined at seven isocenter to detector distances (ISD) ranging from +20 cm to -58 cm, which is referred to as ISD+20 and ISD-58, respectively. This corresponds to air-gaps between the nozzle exit and the detector surface ranging from 84.8 cm to 6.8 cm. The airgap at the isocenter (ISD0cm) is 64.8 cm. The uncertainties on the spot size and position measurements were 0.2 mm [START_REF] Grevillot | Implementation of dosimetry equipment and phantoms at the medaustron light ion beam therapy facility[END_REF].

Integrated Radial Profiles as function of Depth (IRPDs) were acquired in water using the water phantom MP3-PL (PTW, Freiburg) and the Bragg peak chamber TM 34070 (electrode diameter of 81.6±0.2 mm). The uncertainties of the measured range in water were varying from 0.2 mm to 0.6 mm, depending on the depth and were reproduced with a maximum difference of 0.3 mm [START_REF] Grevillot | Implementation of dosimetry equipment and phantoms at the medaustron light ion beam therapy facility[END_REF].

For beam monitor calibration in reference condition, absorbed dose to water was measured with a plane parallel ionization chamber (ROOS, TM34001, PTW) in a uniform single-layer of 12x12 cm 2 following IAEA TRS-398 adapted for scanned ion beams [START_REF] Palmans | Beam monitor calibration in scanned light ion beams[END_REF]. Dose was measured at 1.4 cm depth in water for energies up to 97.4 MeV and at 2.0 cm depth for larger energies, up to 252.7 MeV. DAP to water was derived by multiplying the absorbed dose to water with the product of the constant spot spacing ∆x ∆y as described in Palmans and Vatnitsky [START_REF] Palmans | Beam monitor calibration in scanned light ion beams[END_REF].

The dose was measured at several positions in clinical treatment plans in order to validate the combined beam model. To minimize volume averaging effects, 24 of the smallest available ICs were used in the MP3-PL water phantom (PinPoint TM 31015, volume 0.03 cm 3 and radius 1.45 mm, PTW, Freiburg) characterized in photons and proton beams [START_REF] Carlino | Characterization of PTW-31015 PinPoint ionization chambers in photon and proton beams[END_REF]. The 24 PinPoints were mounted in a special holder (3D Detector Block) such that they do not shield each other in beam's eye view. The same setup was used to measure clinical plans. Table 1 summarizes the measurements performed and their purposes.

2.C. MC beam modeling

All MC simulations were performed including the full nozzle geometry. The MC beam model was developed consecutively. Starting with a preliminary energy tuning using the nominal beam energy. Subsequently, the pencil beam optical properties were fine-tuned and the energy and energy spread were fine-tuned using the measured IRPDs. Finally, the MC simulations were exploited to retrieved the delivered number of particles per spot from the measured DAP. MC beam parameters were optimized for a subset of different energies and a polynomial interpolation was applied for the intermediate energies. GATE v7.2 [START_REF] Sarrut | A review of the use and potential of the GATE Monte Carlo simulation code for radiation therapy and dosimetry applications[END_REF][START_REF] Jan | GATE V6: a major enhancement of the GATE simulation platform enabling modelling of CT and radiotherapy[END_REF] built with Geant4 v10.02 [START_REF] Allison | Geant4 developments and applications[END_REF][START_REF] Allison | Recent developments in geant4[END_REF] was used to carry out the MC simulations. The QBBC_EMZ physics builder [START_REF] Ivantchenko | Geant4 hadronic physics for space radiation environment[END_REF] was selected (Resch et al 2019). A constant maximum step limiter, tracking cut and range cut were set equal to 0.1 mm was applied [START_REF] Grevillot | Optimization of GEANT4 settings for proton pencil beam scanning simulations using GATE[END_REF]. The mean excitation energy of water I w was set to 78 eV following the recommendations in ICRU Report 90 (ICRU report 90 2014). For the beam optics modeling and the IRPDs simulations, 1E5 and 1E6 primaries were simulated.

Beam optic properties

The FWHMs were extracted from the particles phase space using the Gate PhaseSpaceActor attached to planes corresponding to the measured ISDs. A straight pencil beam was considered in our simulation to initially model the optic properties. The beam optics parameters were iteratively adapted to minimize the differences between measured and simulated FWHM at all ISDs for the horizontal and for the vertical plane independently. The tolerance level was defined as 1 mm in absolute deviation or 10% in relative deviation according to clinical practice at MedAustron.

IRDPs

IRPDs were scored in a water cylindrical geometry mimicking the active volume of the Bragg chamber (radius 40.8 mm). The maximum resolution of the measurements was used as MC scoring resolution (0.1 mm) [START_REF] Grevillot | Implementation of dosimetry equipment and phantoms at the medaustron light ion beam therapy facility[END_REF]. Parameters of the MC simulations being tuned were the mean energy E and the energy spread ∆E. Energy properties were estimated by an iterative procedure to match MC simulations with measurements in terms of physical range (R80) and Bragg peak width at the 80% dose level (BPW80). Other parameters were evaluated: clinical range (R90) and practical range (Rp) defined as the 90% and 10% dose level, respectively, in the distal fall-off of a depth dose profile.

Beam model calibration

One possibility to calibrate a MC beam model in absolute dose is to establish the relationship between the delivered number of particles N and DAP in reference conditions [START_REF] Palmans | Beam monitor calibration in scanned light ion beams[END_REF]:

(1)

The numerator of the equation 1 is the experimentally determined DAP in water at the reference depth z ref over the area A. In general, the term represents the "mean stopping power" per incident proton evaluated in water by any dose engine and it is characterized as follow:

(2)

In the expression 2, is the mass electronic stopping power of water for ion species i and is the fluence differential in energy of that ion species at the reference depth z ref [START_REF] Palmans | Beam monitor calibration in scanned light ion beams[END_REF]. The key point of this approach is to calculate by a MC simulation the term for a known number of protons n, as described by equation 2. Therefore, a single pencil beam was simulated and scored in a cylindrical volume of 7 cm radius, sufficiently large to encompass all the charged particles of the radiation field over the full clinical energy range. Electrons were not tracked by applying a production cut of 5 mm meaning that electrons with a CSDA range lower than 5 mm in the current medium are not transported in the simulation; their energy is considered deposited locally. Depending on the energy, two reference depths were considered: 1.4 cm in the energy range between 62.4 MeV and 97.4 MeV, and 2 cm from 97.4 MeV up to the maximum energy (252.7 MeV) (see section 2.B). Mean stopping power values computed with GATE were used during medical commissioning of the MedAustron beam delivery system in order to calibrate the beam monitors in number of particles per MU (N/MU). Hence, measured DAP and simulated DAP should be in perfect agreement.

2.D. Beam model validation

Regular shaped targets in water

The beam model was validated by means of 3D dose distributions in water for four regular shaped dose cubes of different dimensions and depths in water, following a similar methodology as for the commissioning of the Treatment Planning System (TPS) [START_REF] Carlino | Commissioning of pencil beam and Monte Carlo dose engines for non-isocentric treatments in scanned proton beam therapy[END_REF]. The so-called Box6(0,0,6), Box8(0,0,15), and Box10(0,0,25) were delivered at isocenter where the notation BoxN(0,0,Z) stands for a cubic dose distribution with an N cm side-length and its center placed at the water equivalent depth Z. The Box6(0,0,5), was delivered in non-isocentric conditions (ISD-50cm) using a RaShi. Each treatment plan (TP) was simulated using 1.8E8 primaries. The average unsigned and signed local dose deviations were evaluated:

(3) (4) where K is the number of PinPoint ICs used in the measurements, the dose measured for the chamber i and the calculated dose. For the evaluation of and , the 3D holder was placed in a region were at least 95% of the prescribed dose (1 Gy) was reached (target region). Wherever possible, an additional evaluation was performed for the plateau dose region (proximal region), where a dose level within 60% and 95% of the prescribed dose was reached. IC measurements located at a local dose gradient higher than 0.04 Gy/mm were excluded from the analysis as small positioning and range uncertainties or volume averaging effects would cause high deviations [START_REF] Magro | Dosimetric accuracy of a treatment planning system for actively scanned proton beams and small target volumes: Monte carlo and experimental validation[END_REF][START_REF] Carlino | Implementation of advance methologies in the commissioning of a Light Ion Beam Therapy facility[END_REF], Resch et al 2019). Transverse dose profiles were selected at the center of the SOBP and analysed in terms of Field Size at 50% (FS50) and Lateral Penumbra between 80% and 20% (LP80-20) of the dose profile. Furthermore, distal parameters such as R90, R80, Rp, SOBP modulation defined at 90% dose level (MOD90) and distal penumbra between 80% and 20% dose levels (DP80-20) were analysed.

Clinical TP scenarios

For clinical validation, the patient specific quality assurance (PSQA) of seven patients were selected (33 beams in total). The clinical indications were: two clivus chordoma of the brain, a carcinoma of the palate (para-nasal case), a prostate carcinoma, an adenoid cystic carcinoma (head and neck case), a nasopharyngeal carcinoma and a skull-base meningioma. The Planning Target Volume (PTV) ranged from 21.56 up to 302.16 cm 3 located at depths from approximately 3 cm (proximal) to 22 cm (distal) for the non-RaShi fields. Three RaShi beams were included in the evaluation. Each beam was measured up to three different position with 24 PinPoint ICs. In contrast to the regular shaped targets, the planned dose level varied considerably in the beams of the irregular shaped targets. To account for the varying dose level, a relative dose gradient threshold equal to 2%/mm was applied. A total of 1344 dose measurement points were evaluated. In average, 15 points per beam were passing the dose gradient threshold. All the PSQA plans were measured in non-isocentric conditions at ISD-40cm, following the methodology reported in Carlino 2017 [START_REF] Carlino | Implementation of advance methologies in the commissioning of a Light Ion Beam Therapy facility[END_REF].

Results and Discussion

3.A. Beam modeling accuracy

Beam optics properties

The optical beam model parameters found resulted in a 99% pass rate of the FWHM using the clinical 1mm/10% requirement over the entire 86 cm air-gap range. The largest deviations of 1.8 mm (7.2%) and 1.7 mm (6.9%) were obtained for 62.4MeV at ISD+20cm in the vertical and horizontal planes, respectively. The accuracy of the multiple Coulomb scattering model implemented in Geant4 is mainly driving the agreement in FWHM as the scattering in air is the most predominant effect for this energy and position. However, this combination of energy and position is not clinically relevant in proton therapy. With RaShi the residuals in FWHM were systematically shifted in positive direction (see figure 2). This may be due to the multiple Coulomb Scattering uncertainties in the RaShi (PMMA) [START_REF] Fuchs | Evaluation of GATE/Geant4 multiple Coulomb scattering algorithms for a 160MeV proton beam[END_REF], Resch et al 2018). However, observed deviations were still within clinical tolerances. The integrals of the measured and simulated IRPDs normalized at the first measurement point agreed within 1% (see figure 3(b)) in the entire energy range. Simulated ranges deviated less than 0.2 mm with measurements as can be seen in figure 3(c). The BPW80 deviated by 0.1 mm on average which corresponds to a relative deviation of about 3% (see figure 3(d)). Maximum absolute deviation was 0.3 mm (6%) at 169.3 MeV. Deviations of R80, R90 and Rp for depth dose profiles in water at ISD-50cm with RaShi were less than 0.3 mm while mean integral dose deviations as well as the average peak-to-plateau differences were better than 0.4%.

Beam model calibration

We present in table 2 the differences between our calculated "mean stopping power" (full MC simulation) and ICRU stopping power data (mono-energy beam) for protons with the same residual range at the measurement depth. Differences up to 10.4% were found, with an expected energy dependent behavior. Indeed, ICRU stopping power data are provided for mono-energetic pencil beams, while the "mean stopping power" accounts for the loss of primary protons but includes the nuclear secondaries at the measurement depth. This result illustrates the importance of taking the loss of primary protons and the production of secondary charged particle into account.

Table 2: Summary of the "mean stopping power" calculated according to equation 2. The energy is provided at the nozzle entrance. The last two columns shows the ICRU90 stopping power for mono-energetic protons with the same residual range and the difference of those values from the "mean stopping power" per incident proton.

Our results are in agreement with previous publications (Laitano et al 2000, Gomá et al 2016). The number of incident proton at nozzle entrance is used to normalized the DAP values. In our simulations, less than 1.0% and 0.5% proton loss was found at the phantom entrance and at the nozzle exit, respectively. The contribution of the nozzle to the primary fluence loss is in agreement with the 0.58% reported by Grassberger et al. [START_REF] Grassberger | Characterizing a proton beam scanning system for Monte Carlo dose calculation in patients[END_REF].

Agreement in terms of DAP is shown in figure 4. A mean agreement of -0.1% between simulated and measured DAP was found and maximum deviation was -0.5% at 62.4 MeV. The comparison to measured DAP is an internal consistency test of the correct implementation of the MC calculated values of the "mean stopping power" per incident proton for beam monitor calibration. The residuals shown in figure 4 are the result of the fitting of the calibration curve as we have observed before.

3.B. Beam model validation

Regular shaped targets in water

The transverse and longitudinal dose distribution parameters were in satisfactory agreement with measurements and are summarized in table 3. The FS50 and the LP80-20 agreed on average by 0.4 mm and -0.1 mm, respectively. For the FS50, a maximum deviation of -1.0 mm (-11.4% in relative terms) was observed for the shallower target (95% dose level within 3 and 8 cm range) when the RaShi was used. Deviations for R80 and R90, initially found never higher than 0.2 mm, increase up to 0.5 mm and 0.7 mm, respectively, when a RaShi was used. A spread of 1.3 mm was found in the evaluation of the MOD90 over all the boxes, with a maximum deviation of 2.1 mm for the deepest target (95% dose level within 22 and 30 cm range). The DP80-20 was -0.3 mm in average with a maximum deviation of -0.6 mm.

The averaged signed dose deviated on average over all measured points in the regular shaped targets by +2.7%. This dose difference was found consistent for all boxes and independent of the position (proximal or SOBP), air gap, box size and use of range shifter. The reason for the systematic difference of about 2.7% initially observed in the 3D delivery is currently not understood and is still under investigation. A similar off-set in the same direction was also found during the validation of the beam model implemented in the RayStation TPS (RaySearch Laboratories (RSL), Sweden) as reported by Carlino et al. [START_REF] Carlino | Commissioning of pencil beam and Monte Carlo dose engines for non-isocentric treatments in scanned proton beam therapy[END_REF]. Several factors from the absolute dosimetry to the beam modeling process itself can add up and lead to a systematic offset. One of the main contributions could be imprecisions in the description of the shape of the pre-defined mono-energetic Bragg peaks to fit the experimental IRPDS provided by the clinic. The shape of the pre-defined mono-energetic Bragg peaks depends substantially on the choice of the MC code and physics settings (mainly the hadronic models and the nuclear cross sections [START_REF] Grevillot | Optimization of GEANT4 settings for proton pencil beam scanning simulations using GATE[END_REF], Resch et al 2019) and the deficiency of the fitting process may cause the off-set. Therefore, the absolute dose calibration was rescaled by -2.7%. After rescaling, the largest deviation was -1.2% at the proximal position for the Box10 (see figure 5). The target dose in the regular shaped dose distributions was accurate to within +/-1.0% in isocentric as well as non-isocentric conditions (varying air gaps and the use of the RaShi). A maximum average deviation of -0.2% was found in the SOBP for the non-isocentric Box6 with RaShi. All the analyzed regular shaped targets in water are shown in figure 6 and figure 7. The presented analysis suggest the tendency of an underestimation of the dose in the plateau by about -1% and an overestimation towards the end of the SOBP by about +1%.

Clinical TP scenarios

The dose differences in the clinical treatment plans at ISD-40cm were negligible at low dose gradients, whereas up to 3% deviations were found at high dose gradients. Almost zero deviation and a low standard deviation was found where the average local dose gradient at the measurement points was less than 1 %/mm (see figure 8). Fluctuations and standard deviation increased for higher dose gradient indicating a sensitivity to positioning uncertainty. The three beams with RaShi showed a tendency of underestimating the dose, but the significance of those findings is limited due to the high dose gradient at the measurement points which was 1.2 %/mm on average. The RaShi is required to cover targets shallower than 3 cm, which is difficult to measure with the 24 PinPoints system.

Conclusion

A full nozzle beam model for scanned proton beam delivery has been validated with GATE/Geant4. This beam model starts at nozzle entrance tracking all particles through the nozzle elements, which allows to account for physical interactions within the treatment head, including the passive elements.

The beam model was validated in isocentric and non-isocentric treatment conditions, as routinely performed at MedAustron since 2016. The beam spot size evaluated in FWHM in vertical and horizontal direction fulfilled clinical requirements of 1 mm and 10% at seven positions ranging over 86 cm of air gap. Ranges were reproduced within 0.2 and 0.3 mm (max deviation) without and with range shifter, respectively.

The calibration of the beam monitors is based on a new formalism in dose-area-product and has been clinically implemented for the first time.

Three-dimensional delivery was benchmarked using commissioning data and clinical PSQA measurements. Typical parameters such as field size, modulation, lateral and distal penumbra and range were well reproduced. The dose in the 3D regular shaped targets showed a systematic deviation of +2.7% on average, independently of the treatment depth, field size and modulation. Therefore a final rescaling factor of 0.973 was applied to the beam model and used as such in the rest of the work. The reasons for this deviation are currently unclear and will be the subject of further studies.

Final validation of the beam model using clinical TP scenarios agreed within 1% when the average dose gradient at the measurement positions was less than 1 %/mm.

The detailed characterization of the pencil beam offered by GATE and the implementations of the full nozzle design in our simulations makes our MC beam model accurate, even in complex situations such as the non-isocentric patient positioning in presence of range shifter. The MC beam model will be exploited for future research work and will be inserted as a tool in clinical routine for independent dose calculation in isocentric and non-isocentric conditions.

Figure 1 :

 1 Figure 1: Schematic illustration of the MedAustron (MA) nozzle for the horizontal fixed beam line. It contains a double foil vacuum window, the Independent Termination System (ITS) box, Dose Delivery System (DDS) boxes, and the nozzle exit window. The passive elements, two Ripple Filters (RiFis) and the Range Shifter (RaShi) are also shown. Distances are reported in mm.

Table 1 :

 1 Summary of experimental data used for generation of the MC beam model and for its validation.

Figure 2 :

 2 Figure 2: FWHM agreement in relative terms between simulated and measured values with and without RaShi for the horizontal plane. For clarity purpose, only three representative energies are shown. Similar results were obtained for the vertical plane.

Figure 3 :

 3 Figure 3: In figure 3(a) depth dose profiles simulated in water compared to measured profiles for 20 energies are shown. For this comparison, energy values deposited in the cylindrical volume representing the Bragg peak chamber of 4.08 cm radius was used. In figure 3(b) relative deviations between simulated and measured integral doses are shown. In figure 3(c) absolute deviations between simulated and measured ranges are shown while in figure 3(d) absolute deviations of Bragg peak width values at 80% dose level are reported. Error bars always correspond to the standard deviation. For sake of clarity, in figure 3(c) error bars are reported for only one set of deviations.

Figure 4 :

 4 Figure 4: Relative deviations in terms of DAP per incident proton for a 2D scanned field 12x12 cm 2 and measured at the reference position of 1.4 and 2 cm in depth (depending on the energy).

Table 3 :

 3 Deviations in water for the 3D dose distributions analysed in terms of transverse dose distribution parameters (above) and distal dose distribution parameters (below) in absolute and relative terms. The LP80-20-and the LP80-20+ are referring to the LP80-20 evaluated for the negative and the positive distance, respectively, from the center of the transverse dose profile.

Figure 5 :

 5 Figure 5: Local dose deviations evaluated according to equation 3 and equation 4 after rescaling of -2.7%.

Figure 6 :

 6 Figure 6: Longitudinal (a) and transverse (b) dose profiles for Box6(0,0,6) and longitudinal (c) and transverse (d) dose profiles for the Box6(0,0,5) placed at ISD-50cm with RaShi compared to measurements.

Figure 7 :

 7 Figure 7: Longitudinal (a) and transverse (b) dose profiles for Box8(0,0,15) and longitudinal (c) and transverse (d) dose profiles for the Box10(0,0,25) compared to measurements.

Figure 8 :

 8 Figure 8: Local average dose deviation per beam in clinical cases as a function of the averaged dose gradient. Error bars represents the standard deviation of the relative deviations. The 3% and 5% level are indicated with dashed and solid lines, respectively.

Acknowledgment

The authors are grateful to all the colleagues at the MedAustron medical physics department (Wiener Neustadt, Austria), in particular to Stanislav Vatnitsky for remarkable discussions. The entire CREATIS group at the Léon Bérard hospital (Lyon, France) is fully acknowledged. This work was performed within the framework of the SIRIC LYriCAN Grant INCa_INSERM_DGOS_12563, and the LABEX PRIMES (ANR-11-LABX-345 0063) of Université de Lyon, within the program ANR-11-IDEX-0007. The financial support by the (Austrian) Federal Ministry for Digital, Business and Enterprise and the National Foundation for Research, Technology and Development is gratefully acknowledged.