Jérémy Bourgalais

Nathalie Carrasco

Quentin Changeat

Olivia Venot

Lora Jovanović

Pascal Pernot

Jonathan Tennyson

Katy L Chubb

Sergey N Yurchenko

Giovanna Tinetti

Typeset using

Keywords: planets and satellites, atmospheres -methods, laboratory -ultraviolet, planetary systems -molecular processes

With the upcoming launch of space telescopes dedicated to the study of exoplanets, the Atmospheric Remote-Sensing Infrared Exoplanet Large-survey (ARIEL) and the James Webb Space Telescope (JWST), a new era is opening in the exoplanetary atmospheric explorations. However, especially in relatively cold planets around later-type stars, photochemical hazes and clouds may mask the composition of the lower part of the atmosphere, making it difficult to detect any chemical species in the troposphere or to

INTRODUCTION

Thanks to the Kepler mission (Borucki 2016), many exoplanets of a new kind, bigger than the Earth but smaller than Neptune (1R ⊕ < R p < 4R ⊕ ) have been observed and constitute the majority of the worlds detected outside the solar system (Coughlin et al. 2016). Based on observations and theoretical considerations, many of these low-gravity planets seem to have thick atmospheres which retain light species (H 2 , He) (Jontof-Hutter 2019). However, these small planets show a great diversity of densities, suggesting a diverse range of bulk and atmospheric compositions including rocky planets with thin atmospheres mainly made of heavy molecules such as N 2 , CO, or H 2 O (Rimmer & Rugheimer 2019;[START_REF] Hu | [END_REF]Moses et al. 2013;Elkins-Tanton & Seager 2008;Schaefer et al. 2012).

These atmospheric constituents have a significant impact on the evolution of the planet, including its potential habitability. Many low-gravity planets are detected mainly in orbit around low-mass M-type stars which are the most abundant stars in our galaxy [START_REF] Laughlin | [END_REF]Chabrier 2003). These stars emit strongly in the UltraViolet (UV) and the atmospheric properties as well as the dynamics and fate of their associated exoplanets are thus strongly impacted by stellar radiation (Shkolnik & Barman 2014). Based on our current knowledge of the solar system's bodies, such as Titan (Hörst 2017), high-energy radiation, in particular in the Extreme-UV (EUV) range (10 -120 nm), has a strong effect on the chemical composition of the atmosphere (Chadney et al. 2015). Similar effects are therefore also expected in the atmospheres of exoplanets especially those orbiting M dwarfs which remain active in the EUV much longer than Sun-like G-stars (luminosity factors of 10 -100 times higher than the solar flux (Ribas et al. 2016)) (Lammer et al. 2011). EUV radiation significantly alters the chemical composition of the upper atmosphere of exoplanets during the early phase of their evolution through interaction with the most abundant molecular species, but its long-term effect on the habitability of planets is complex and remains an open question (Luger & Barnes 2015;Wordsworth & Pierrehumbert 2014;[START_REF] Bolmont | [END_REF]Owen & Mohanty 2016). The thermosphere is a complex chemical environment where the main atmospheric components are dissociated and ionized, leading to efficient ion-molecule reactions involving their primary volatile elements (Moses et al. 2011). Many atmospheric chemical models, adapted mostly to hot exoplanets, are available in the literature, some of which considering the effects of metals and simple ions, like H + 3 , in nonequilibrium models [START_REF] Muñoz | [END_REF]Helling & Rimmer 2019;Drummond et al. 2016;Yelle 2004;Lavvas et al. 2014;Koskinen et al. 2013;Rimmer & Helling 2016). However, there is no existing study focused on the characterization of the main positive ions that may form within the thermosphere of warm exoplanets and their expected important effects on the chemistry are yet to be quantified.

In particular, any effect on the distribution of the IR-radiating species which drive thermal and non-thermal atmospheric escape processes will influence the retention of planetary atmospheres.

Determining whether a low-gravity planet around a small active star is able to retain its volatile components is a critical question to understand its nature and evolution. For rocky temperate planets, in particular, this question is pivotal to understand their potential habitability (Luger & Barnes 2015;Bourrier et al. 2017) since UV irradiation may also stimulate the synthesis of prebiotic molecules (Rimmer et al. 2018).

We are at the dawn of a new era in space exploration with the upcoming launch of new telescopes like the Atmospheric Remote-Sensing Infrared Exoplanet Large-survey (ARIEL) [START_REF] Tinetti | [END_REF] and the James Webb Space Telescope (JWST) (Gardner et al. 2006). They will provide a unique insight into the thermal structure and chemical composition of the atmosphere of exoplanets, due to their spectral range and array of instruments that will be used for transmission, emission and phase-curve spectroscopy. To make the most of these future observations, it is important to learn as much as possible about their atmospheric chemistry. Among the different approaches, experimental laboratory simulation is a powerful tool for understanding the chemical evolution of an exoplanet atmosphere and also for guiding observations and the development of atmospheric models. However, so far in the literature, previous laboratory experiments have focused on solid particles (production rates, size, chemical composition) that may constitute photochemical hazes within the cold and warm atmospheres of exoplanets (Hörst et al. 2018;He et al. 2018a,b;Wolters et al. 2019;Moran et al. 2020;Vuitton et al. 2019a).

In order to complete this crucial phase of experimental simulation, it is also necessary to study the chemical evolution of the gas-phase in the primary upper atmosphere of exoplanets following UV irradiation of host stars and to identify the abundant species that may form. While He et al. (2019) studied the neutral species involved in the gas-phase chemistry of simulated exoplanet atmospheres, the implication of the positive ions has never been studied. Photoproducts drive the formation of hazes lower in the atmosphere but are also the starting point for the formation of prebiotic molecules.

To be able to detect biosignatures that will be critical to the search for life, this work paves the first milestones with innovative photochemical-driven molecular growth laboratory experiments and shows the direct relationship with space observations.

The approach in this article aims at providing information on chemical species, in particular, charged species that are expected in the thermosphere of rocky and gas-rich exoplanets and could detectable with future ARIEL and JWST observations. To achieve this goal and to show how the combination of laboratory simulation and numerical modeling constrains the inputs in the modeling of the observations, this work has been divided into three steps:

-In a first step, we have experimentally simulated the effect of UV radiation on simplified gas phase samples (H 2 -CO-N 2 -H 2 O) whose variation of the hydrogen mixing ratio allows us to mimic cases of super-Earths and mini-Neptunes. This first experimental study allowed the measurements of both neutral and cations and to identify the main positive ions formed in H 2 -poor and H 2 -rich environments.

-In the second step, we analyzed the main pathways of formation of these ions using a photochemical model. The characterization of the photoproducts and the identification of key reactions allowed us to confirm that these experimental results were transposable to planetary environments and that these species could also be formed in significant abundance in the thermosphere of exoplanets.

-In the third step, we carried out simulations of observations using the resolution of the future JWST and ARIEL space telescopes to determine whether it was possible to detect in a transmission spectrum of an exoplanet two of the major ions observed in laboratory experiments with realistic mixing ratios. Due to the lack of chemical information on the atmospheres of exoplanets, we chose a planet whose composition most closely resembles one of the gaseous mixtures used in the laboratory (mini-Neptune GJ1214b) where CO replaced by methane is the only difference, but whose impact on the formation of ions is nil.

2. METHODS

Photochemical Reactor Coupled to an EUV Photon Source

Simplified super-Earth and mini-Neptune atmospheres are reproduced in the laboratory using a photochemical reactor that has already been used in the past for photochemical experiments simulating complex environments similar to the upper layers of planetary atmospheres (Carrasco et al. 2013;Peng et al. 2014;Bourgalais et al. 2019). The cell is filled with a gas mixture representative of the environments to be reproduced. The species of the gas mixture are then irradiated by an EUV photon source windowless coupled to the photochemical reactor. The energetic particles are produced by a microwave plasma discharge of a rare gas at low pressure providing radiation at different wavelengths depending on the gas used (Tigrine et al. 2016). The chemical species present in the reactor are probed at steady-state using a quadrupole mass spectrometer (Hiden Analytical EQP 200 QMS ) to measure neutral compounds and positive ions with a high sensitivity below parts-per-million (ppm) level (Bourgalais et al. 2019;Dubois et al. 2019).

In this work, we use a surfatron-type discharge with a neon gas flow that allows us to irradiate the gas mixture of the photochemical cell at an overall pressure of 0.9 mbar and at room temperature.

The pressure is adjusted to represent atmospheric pressures of exoplanets but is experimentally constraint in the simulated environments. Reaching the low pressures found in the higher layer of exoplanet environments is not possible due to experimental constraints. The configuration of the device shown in Bourgalais et al. (2019), shows that the inlet of the MS is located very close from the outlet of the VUV lamp in order to limit wall effects and as a consequence the working pressure adjusted to ensure that reactions would take place before detection with the mass spectrometer. The working pressure is barely higher than the pressure (ca. 0.01-0.1 mbar) at which EUV radiation penetrates (Lavvas et al. 2014) and is low enough to avoid termolecular effects as in the actual lowdensity exoplanet thermospheres. Thus, the higher pressure increases the kinetics in the experiments but the results obtained on the chemistry of the ions remain transposable to the atmospheres of exoplanets. The temperature of the thermospheres of exoplanets can reach tens of thousands of Kelvin, whereas our experiments were conducted at room temperature. However, in the lower part of the thermospheres, the deposition of stellar energy can be compensated mainly by the radiative emission of cooling ions and the subsequent ion chemistry, so that the temperature can drop to 1,000 or even a few hundred of Kelvin (Yelle 2004;Koskinen et al. 2007). At these relatively high temperatures, the temperature dependence of the rate coefficients and branching ratios can be very important. Therefore, it should be noted that the temperature dependence of the formation pathways of the main ions highlighted in this work, was studied with a 0D-photochemical model described in the following section. Numerical simulations performed at room temperature and at 1,000 K, show the same major ions for the two extreme cases discussed in this article. As a conclusion, the simulated environments are relevant to study photochemical mechanisms in upper layer of cold and potentiallly warm exoplanetary atmospheres (< 1,000 K) but these extrapolations at high temperatures must be taken with caution as there is little data on rate coefficients and branching ratios in the literature.

A wavelength of 73.6 nm corresponding to a photon energy of 16.8 eV is used with a flux of about 10 14 ph cm -2 s -1 leading to an ionization rate of ca. 10 -6 . This wavelength is chosen to reach the ionization threshold of the chemical species of the gas mixture that is composed of H 2 , CO, and N 2 .

With this wavelength and the pressure in the reactor, we mimic the stellar EUV field of a star in an upper atmospheric layer (France et al. 2016;Youngblood et al. 2016;Loyd et al. 2016). Transitional planets with a size between Earth and Neptune appear to have a diversity of compositions that is not well-explained but everyone agrees that low-density planets are mainly composed of molecular hydrogen but some super-Earths and mini-Neptunes will likely have thick atmospheres that are not H 2 -dominated. Bulk atmospheric composition could have a wide variety of elemental abundance ratios in which the amount of H 2 has profound effects on the evolution of the atmosphere (Pierrehumbert & Gaidos 2011). Our approach is to use gas mixtures containing H 2 along with other expected important molecules in upper layers of cool super-Earth atmospheres. We vary only the relative proportions of H 2 within a large range to test the contribution of H 2 in the formation of photoproducts in the cases of H 2 -poor and H 2 -rich environments which are representative of simplified super-Earth and mini-Neptune upper atmospheres (Rimmer & Rugheimer 2019). We focus on Earth-like planets by adding N 2 to the gas mixtures (Wordsworth & Pierrehumbert 2014) containing H 2 and since the bulk atmospheric composition of such environment is defined and affected largely in term of the C/O ratio, CO is adding as the only major source of carbon and oxygen to keep a ratio of 1 (Rimmer & Rugheimer 2019). Three gas mixtures have been used with different initial fractional abundances of H 2 . The first is a gas mixture at 1% H 2 , the second is an intermediate at 33% H 2 and finally the third is a mixture at 96% H 2 . The mixing ratio of N 2 was arbitrarily chosen equal to that of CO. Finally, water vapor is the major trace species at ppm level in the reactor coming from residual adsorption on the reactor walls. Contrary to a previous study on Titan conducted with the same device in Bourgalais et al. (2019), here the reactor was not heated to minimize the presence of residual water and the pressure in the reactor is higher (0.9 mbar in this work versus 0.01 mbar in Bourgalais et al. (2019). Gao et al. (2015) showed that water vapor can have profound effects on the evolution of the planet's atmosphere. Thus, the presence of water vapor traces in the reactor provides a study of the impact of trace level of water vapor in the atmosphere of exoplanets.

Mass spectra reported in this study are the average of 10 scans obtained at 2 s/amu with a channeltron-type detection, over the 1-50 amu mass range. In all the spectra presented below, signals with intensities lower than 1,000 cps are considered as background noise. The reported mass spectra show the direct signal recorded without background subtraction. The duration of the irradiation is only a few minutes because, as the numerical simulations in this manuscript show, the state of equilibrium is reached in less than 1 ms for the main species because their abundance is constant.

The stability of the ion signal over the 10 mass spectra recorded during the irradiation supports the claim that the steady state is reached very quickly.

0D-Photochemical Model

To interpret the experimental mass spectra, a coupled ion-neutral photochemical model was used to reproduce the chemistry occurring in the reactor. The model has been detailed in the literature (Peng et al. 2014). We provide here a short description, mentioning the main updates since the initial paper.

The Reactor Model

We consider a simplified 1-cell (0D) geometry of the reactor, assuming a uniform spatial distribution of the species. The gas inflow and outflow are taken into account to ensure a constant pressure, and radiative transfer is treated through Beer-Lambert-type photoabsorption, in the assumption of a uniform gas. The chemical reaction are treated as a system of ordinary differential equations (ODEs).

The implicit-explicit Runge-Kutta-Chebyshev (IRKC) method (Shampine et al. 2006) is used to integrate the system of partial differential equations. The photolysis and transport are treated ex-plicitely (complex evaluation of jacobian elements for radiative-transfer makes it impractical to regard the photolysis rate equations as ordinary differential equations). The chemical ODEs are treated implicitely. The IRKC method handles stiff ODE systems (Verwer & Sommeijer 2004).

The fortran Reactor code is publicly available at https://github.com/ppernot/Reactor, with bash scripts for execution of parallel Monte Carlo runs on an OpenStack cloud (CloudReactor). A graphical interface (ReactorUI) is provided for the analysis of the results.

Chemistry

The MC-ChemDB Database -The MC-ChemDB database (MC-ChemDB) has been designed to handle the uncertainty on the rate parameters of chemical reactions (Hébrard et al. 2006(Hébrard et al. , 2009;;Carrasco et al. 2008;Plessis et al. 2012) and of photolysis processes (Gans et al. 2010;Peng et al. 2014).

The principle is based on a server-client architecture. The server generates and hosts random sets of databases for Monte Carlo uncertainty propagation [START_REF] Bipm | Evaluation of measurement data -Supplement 1 to the "Guide to the expression of uncertainty in measurement" -Propagation of distributions using a Monte Carlo method[END_REF]. The client gathers the required number of samples for adapted chemical schemes: depending on the gas mixture in the reactor, only the pertinent subset of the database is used.

The database is split in three modules (neutral reactions, ionic reactions and photolysis) that require different treatments to generate random reaction rate parameters:

-Neutral chemistry (bimolecular thermal reactions and termolecular thermal recombinations) is built from Hébrard et al. (2006Hébrard et al. ( , 2009)); Dobrijevic et al. (2016). Bimolecular rate constants are represented with the Kooij (modified Arrhenius) rate law, and termolecular association with a pressure-dependent Lindeman-Troe rate law (Dobrijevic et al. 2016;Vuitton et al. 2019b). The scheme is based on partial rate constants.

- A global rate analysis is performed to identify key reactions and dominant formation pathways (Hébrard et al. 2009).

Simulations of Observations with ARIEL and JWST

We investigate the feasibility to detect the ions H + 3 and H 3 O + in sub-Neptune type planets using the radiative transfer code TauREx3 (Waldmann et al. 2015b,a). Our approach consists in simulating a high resolution spectrum using TauREx in forward mode. Then we use our instrument simulators ArielRad [START_REF] Mugnai | EPSC-DPS Joint Meeting[END_REF] to bin down the observations and estimate the instrument noise. As we are investigating the theoretical detection biases through retrieval techniques, we do not scatter the observed spectra (Changeat et al. 2019;Feng et al. 2018). Finally, we retrieve the simulated observations using TauREx3 to check whether the ion signals can be statistically recovered.

We base our simulations on the bulk parameters for GJ1214b from Harpsøe et al. (2013): we use a planet radius of 0.254 R j , a planet mass of 0.0197 M j and a star radius of 0.216 R s . In the case of ARIEL, we consider that all the Tier 3 observations for this planet are stacked together (total of 10 transits) [START_REF] Tinetti | Space Telescopes and Instrumentation 2016: Optical, Infrared, and Millimeter Wave[END_REF], while for JWST we combine single observations from NIRISS and NIRSpec (total of 2 transits) to roughly match the wavelength coverage. The atmosphere is assumed to be isothermal at 700 K and filled with H 2 and He. We add the trace gases H 2 O and CH 4 along with the two ions H + 3 and H 3 O + . For our simplified model, the abundances for H 2 O and CH 4 are set to 10 -5 and are constant with altitude. For the ions, we base our input abundances on profiles resembling those from Helling & Rimmer (2019). As it is expected for sub-Neptune type planets, we enhance their abundance findings by approximately one order of magnitude and use the 2-layers model (Changeat et al. 2019) to describe the abundance decrease deeper in the atmosphere. For H + 3 , the top abundance is 10 -5 and the surface abundance is 10 -14 with a layer pressure change at 10 -5 bar. For H 3 O + , the top abundance is also 10 -6 and the surface abundance is 10 -14 . The layer change is done deeper at 10 -3 bar. On top of this, we add a grey cloud cover which is fully opaque below 5×10 -3 bar as expected for this planet (Kreidberg et al. 2014;Kempton et al. 2011;Morley et al. 2013). For the molecular absorptions, we use the most up-to-date molecular line lists from the ExoMol project (Tennyson et al. 2016), HITEMP (Rothman et al. 2010) and HITRAN [START_REF] Gordon | AAS/Division for Planetary Sciences Meeting Abstracts#[END_REF]) and build the cross-section tables at a resolution of 10, 000. We also consider absorptions from Rayleigh scattering and Collision Induced Absorption (only the pairs H 2 -H 2 and H 2 -He). The list of opacities used in this paper is summarised in Table 1.

H 3 O + Log 0.1-10 -12
Temperature Linear 300-1500 K Clouds Pressure Log 10-10 -7 bar

For our retrieval step, the parameter space is explored using the Multinest algorithm [START_REF] Feroz | Classical and Quantum Gravity[END_REF]) with 500 live points and an evidence tolerance of 0.5. The list of retrieved parameters and their corresponding uniform prior bounds are shown in Table 2.

3. RESULTS We are conducting our experiments in an open system, where the gas mixture is introduced and pumped continuously. In the reactor, the residence time is short, making it difficult to detect neutral molecules larger than 40 amu. Neutral species beside initial gas mixture contain only atomic species H, C, N and O resulting from the fragmentation of their parental molecules upon photodissociation, neon which is the non-reactive gas owed continuously in the photochemical reactor and small molecules formed during the first reaction such as CO 2 .

Identification by mass spectrometry of the neutral species formed in the experiments

Thus, the experimental results obtained for neutral species in this work do not provide as much information and the rest of the study will focus on positive ions, about which our knowledge in this type of environment is much more limited. It should also be noted that the rate of production of ionic products relative to neutral products is similar, with molar fractions of the order of ppm relative to neutral species in the initial gas mixtures.

Identification by mass spectrometry of the main ions formed in the experiments

Figure 2. Positive ions measured after irradiation at 73.6 nm of H 2 /CO/N 2 gas mixtures. H 2 mixing ratios was chosen at 1% (orange), 33% (grey) and 96% (blue) at an overall pressure in the reactor of 0.9 mbar.

Twice the standard deviation, 2σ, in 10 measurements of each mass is displayed.

Fig. 2 shows in situ non-normalized gas-phase mass spectra of cationic species in their stationary point in the reactor after 15 min EUV irradiation at an overall pressure of 0.9 mbar. The species are detected at m/z 1, 2, 3, 18, 19, 20, 21, 28, 29, 30, 31 and 32 (named H x ions) does not follow the same trend as the primary ions N + 2 and CO + at mass 28 which are formed by the ionization of N 2 and CO. In the case of H 2 -rich gas mixtures, H x ions, and in particular H 3 O + , dominate while the primary ions N + 2 and CO + occur in significant amounts in the case of gas mixtures where H 2 is not in majority.

3.3.

Unraveling the chemical processes through ion-neutral modelling

Hydrogen-poor atmospheres

Fig. 3 displays the simulated evolution of the mole fractions of the main ionic photo-products after irradiation at 73.6 nm of H 2 /CO/N 2 gas mixture with a H 2 mixing ratio of 1%. The two most abundant predicted ions using the 0D-photochemical model are HCO + and H 3 O + as observed in the experimental mass spectrum of Fig. 2. The main HCO + formation pathways are found to be through the reactions R1, R2, and R3. HOC + is formed in a second exit channel of R1 but is rapidly converted into HCO + through proton transfer (R3). The main loss pathways of HCO + are found to be dissociative recombination with electrons, along with the proton transfer reaction with water (R4 and R5):

CO + + H 2 --→ HCO + /HOC + + H N 2 H + + CO --→ HCO + + N 2 HOC + + CO --→ HCO + + CO (R1) (R2) (R3)
HCO + + e ---→ H + CO HCO + + H 2 O --→ H 3 O + + CO (R4) (R5)
R5 is the only formation pathway for H 3 O + , using trace water vapor in the reactor. H 3 O + is consumed through dissociative recombination leading to water and OH radicals. Modeling of chemistry within the reactor also shows that the formation of N 2 H + is negligible compared to that of HCO + .

Its main production pathways are:

HOC + + N 2 --→ N 2 H + + CO N 2 + + H 2 --→ N 2 H + + H (R6) (R7)
The low abundance of N 2 H + compared to HCO + is due to the efficient reaction R2 showing the propensity of CO to consume N 2 H + to produce HCO + . As observed in the experimental mass spectrum, the H x ions are expected to be negligible through the main reactions (R8, R9, and R10).

CO + + H --→ H + + CO H 2 + hν --→ H + 2 + e - HOC + + H 2 --→ H + 3 + CO (R8) (R9) (R10)
HOC + indeed mainly reacts with N 2 and CO (R3 and R6) to produce N 2 H + and HCO + respectively rather than with H 2 .

Hydrogen-rich atmospheres

Fig. 4 displays the simulated evolution of the mole fractions of the main ionic photo-products after irradiation at 73.6 nm of H 2 /CO/N 2 gas mixture with a H 2 mixing ratio of 96%. The two most abundant predicted ions using the 0D-photochemical model remain HCO + and H 3 O + as observed in the experimental mass spectrum of Fig. 2. However, the main formation pathways for HCO + is still via R3 but with an additional parallel reaction involving H + 3 :

H + 3 + CO --→ HCO + + H 2 (R11)
Independent of the mixing ratio of H 2 , the main destruction pathways of HCO + remain the same in both H 2 environments (R4 and R5) which includes the main production pathway for H 3 O + (R5). Compared to the case of H 2 -poor atmospheres, the mole fraction of N 2 H + increases significantly in the presence of a large abundance of H 2 due to a production pathway more efficient than R6 and R7:

H + 3 + N 2 --→ N 2 H + + H 2 (R12)
In this H 2 -rich environment, H + 3 starts to be the dominant ion due to a more efficient formation pathway:

H 2 + + H 2 --→ H + 3 + H (R13)
which supersedes R10. Once formed H + 3 leads to the formation of the dominant ions HCO + and N 2 H + by R11 and R12 respectively.

Comparison of the two cases

Laboratory simulations show that the same main ions (HCO + , H 3 O + ) are observed in H 2 -poor and rich environments. This implies that we can expect to observe similar contributions of these ions in the thermosphere of super-Earths and mini-Neptunes. However, they are formed through different chemical pathways by H-transfer reactions leading to the observed stable protonated ions. Another remarkable result is that even with trace amounts of water, the formation of the hydronium ion (H 3 O + ) is very effective. These results highlight the need to accurately model ion-neutral chemistry in atmospheric models of exoplanets. In the case of the low H 2 environment, the primary ions N + 2 and CO + formed at 73.6 nm trigger the molecular growth reacting with H 2 . In the specific case of environments rich in H 2 , the triatomic hydrogen ion H + 3 is one of the most abundant ions and is

shown to be very reactive, leading to the formation of heavier ions (N 2 H + and HCO + ).

Simulations of Observations with ARIEL and JWST

Choice of the planet type and the ion targets

The fact that the H + 3 ion is predicted to be abundant and highly reactive leading to the formation of heavier ions in H 2 -rich environments is of great interest in the context of exoplanetary observations considering the upcoming launches of the ARIEL and JWST telescopes. The abundance and distribution of H + 3 , the stable ionic form of H 2 is essential to understand the chemistry and to get information about the thermal structure, the dynamic, and the energy balance of exoplanet atmospheres (Miller et al. 2000). H + 3 has been thoroughly studied in the mid-InfraRed (IR) on Jupiter (Stallard et al. 2001), Saturn (Geballe et al. 1993) and Uranus (Trafton et al. 1993) but detection in the thermosphere of hot Jupiter exoplanets, with a predicted mixing ratio ranging from 10 -6 up to 10 -4 , so far remains illusive (Lenz et al. 2016;Shkolnik et al. 2006). Helling & Rimmer (2019) assume that the formation pathway of H + 3 through the ionization of H 2 (R11), although very efficient, is difficult in hot exoplanet atmospheres due to the thermal decomposition of H 2 which becomes important when temperature reaches 1,000 K (Yelle 2004). Indeed, the very hot thermosphere of these exoplanets (> 8,000 K) combined with very intense stellar irradiations promotes the formation of atomic neutral and H + making the detection of H + 3 difficult (Koskinen et al. 2010). However, H + 3 is an efficient coolant up to 10,000 K, temperature above which thermal dissociation becomes significant and this thermostat effect is important in controlling the atmospheric stability of giant exoplanets (Neale et al. 1996;Miller et al. 2010Miller et al. , 2013)). By lowering the temperature of the exosphere and the related pressure, H + 3 helps to counteract the atmospheric escape of atomic and molecular species. Theory predicts that the cooling potential of H + 3 would allow one to observe H + 3 in giant exoplanets inside 0.2 to 1 AU orbits (Yelle 2004;Miller et al. 2013;Koskinen et al. 2007). H + 3 is able to offset the increased heating due to EUV radiation which allows to produce increased ion densities by photoionization. Thus, to get sufficient H + 3 an atmosphere relatively cold or far-enough from its star to have a low-temperature thermosphere is needed. In order to evaluate the ability of future telescopes to detect these two ions, we simulate a simplified atmosphere for a GJ1214b-type warm Neptune planet as observed by JWST and ARIEL with respectively 2 and 10 stacked transits. The chemical profiles we use for the H + 3 and H 3 O + ions are inspired by the work of Helling & Rimmer (2019); we also include the absorption of H 2 O and CH 4 , which are the main absorbers. For our simplified model, the abundances for water vapor and methane are set to 10 ppm and are constant with altitude. Those two molecules are predicted to be abundant in sub-Neptune planets (Kempton et al. 2011;Venot et al. 2014;[START_REF] Hu | Planetary Exploration and Science: Recent Results and Advances[END_REF][START_REF] Tinetti | [END_REF]Changeat et al. 2019) and have very strong features in the wavelength range covered by Ariel and JWST. As there are still few models of exoplanets to incorporate reactions with ions, the abundances of H + 3 and H 3 O + that can be expected remain unclear. We present the results with realistic mixing ratios inspired by the work of Helling & Rimmer (2019) though. The purpose of these simulations was to determine with the resolution of future telescopes if these mixing ratios would allow the cations to be detected in the atmospheres of exoplanets and through which transitions. While their profiles may variate with altitude, for the purpose of our simplified example we consider that they are vertically mixed, hence constant with altitude (see Fig. 5).

It should be noted that no object exactly similar to our experimental gas mixtures has yet been identified, but the fact that CH 4 replaces CO/N 2 among the secondary species does not change the final result, which is that H 3 O + and H + 3 will be among the major ions formed in a H 2 -dominated environment where there is the presence of oxygen, even if only in trace amounts (Hollenbach et al. 2012;[START_REF] Beuther | Protostars and planets[END_REF]Gerin et al. 2016) In both ARIEL and JWST cases, our retrieval analysis manages to reproduce the observed spectra (see Fig. 6). We can see that H 3 O + is detectable by both instruments (see also posterior distributions Fig. 7). Indeed, H 3 O + presents a wide feature around 2.8 µm, which is easily captured. H + 3 , however, exhibits much weaker features that are largely masked by the clouds. In general, the grey cloud cover, which represents a pessimistic assumption, introduces a large degeneracy with the planet radius and the retrieved abundances (see Fig. 8). This, along with the assumption of constant chemistry in the retrievals, explain the slight difference from the true and retrieved chemistry for H 2 O and CH 4 .

H + 3 is only successfully detected in the JWST simulated spectrum thanks to the presence of small features, only seen at high resolution. In the ARIEL case, an upper limit of 10 -5 on H + 3 abundance can be deduced but the detection remains unclear. In order to detect this ion with ARIEL, higher abundances or a slightly higher number of transits will be necessary. tures with different mixing ratios of H 2 and important exoplanetary trace compounds (CO and N 2 ) were irradiated at a pressure of 0.9 mbar, representative of high altitude in an exoplanetary atmosphere, and a wavelength of 73.6 nm. The photoproduced ions were measured using a quadrupole mass spectrometer and the experimental data were compared to the expected abundances of a 0Dphotochemical model. A sensitivity analysis based on the Monte Carlo method allows us to confirm H + 3 as one of the main ions in environments rich in H 2 . Our retrieval analysis based on these sim- ulations show that in realistic conditions, the next generation of space telescopes should be able to accurately constrain H 3 O + abundances in sub-Neptunes due to the large broadband feature at around 2.8 µm. For H + 3 , the strongest features seem to require high resolution, which can easily be achieved by JWST but would require the stacking of a large number of observations in the case of ARIEL. Even if the presence of other molecules make detections more difficult, spectral retrieval models provide a powerful technique to disentangle the contribution of the different species from the data (Lammer et al. 2019;Lavie et al. 2017;[START_REF] Goyal | [END_REF]Gandhi & Madhusudhan 2017).

Additionally, H + 3 is very reactive and tends to donate a proton to all the molecules it encounters, enriching chemistry through chain reactions to form more complex and diverse molecules. H + 3 will thus be difficult to observe if it reacts very quickly with neighboring molecules such as CO and N 2 , as shown in this work. Given the results of this work, we suggest that it should be fruitful to search upper exoplanetary atmospheres for abundant protonated ions such as H 3 O + produced by reactions with H + 3 . The simultaneous detection of H + 3 and H 3 O + ions would serve to constrain parameters when characterizing exoplanets.

Several future ground-based and space telescopes (JWST, ARIEL, E-ELT, TMT) will allow the characterization of the atmospheres of exoplanets using infrared spectroscopy, in particular to search for biosignatures within the atmospheres of rocky planets. Cool stars such as M dwarfs are preferred targets for the search of rocky temperate planets, as they have a closer habitable zone due to their size and temperature, which favors the detection and characterization of terrestrial exoplanets (Tarter et al. 2007;Scalo et al. 2007). However, the exposure of planetary atmospheres to the strong UV radiation from these active stars may have drastically negative effect on the habitability of planets (Luger & Barnes 2015;Tian 2009;Erkaev et al. 2013;Miguel et al. 2014;Chadney et al. 2016). Highenergy radiation drives the chemical composition, including the formation of IR-radiating species that will influence the evolution of the atmosphere by controlling its thermal structure. Species like H + 3 and H 3 O + with different cooling rates will allow different EUV-response to thermospheric expansion and hinder the escape of atoms and other volatiles that could be required for life. Thus, the detection of these ions can provide information on the physical parameters of atmospheres and the formation of atmospheric constituents. For instance the saturated molecular hydronium ion is known to play an important role in planetary and interstellar oxygen chemistry network in a way that dissociative recombination of hydronium ions with electrons could be a source of target molecules related to habitability like water vapor and O 2 [START_REF] Larsson | Dissociative recombination of molecular ions[END_REF]. Thus H + 3 or H 3 O + may be used to estimate water vapor and O 2 abundances when direct detection is unfeasible.

As a conclusion, this work proposed H + 3 and H 3 O + as potential tracers of the nature and evolution of low-gravity exoplanets. H + 3 formation is negligible in low H 2 environments but important in H 2 -rich gas mixtures, while H 3 O + should remain abundant in both environments if water vapor is present in the atmosphere even at trace level. Thus, non-detection of H + 3 associated with the observation of H 3 O + in an exoplanetary atmosphere would provide additional information for the classification of planets in the transition between super-Earths and mini-Neptunes.

Future studies should focus on the potential detection of the other major ions highlighted in this work, HCO + and N 2 H + , which may be relatively abundant within the atmospheres of exoplanets where CO and N 2 would be present. However, theoretical calculations are required to provide a list of absorption transitions and their related intensities for these species. One strategy to support the theoretical calculations would be to develop experimental devices to measure ion absorption cross sections in IR. Experimental devices must begin to focus on the evolution of thermospheric chemistry as a function of temperature, which requires the development of simulation chambers that can be temperature-controlled between 300 and 1,000 K. Finally, photochemical models of exoplanets should begin to take into account the coupling between the chemistry of neutral species and ions for the simulation of atmospheres of cold and warm exoplanets.

Figure 1 .

 1 Figure 1. Typical in situ gas-phase mass spectrum of neutral species present in the reactor with a H 2 -CO-N 2 gas mixture after EUV irradiation.

Figure 3 .

 3 Figure 3. Simulated evolution of the mole fractions of the main ionic photo-products after irradiation at 73.6 nm of H 2 /CO/N 2 gas mixture with a H 2 mixing ratio of 1%. All species with a mole fraction above 10 -11 are represented.

Figure 4 .

 4 Figure 4. Simulated evolution of the mole fractions of the main ionic photo-products after irradiation at 73.6 nm of H 2 /CO/N 2 gas mixture with a H 2 mixing ratio of 96%. All species with a mole fraction above 10 -11 are represented.

Figure 5 .

 5 Figure 5. Atmospheric profiles used as input for foward simulations of exoplanet spectrum, as well as for fitting data in retrievals.

Figure 6 .

 6 Figure 6. Retrieved spectra and volume mixing ratios for our JWST (top) and ARIEL (bottom) simulations of a GJ1214b-like planet with H 3 O + and H + 3 ions. For the abundances the solid line is the true value, the dashed line the mean retrieved value and the shaded area represents the 1σ retrieved value.

  4. CONCLUSIONLaboratory experiments simulating thermospheric chemistry in simplified atmospheres of super-Earth and mini-Neptunes exoplanets were performed. A photochemical reactor containing gas mix-

Figure 7 .

 7 Figure 7. Posterior distributions for ARIEL (blue) and JWST (green) retrievals of our GJ1214b planet case.

Figure 8 .

 8 Figure 8. Contribution of the different absorbers to the transmission spectra of our best fit retrievals for ARIEL (top) and JWST (bottom).

  Ionic chemistry (ion-molecule reactions and dissociative recombinations) is built from Carrasco et al. the Leiden database (http://www.strw.leidenuniv.nl/˜ewine/photo)(Heays et al. 2017). If BR are not provided by Leiden, they are extracted from SWRI (https://phidrates.space. swri.edu) (Huebner & Mukherjee 2015), except for CH 4 , where the representation by Peng et al. -The model simulates chemistry under conditions similar to those of the experiments (Bourgalais et al. 2019). It is run for a time long enough to reach a stationary state (1 s) and the stationary mole fractions of the products are compared with the experimental data. Complementary Monte Carlo simulations are performed to evaluate the uncertainty on the model predictions (500 runs).

	(2014) is used. Cross-sections are provided with a 1 nm and 0.1 nm resolutions. At this stage,
	temperature-dependence of the cross-sections and BRs is ignored.
	For a description of the parameters uncertainty representations, see Peng et al. (2014). The database
	was originally oriented on the chemistry of H, C, and N-bearing species, and it has recently been
	upgraded to account for oxygenated species from Vuitton et al. (2019b) and the KIDA database
	(Wakelam et al. 2012).
	(2008); Plessis et al. (2012). Rate laws for ion-neutral reactions are mostly based on the Langevin law,
	with temperature-dependent laws of type ionpol (Wakelam et al. 2012) for molecules with non-zero

dipole moments. Dissociative recombinations are based on a modified Arrhenius rate law

(Plessis et al. 2012)

. The global reaction rates and branching ratios are treated separately.

-Photo-processes (photo-dissociation and photo-ionization), cross-sections and branching ratios (BR) are extracted from

The photochemical scheme -According to the following procedure a consistent set of reactions is iteratively generated from the irradiated initial mixture (H 2 /N 2 /CO/H 2 O):

1. Select all the reactions involving the list of species; 2. Update the list of species with the generated products; 3. Iterate to (1) until no new species is produced. The resulting chemical model contains 54 photo-processes (photolysis of N 2 , H 2 , CH 4 , C 2 H 2 , C 2 H 4 , C 2 H 6 , HCN, NH 3 , CO & H 2 O), 903 neutral reactions (811 bimolecular and 92 termolecular), 1941 ion processes (1314 bimolecular and 687 dissociative recombinations), involving 177 neutral species, 190 positive ions with masses up to 130 and electrons.

Simulations

Table 1 .

 1 List of IR opacities used in this work.

	Species	References
	H 2 -H 2	(Abel et al. 2011; Fletcher et al. 2018)
	H 2 -He	(Abel et al. 2012)
	H 2 O	(Barton et al. 2017; Polyansky et al. 2018)
	CH 4	(Hill et al. 2013; Yurchenko et al. 2014)
	H + 3	(Mizus et al. 2017)
	H 3 O +	Yurchenko et al. (to be published)

Table 2 .

 2 List of retrieved parameters and their corresponding uniform prior bounds.

	Parameter	Mode	Bounds
	Radius	Linear 0.1-0.4 R j
	H 2 O	Log	0.1-10 -12
	CH 4	Log	0.1-10 -12
	H + 3	Log	0.1-10 -12

  above the limit of the background noise. m/z 1, 2 and 3 are attributed to H + , H + 2 , and H + 3 respectively. Due to the high sensitivity of the mass spectrometer, H 2 O + and its protonated form H 3 O + are detected along with first H 3 O + isotopologues (molecules that only differ in their isotopic composition) at m/z 18, 19, 20 and 21. Finally, m/z 28, 29, 30, 31 and 32 are attributed to N + 2 and CO + , HCO + and/or N 2 H + along with the first HCO + isotopologues. The first notable result is that the most abundant ions at m/z 19 and m/z 29 remain the same with all gas mixtures. The formation of ions H + , H + 2 , and H + 3 ions

  Simulation of the observation of ions in a GJ1214b-like planet by ARIEL and JWST

	3.4.2.
	The results of the laboratory experiments in
	this work suggest that colder atmospheres dominated by H 2 would be more conducive to H 2 stability
	leading to a higher H + 3 abundance. Warm Neptunes (ca. 400 to 800 K) would therefore potentially
	be better candidates than hot Jupiters for the detection of H + 3 .
	Like H + 3 , long lifetimes associated with certain excited states can lead to population trapping and
	unexpected state distributions in non-thermalized environments. A recent theoretical study by Mel-
	nikov et al. (2016) calculated the stability of the ro-vibrational states of H 3 O + , the lifetimes of
	individual states and the overall cooling rates. H 3 O + is present in abundance in diffuse and dense
	regions of the interstellar medium, such as comets and molecular clouds (Goicoechea & Cernicharo
	2001; Barber et al. 2009). So far, no attempt has been made to observe H 3 O + in the thermosphere
	of exoplanets, although it may, like H + 3 , act as a cooling agent and bring some constraints on the
	physical parameters of the atmosphere.

  . Efficient series of proton transfer reactions from H +

	3
	or other abundant ions, depending on the secondary species in the environments, will easily lead to
	the formation of H 3 O + (Indriolo & McCall 2013; Van Dishoeck et al. 2013), making the detection of
	H 3 O + and H + 3 in an object like GJ1214b quite legitimate.

ACKNOWLEDGMENTS

This research was supported by the European Research Council Starting Grant PRIMCHEM 636829 to N.C.. O.V. thanks the CNRS/INSU Programme National de Planétologie (PNP) and the CNES for funding support. We acknowledge the help and the fruitful discussions with Dr Billy Edwards, who provided access to his JWST noise simulator ExoWebb.


Jean-François F Ritz
email: jean-francois.ritz@gm.univ-montp2.fr

Regis Braucher
email: braucher@cerege.fr

Didier Bourlès
email: bourlès@cerege.fr

Bertrand Delouis
email: delouis@geoazur.unice.fr

Carlos Marquardt
email: cmarquardt@ing.puc.ci

Alain Lavenu
email: alavenu2@orange.fr

Hervé Philip
email: philip@gm.univ-montp2.fr

Luc Ortlieb
email: luc.ortlieb@bondy.ird.fr

Régis Braucher

Slip rate of trench-parallel normal faulting along the Mejillones Fault (Atacama Fault System): Relationships with the northern Chile subduction and implications for seismic hazards

des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

INTRODUCTION

The northern part of the arid Chile Andean forearc is known to have undergone major interplate underthrusting seismic events over the past 150 years (e.g. [START_REF] Comte | Stress distribution and geometry of the subducting Nazca plate in northern Chile using teleseismically recorded earthquakes[END_REF][START_REF] Madariaga | Central Chile finally breaks[END_REF] (Figure 1A). It is described as a strongly coupled plate boundary between oceanic and continental crusts, with a strong interseismic coupling (i.e. [START_REF] Madariaga | Central Chile finally breaks[END_REF][START_REF] Chlieh | Interseismic coupling and seismic potential along the Central Andes subduction zone[END_REF]. Since the 1877-M8.6 mega earthquake, only four events of lesser magnitude occurred, in 1995 (Mw 8.1 Antofagasta), in 2007 (Mw 7.7, Tocopilla), and in 2014 (Mw 7.6 and Mw 8.1, Iquique), and a large seismic moment deficit remains in this region (e.g. [START_REF] Ruegg | The Mw=8.1 Antofagasta (North Chile) earthquake of July 30, 1995: First results from teleseismic and geodetic data[END_REF][START_REF] Delouis | Recent crustal deformation in the Antofagasta region (northern Chile) and the subduction process[END_REF][START_REF] Delouis | The Mw 7.7 Tocopilla earthquake of 14 November 2007 at the southern edge of the northern Chile seismic Gap: Rupture in the deep part of the coupled plate interface[END_REF][START_REF] Peyrat | Kinematic rupture process of the 2007 Tocopilla earthquake and its main aftershocks from teleseismic and strong-motion data[END_REF]Béjar-Pizzaro et al., 2010[START_REF] Ruiz | Intense foreshocks and a slow slip event preceded the 2014 Iquique Mw 8.1 earthquake[END_REF] (Figure 1A).

The whole northern Chile region is also characterized by a set of N-S sub-parallel faults that stand out clearly in the morphology (Figure 1B). Called the Atacama Fault System (AFS) [START_REF] Arabasz | Geological and geophysical studies of the Atacama Fault Zone in Northern Chile[END_REF], these eastward dipping faults, inherited from the Mesozoic tectonics (e.g., [START_REF] Brown | Displacement history of the Atacama fault system 25°00'S-27° 00'S, northern Chile[END_REF][START_REF] Scheuber | Tectonics of the Jurassic-Early Cretaceous magmatic arc of the north Chilean Coastal Cordillera (22º-26ºS): a story of crustal deformation along a convergent plate boundary[END_REF][START_REF] Cembrano | Fault zone development and strain partitioning in an extensional strike-slip duplex: A case study from the Mesozoic Atacama fault system, Northern Chile[END_REF] show evidences of recent activity (e.g., [START_REF] Okada | On the neotectonics of the Atacama fault zone region. Preliminary notes on late Cenozoic faulting and geomorphic developement of Coast Range of northern Chile[END_REF][START_REF] Armijo | Active faulting in northern Chile: ramp stacking and lateral decoupling along a subduction plate boundary?[END_REF]Gonzales and Carrizo, 2003;Gonzales andCarrizo, 2003, 2006;[START_REF] Vargas | U-series dating of co-seismic gypsum and submarine paleoseismology of active faults in Northern Chile (23°S)[END_REF][START_REF] Cortés | Paleoseismology of the Mejillones Fault, northern Chile: Insights from cosmogenic 10Be and optically stimulated luminescence determinations[END_REF]. [START_REF] Armijo | Active faulting in northern Chile: ramp stacking and lateral decoupling along a subduction plate boundary?[END_REF] started describing the Atacama fault as a left-lateral strike-slip fault, and invoked a local consequence of an E-W extension affecting the coast line (i.e. Coastal Scarp).. But, the following works showed that the Atacama fault system, including the Mejillones peninsula, corresponded to a system of normal faults [START_REF] Delouis | Recent crustal deformation in the Antofagasta region (northern Chile) and the subduction process[END_REF]Gonzales and Carrizo, 2003;Gonzales et al., 2003;2006;[START_REF] Marquardt | Déformations néogènes le long de la côte nord du Chili (23°-27° s), avantarc des Andes Centrales[END_REF][START_REF] Vargas | U-series dating of co-seismic gypsum and submarine paleoseismology of active faults in Northern Chile (23°S)[END_REF][START_REF] Cortés | Paleoseismology of the Mejillones Fault, northern Chile: Insights from cosmogenic 10Be and optically stimulated luminescence determinations[END_REF][START_REF] Cortés-Aranda | Normal upper plate fault reactivation in northern Chile and the subduction earthquake cycle: From geological observations and static Coulomb Failure Stress (CFS) change[END_REF].

Observing fresh normal displacement along the Paposo fault that did not exist before the 1995 Antofagasta earthquake (see location "Fig. 8A" in Figure 1C), [START_REF] Delouis | Recent crustal deformation in the Antofagasta region (northern Chile) and the subduction process[END_REF] proposed that the Atacama fault system was strongly coupled with the seismic cycle along the 4 subduction zone, with normal faulting events contemporaneous with strong co-seismic subduction events. The same kind of observations were made during the M w 8.8, 2010 Maule earthquake near Santiago, with the occurrence of large aftershocks nucleated on normal faults in the upper plate [START_REF] Allmendinger | Permanent fore-arc extension and seismic segmentation: Insights from the 2010 Maule earthquake, Chile[END_REF]. This interpretation is consistent with the geodetic data recorded during the 1995 Antofagasta earthquake that evidence a rapid increase of the horizontal coseismic displacement towards the trench when moving from points located inside the forearc to points located near the coast [START_REF] Ruegg | The Mw=8.1 Antofagasta (North Chile) earthquake of July 30, 1995: First results from teleseismic and geodetic data[END_REF][START_REF] Klotz | GPS-derived 23 deformation of the Central Andes including the[END_REF]). On the contrary, the horizontal displacement of the forearc during the interseismic period is characterized by vectors directed landward that decrease progressively from west to east [START_REF] Chlieh | Crustal deformation and fault slip during the seismic cycle in the North Chile subduction zone, from GPS and InSAR observations[END_REF]. Analysing Insar data, Loveless and Prichard (2008) obtained ambiguous results with some images showing offset consistent with coseismic faulting on the Paposo segment, while others lack such signal.

For Gonzales and Carrizzo (2003) and Gonzales et al. (2003Gonzales et al. ( , 2006) ) the normal faulting process along the AFS is associated with a long-term extensional strain due to the buckling of the margin with no direct co-seismic coupling with the subduction events.

In [START_REF] Delouis | Recent crustal deformation in the Antofagasta region (northern Chile) and the subduction process[END_REF]'s interpretation, large, shallow underthrusting events at the plate interface induce an increment of ~E-W extension within the overlying upper crust, triggering ruptures along the Atacama normal fault system. According to [START_REF] Delouis | Recent crustal deformation in the Antofagasta region (northern Chile) and the subduction process[END_REF], this increment of extension is associated with an increase of the Coulomb stress for normal failure on the AFS. A comparison of the distributions of coseismic and interseismic horizontal displacements on the overriding plate [START_REF] Ruegg | The Mw=8.1 Antofagasta (North Chile) earthquake of July 30, 1995: First results from teleseismic and geodetic data[END_REF][START_REF] Klotz | GPS-derived 23 deformation of the Central Andes including the[END_REF][START_REF] Chlieh | Crustal deformation and fault slip during the seismic cycle in the North Chile subduction zone, from GPS and InSAR observations[END_REF] indicates that coseismic extension is concentrated on the coastal area whereas interseismic contraction is more broadly distributed in the whole forearc area and beyond towards the East.

5

This may explain why extension across the AFS is not cancelled at the end of a seismic cycle, and may slowly increase from one cycle to the other. [START_REF] Delouis | Recent crustal deformation in the Antofagasta region (northern Chile) and the subduction process[END_REF]'s interpretation has three implications: 1) each complete seismic cycle at the subduction interface produces a net increase of extension in the coastal area, explaining the long-lasting activity and persistent morphological signature of the normal faults;

2) during the interseismic period, extension in the coastal area is only partially reduced, which implies a partial reduction of the Coulomb Failure Stress (CFS) on the AFS. This is supported by the absence of present-day seismicity along the AFS; and 3) movements along the AFS are expected to coincide with the abrupt increase of the CFS on the normal faults, i.e. when large subduction events occur. As confirmed by [START_REF] Allmendinger | Permanent fore-arc extension and seismic segmentation: Insights from the 2010 Maule earthquake, Chile[END_REF], the static deformation field produced by a great subduction earthquake is an effective mechanism for generating permanent extension above the seismogenic zone, reactivating suitably oriented, long-lived normal faults. In the case of the Maule earthquake, [START_REF] Allmendinger | Permanent fore-arc extension and seismic segmentation: Insights from the 2010 Maule earthquake, Chile[END_REF] calculated that the earthquake produced a radial pattern of static extension and deviatoric tension, with a CFS increment larger than 1.5 MPa, along the Coastal Cordillera enclosing the rupture area. [START_REF] Cortés-Aranda | Normal upper plate fault reactivation in northern Chile and the subduction earthquake cycle: From geological observations and static Coulomb Failure Stress (CFS) change[END_REF] explored further the idea of a CFS link between the subduction, the Mejillones fault system (MFS) within the Mejillones Peninsula, and the rest of the Atacama fault system further east, with a more thorough modelling of the CFS associated to both the co-seismic and inter-seismic periods. In particular, they show that the lower seismic coupling below the Mejillones peninsula [START_REF] Métois | Revisiting the North Chile seismic gap segmentation using GPS derived interseismic coupling[END_REF] may even be associated with a positive increase of CFS on the MFS in the interseismic stage. They hence confirm that the MFS activation can be explained in terms of CFS changes generated by slip on the plate interface.

In the present study, we do not discuss the role and origin of the west-dipping coastal scarp, which could act differently than the AFS [START_REF] Armijo | Active faulting in northern Chile: ramp stacking and lateral decoupling along a subduction plate boundary?[END_REF]Contreras-Reyes et al., 6 2012). We consider, as in [START_REF] Delouis | Recent crustal deformation in the Antofagasta region (northern Chile) and the subduction process[END_REF] and [START_REF] Cortés-Aranda | Normal upper plate fault reactivation in northern Chile and the subduction earthquake cycle: From geological observations and static Coulomb Failure Stress (CFS) change[END_REF], that the MFS is part of the AFS and respond to the same kind of CFS mechanical interaction with the subduction plate interface.

If the mechanism driving this trench parallel normal faulting is associated with the seismic cycle within the northern Chile subduction zone, analyzing geomorphology along these faults may provide information on the seismic potential associated with the subduction zone. In this paper, we present a study carried out in February 1999 and April 2002 (Marquardt, 2005 PhD thesis), aiming at quantifying precisely the long-term average slip rate along the Mejillones fault, one of the most prominent features within the AFS (Figures1B and 2). The precise quantification of the offsets affecting the alluvial surfaces at three sites along the fault, and the determination of their ages taking into account the recent developments and knowledge enhancements concerning 10 Be cosmic-ray exposure dating, allow us to determine the slip rate along the fault for the past 60 ka. Combining our results with new observations about the size of potential coseismic displacements along the AFS, we bring new insights about the seismic hazard in the region.

TECTONIC SETTING

The AFS extends along the northern Chilean coast over more than 1,000 km. It corresponds to an eastward dipping fault system with a mean N-S direction. Between 23° and 23.5°S, at the latitude of the Mejillones peninsula, the AFS distributes over a 40 km wide fault zone (Fig. 1B) with fault directions evolving from NNE-SSW (i.e. Salar del Carmen, Cerro Gordo, Cerro Moreno and northern Cerro Fortuna faults) to NNW-SSE (i.e. Caleta Herradura fault) with N-S trending faults or fault sections in between (i.e. Mejillones and southern Cerro Fortuna faults). Bathymetry and seismic imagery show that this set of normal faults extends offshore and forms a seaward verging domino-type structural domain [START_REF] Von-Huene | Subduction erosion and basal friction along the sedimentstarved convergent margin off Antofagasta, Chile[END_REF].

7

On land Neogene to Quaternary reactivated segments of the fault system are recognized north of 27ºS, where fault scarps and open cracks are found within the youngest deposits [START_REF] Armijo | Active faulting in northern Chile: ramp stacking and lateral decoupling along a subduction plate boundary?[END_REF][START_REF] Riquelme | A geomorphological approach to determining the Neogene to Recent tectonic deformation in the Coastal Cordillera of northern Chile (Atacama)[END_REF][START_REF] González | Segmentación, cinemática y cronología relativa de la deformación tardía de la Falla Salar del Carmen[END_REF][START_REF] Vargas | U-series dating of co-seismic gypsum and submarine paleoseismology of active faults in Northern Chile (23°S)[END_REF][START_REF] Cortés | Paleoseismology of the Mejillones Fault, northern Chile: Insights from cosmogenic 10Be and optically stimulated luminescence determinations[END_REF]. Age constraints on neotectonic features are scarce because of the hyperarid climate [START_REF] Vargas | Aluviones históricos en Antofagasta y su relación con eventos El Niño/Oscilación del Sur[END_REF], which precludes the use of radiocarbon dating. K-Ar dating of biotites from reworked ashes interbedded with alluvial deposits cut by the AFS leads a maximum age for the onset of neotectonic activity between 2.9 to 6.1 Ma (González et al., 2003b). Morphological modelling based on a diffusional model of erosion allowed [START_REF] González | Segmentación, cinemática y cronología relativa de la deformación tardía de la Falla Salar del Carmen[END_REF] estimating a post-400 ka activity along the Salar del Carmen fault, the easternmost fault of the AFS (Fig. 1B). During the past 500 years, the reactivation of the AFS has been observed only once, during the 1995, Mw8.0 Antofagasta earthquake. Surface breaks were found along several segments of the southern Salar del Carmen fault, sometimes associated with a slight normal offset [START_REF] Delouis | Recent crustal deformation in the Antofagasta region (northern Chile) and the subduction process[END_REF][START_REF] González | Segmentación, cinemática y cronología relativa de la deformación tardía de la Falla Salar del Carmen[END_REF].

Within the Atacama Fault System, the Mejillones Peninsula appears as a prominent geomorphic feature. Interrupting the N-S linear trend of the coastline, the peninsula displays uplifted marine terraces, beach ridges and alluvial fans, which are crosscut by several normal faults and large-scale tension gashes (Figures 1B,1C). One of the main fault scarps is observed at the north-eastern termination of the peninsula where the Mejillones fault cuts through Quaternary alluvial fans. There, cumulative offsets of several tens of meters are observed (Figure 2). The low denudation rates due to the regional pervasive hyper-aridity, and the lithology makes this site suitable to apply the cosmic-ray exposure age dating method using in situ-produced 10 Be (e.g. [START_REF] Ritz | Slip rates along active faults estimated with cosmic-ray-exposure dates : Application to the Bogd fault, Gobi-Altaï, Mongolia[END_REF][START_REF] Bierman | Cosmogenic ages for earthquakes recurrence intervals and debris fl ow fan deposition[END_REF][START_REF] Marquardt | Déformations néogènes le long de la côte nord du Chili (23°-27° s), avantarc des Andes Centrales[END_REF].

Using both 10 Be and OSL dating techniques and detailed analysis of paleoseismic trenches, [START_REF] Cortés | Paleoseismology of the Mejillones Fault, northern Chile: Insights from cosmogenic 10Be and optically stimulated luminescence determinations[END_REF] determined a steady slip rate of 0.61± 0.26 mm/yr along the 8 Mejillones fault between 35 and 3 ka. They also estimated that the fault was capable of producing Mw ~7 earthquakes with mean recurrence intervals of 5 ± 3.5 ka. This is consistent with Vargas et al., (2011)'s paleoseismic investigations, which used the U-series disequilibrium method to constrain the timing of earthquake-induced fissures in gypsum precipitated in extensional fissures along the fault. Both studies agree to date the last significant earthquake around ~3, and concluded that the cycle of the large earthquakes along the Mejillones Fault was not in phase with the cycle of large (Mw > 8.5) subduction earthquakes.

MORPHOTECTONIC STUDY

The 20 km long Mejillones normal fault bounds to the east an uplifted faulted block (Morro Mejillones) which displays a series of Pliocene and Pleistocene wave-cut platforms abraded in basement rocks (Figure 3). Easwards the Mejillones normal fault, in the hanging-wall, the Pampa Mejillones is a sedimentary basin filled by Neogene marine deposits and covered by E-W trending Quaternary sequences of beach-ridges formed over the last 400 ka [START_REF] Ortlieb | A warm interglacial episode during isotope stage 11 in Northern Chile[END_REF] (Figure 3). Near the fault, the marine deposits are interbedded with thick alluvial sequences deposited at the outlets of the drainage basins incising the eastern flanks of the Morro Mejillones.

The obvious cumulative fault scarp cutting through the latest alluvial deposits demonstrates the recent activity on the Mejillones fault (Figure 2).

We carried out a morphological study at three sites located at the outlet of three different drainage basins (Figure 3). At all three sites (O1, O2 and O3), two main alluvial surfaces are displaced vertically above the youngest fan surface S0 incised by intermittent stream channels (Figure 4A). These features reflect the interplay between climate-forced alluviation cycles and co-seismic movements along the fault leading to the abandonment of the surfaces. S2, the most obvious is the landscape, corresponds to the oldest preserved fan surface. S1, inset in S2, 9 corresponds to a younger surface. The youngest surface S0 does not show any offset. South of site O3, in the direct extension of the easternmost normal fault strand (see dashed line in Figure 3B), we noticed that S0 was affected by a centimetric width fissure (Figure 4B). This feature could be related to the 1995 Antofagasta earthquake, like the surfaces ruptures observed by [START_REF] Delouis | Recent crustal deformation in the Antofagasta region (northern Chile) and the subduction process[END_REF] along the Salar del Carmen fault after this event.

The surfaces S1 and S2 are present on both sides of the Mejillones fault. In the hangingwall, eastwards of the Mejillones fault, S1 is slightly inset in S2 (Figure 4C) showing that there is a general uplift over the entire area as also shown by the uplifted beach-ridges within the Pampa Mejillones [START_REF] Ortlieb | Late Quaternary Coastal Changes in northern Chile: Guidebook for a fieldtrip[END_REF]. Fan material contains both clast-supported stream flow deposits and poorly consolidated coarse-grained debris flow deposits found at the top of alluvial sequences.

These debris flows contain angular boulders embedded in a sand-gravel matrix with dimensions scaling with the size of the feeding mountain streams. The largest boulders have diameters up to 1 m, with the exposed faces of the boulders displaying significant desert varnish.

Cumulative offsets and cosmic-ray exposure ages

Surfaces S2 and S1 found on both sides of the fault scarp enable accurate estimates of the tectonic cumulative offsets. A differential GPS survey, with 2 Trimble 5700 GPS stations was carried out within three outlets O1, O2 and O3. Digital Elevation Models (DEM) were generated from the interpolation of the topographic data using Surfer software and a Kriging gridding method, from which fault offsets were measured (Figure 5, see also Figures Sup1, Sup2 andSup3 in Supplementary Material). The vertical offsets affecting S1 surface at the three studied outlets O1, O2 and O3, are 5.2 ± 0.4 m, 4.2 ± 0.4 m and 5.8 ± 0.2 m, respectively (these values correspond to the mean value between the maximum and minimum vertical offsets affecting the surface (see Figure 5); the uncertainties correspond to the differences between the mean value 10 and the extreme ones). The vertical offsets affecting S2 surface at the three studied outlets are 12.7 ± 0.7 m (Site O1), 11.4 ± 0.8 m and 14.7 ± 0.4 m (Site O2), 11.8 ± 0.5 m and 13.2 ± 0.8 m (Site O3). These measurements yield average cumulative displacements of 5.1  0.8 m for surface S1 and 12.8  1.3 m for surface S2.

The minimum cosmic ray exposure ages of boulders sampled from the surfaces of offset S1 and S2 alluvial fans were determined measuring in situ-produced cosmogenic 10 Be concentrations in quartz (e.g., [START_REF] Brown | Examination of surface exposure ages of Antarctic moraines using in situ-produced 10Be and 26Al[END_REF][START_REF] Gosse | Terrestrial in situ cosmogenic nuclides: theory and applications[END_REF]. Our sampling strategy aimed to minimize the effects of exposure prior to deposition and of denudation following deposition. Debris flows related to surfaces S1 and S2 correspond to alluvial pulses reworking massive quantities of slope material, which suggests rapid transport from the upstream drainage basins to their outlets. However, it has been shown that pre-exposure due to complex exhumation-transport history may occur even in the case of major alluviation (e.g. [START_REF] Ritz | Late Pleistocene to Holocene slip rates for the Gurvan Bulag thrust fault (Gobi-Altay, Mongolia) estimated with 10Be dates[END_REF]2006;[START_REF] Vassallo | Dating faulted alluvial fans with cosmogenic 10Be in the Gurvan Bogd mountain range (Gobi-Altay, Mongolia) : climatic and tectonic implications[END_REF]2011;[START_REF] Rizza | Slip rate and slip magnitudes of past earthquakes along the Bogd left-lateral strike-slip fault (Mongolia)[END_REF]. Since 10 Be production decreases exponentially with depth, any deviation from this theoretical distribution indicates pre-deposition inheritance or post-deposition displacement. When possible, two or three samples (top-bottom or top-middle-bottom) have been collected from granite boulders embedded in S1 and S2 surfaces to identify any such previous complex exposure histories (Figure 6).

Within both S1 and S2 surfaces, five boulders distributed on both sides of the fault scarp were sampled (Table 1). To minimize the effects of denudation after deposition, the biggest boulders showing desert varnished and standing above the surface, but still encased in it, were sampled (Figure 6). Surface S2 was sampled at Site O3 (Figures 3B,5). Five boulders were sampled; four in the footwall and one in the hanging wall. Surface S1 was sampled at Site O2 11 (Figures 3B,5). Three samples were we also collected along a 2 m deep soil pit excavated in surface S1 to determine the evolution of the 10 Be concentration as a function of depth.

Quartz was isolated from crushed and sieved samples by dissolving all other minerals in a mixture of 1/3 HCl -2/3 H 2 SiF 6 acidic solution. Atmospheric 10 Be was then eliminated by successive HF sequential dissolutions (e.g., [START_REF] Brown | Examination of surface exposure ages of Antarctic moraines using in situ-produced 10Be and 26Al[END_REF]. After complete dissolution of the purified quartz in Suprapur HF, a reference 300 µg 9 Be spike (Merck Certipur) was added (e.g., [START_REF] Bourlès | Etude de la géochimie de l'isotope cosmogénique 10Be et de son isotope stable 9Be en milieu océanique. Application à la datation des sédiments marins[END_REF]. Beryllium was separated from these solutions by successive solvent extractions and alkaline precipitations. All 10 Be measurements were performed in June 2002 by accelerator mass spectrometry at the Tandétron facility, Gif-sur-Yvette, France [START_REF] Raisbeck | The AMS facility at Gif-sur-Yvette: progress, perturbations and projects[END_REF] against the NIST SRM4325 standard of certified 10 Be/ 9 Be value of 2.68.10 -11 . However, together with the absolute measurement of this ratio by Nishiizumi et al. (2007), the half live of 10 Be has also been more precisely measured in 2010 by two independent teams (Korschineck et al; 2010, [START_REF] Chmeleff | Determination of the 10 Be halflive by multicollector ICP-MS and liquid scintillation counting[END_REF]. Therefore, to incorporate these two recent advances, the original ratios have been recalculated, and the measured 10 Be concentrations presented in Table 1 are now consistent with an assigned NIST value of 2.79.10 -11 and a 10 Be half-life of 1.387 Ma.

-Minimum exposure ages have been calculated from: the measured concentrations extrapolated to surface assuming neutron attenuation only (see Table 2), the CREp (Cosmic Ray Exposure program) [START_REF] Martin | The CREp program and the ICE-D production rate calibration database: A fully parameterizable and updated online tool to compute cosmic-ray exposure ages[END_REF] with a production rate at sea level high altitude of 4.02 ± 0.16at/g/a [START_REF] Borchers | Geological calibration of spallation production rates in the CRONUS-Earth project[END_REF], and using the polynomials of [START_REF] Stone | Air pressure and cosmogenic isotope production[END_REF], a standard atmosphere and a custom geomagnetic database based on authigenic 10 Be/ 9 Be reconstruction from [START_REF] Simon | Authigenic 10Be/9Be ratio signatures of the cosmogenic nuclide production linked to geomagnetic dipole moment variation since the Brunhes/Matuyama boundary[END_REF]. The muon scheme is based on [START_REF] Braucher | Production of cosmogenic radionuclides at great depth: a multi element approach[END_REF], which uses an altitudinal scaling only and sea level production rates of 0.012 at/g/a and 0.039 at/g/a for slow and fast muons, respectively. Maximum denudation rates 12 have been similarly calculated by setting the time to infinite. All necessary equations (eq.1 and eq. 2) are presented in the Supplementary Material.

To combine and compare the different sample concentrations, the method proposed by [START_REF] Ward | Procedures for comparing and combining radiocarbon age determinations: a critique[END_REF] was then applied (Table 2, Figure 7, andFigure Sup4 in Supplementary Material): -This method is based on Chi-square analysis. To obtain the most adequate number of samples (n) per surface (i.e. samples that are statistically undistinguishable), the 0.05 critical value for a Chi-square with (n-1) degrees of freedom is calculated and compared with the theoretical value given by the chi-square table. If the calculated value for these n samples, is lower than the theoretical one, then all samples are used to calculate a weighted mean concentration; otherwise outliers are rejected until the distribution passes the test and the weighted mean concentration is calculated with the remaining samples (i.e. the n samples less the outliers ) (Figure 7). We applied the method described above, on the oldest alluvial surface S2, considering the 12 samples, and found a chi2 value of 98.2. This is much higher than the theoretical one of 19.68. It appears that for boulder B, the concentration measured in the lowest sample B3 is larger when extrapolated to the surface than that measured in the sample B1 collected at the surface. Moreover, sample B2, collected in between leads to a depleted concentration once extrapolated to the surface. The decrease of the in situ-produced 10 Be concentration at the core of the boulder and its increase at the bottom suggests that the boulder has been flipped and exposed from one face to another at least once before its present exposition (e.g. [START_REF] Carretier | Is it possible to quantify pebble abrasion and velocity in rivers using terrestrial cosmogenic nuclides[END_REF]. We therefore consider boulder B as an outlier. Similarly, the wide dispersion of the in situ-produced 10 Be concentrations measured in the samples from boulder E most likely indicates a complex exposure 13 13 history. Boulder E is therefore considered as an outlier. Thus, considering these concentrations extrapolated to the surface for the samples A, C and D, a chi2 value of 6.11 is calculated, yielding to a minimum surface cosmic ray exposure age of 61.0 ± 3.1 ka for surface S2 assuming negligible denudation rate (Table 2; note that because 10 Be surface concentrations of samples A, C and D are similar, inheritance can be neglected).

The same approach was applied to the concentrations obtained from the younger surface S1. The surface concentrations (F1, G1, I1) and the concentrations extrapolated to the surface of the samples collected at depth (F2, G2, and I2) are internally consistent considering all the samples (F, G, H, I) excluding here also inherited 10 Be from previous history. A mean surface concentration of 121.7 ± 7.4 kat/g is thus calculated, yielding to a minimum surface exposure age of 35.3 ± 2.0 ka. (Table 2).

However, one can argue that using a negligible denudation rate may not be right in such environment. Therefore, taking benefit of the multiple sampling within the same boulder, a maximum denudation rate can be determined for each sample from its measured 10 Be concentration. This is presented in Table Sup1 in Supplementary Material. One can observe that top and bottom denudations rates are similar within uncertainties. This may suggest that boulders are closed to steady state. To estimate their minimum exposure age, a so called integration times [START_REF] Lal | Cosmic ray labelling of erosion surfaces: In-situ nuclide production rates and erosion models[END_REF][START_REF] Rizza | Slip rate and slip magnitudes of past earthquakes along the Bogd left-lateral strike-slip fault (Mongolia)[END_REF] can be determined from each sample concentration, as well as mean maximum denudation rates and their related integration time for surfaces S1 and S2, from the 10 Be mean surface concentrations of S1 and S2 (Table 2). This yields denudation rates of 19.59 ± 1.18 m/Myr and 10.18 ± 0.78 m/Myr, and integration times of 37.1 ± 2.2 ka and 67.4 ± 5.1 ka for S1 and S2, respectively.

We also modeled the depth profile 10 Be concentrations obtained in the soil pit dug in S1 surface following Braucher et al (2009) and using the approach of [START_REF] Hidy | A geologically constrained Monte Carlo approach to modeling exposure ages from profiles of cosmogenic nuclides: An example from Lees Ferry, Arizona[END_REF]. The best fit (associated to the lower chi2) yields a minimum age of 38.49 (+ 14.43 / -3.45) ka (Figure Sup5 and Table Sup2 in Supplementary Material). This is consistent with the integration time calculated for S1. Note that the inheritance is totally negligible within the soil pit material suggesting that S1 emplaced rapidly with short transport time. This is consistent with Vargas et al. (2006)'s observation indicating that the period 35-15 ka corresponded to a period of strong alluviations.

Slip rates

To determine the slip rate along the Mejillones fault, we used the minimum surface exposure ages calculated above without denudation, which are very close to the integration ages. Dividing S2 and S1 mean offsets (12.8  1.3 and 5.1  0.8 m, respectively) by their minimum surface cosmic ray exposure ages (61.0 ± 3.1 ka and 35.3 ± 2.0 ka, respectively) yield maximum vertical slip rates of 0.21 ± 0.03 mm/yr and 0.15 ± 0.03 mm/yr, respectively. These results suggest that the vertical slip rate along the Mejillones fault has slightly decreased through time during the past 60 ka. However, given the uncertainties, these two estimates are not significantly different, and a mean 0.18 mm/yr vertical slip rate can be proposed for the Mejillones fault since the mid-upper Pleistocene.Where the fault cross cut the present shoreline, in the northern part of the peninsula, we could observe the main fault plane outcropping (see Figure 3B for location). There, we measured a fault azimuth of N008°E, a dip of 63°E, and a pitch for the slickensides of 84°N. This allows calculating an actual mean slip rate along the Mejillones fault of 0.20 mm/yr. 15

Coseismic displacements along the AFS

Along the Mejillones fault, we observed few subtle cracks within S0 surface that can be interpreted as surface breaks associated with the Mw 8.1, 1995 Antofagasta earthquake (See section 3.1; Figure 4B). However, the cumulative displacements observed along the fault within S1 and S2 surfaces are clearly the result of larger dislocations. After the Antofagasta earthquake, [START_REF] Delouis | Recent crustal deformation in the Antofagasta region (northern Chile) and the subduction process[END_REF] observed surface ruptures with decimetric vertical displacements along the Paposo fault (the southern extension of the Salar Del Carmen fault) 100 km southwards the Mejillones peninsula (Figures 8A,B), within a site situated in the middle of the epicentral area of the 1995 Mw8.1 Antofagasta earthquake (see Figure 1A for location). Looking for the extension of this surface rupture further north along the Paposo fault, and the northern Salar del Carmen fault, we observed an outstanding fault scarp predating the 1995 earthquake (Figure 8C).

Given the absence of beveled surfaces above and beneath the scarp, and given also the fact that we did not observed stepped surfaces, we interpreted this fault scarp as the result of a single coseismic surface dislocation, and not the result of cumulative offsets. A detailed topographic survey of the fault scarp allowed us determining a co-seismic vertical offset of 2.5 ± 0.5 m (Figure 9), a value which is of the same order than the offsets estimated along the Mejillones fault by [START_REF] Cortés | Paleoseismology of the Mejillones Fault, northern Chile: Insights from cosmogenic 10Be and optically stimulated luminescence determinations[END_REF] from paleoseismic investigations. Therefore, the question of the magnitude of the event during which such a vertical displacement occurred along the Salar del Carmen fault arises if we assume that it is contemporary with a subduction earthquake, as is the 1995 fault rupture described by [START_REF] Delouis | Recent crustal deformation in the Antofagasta region (northern Chile) and the subduction process[END_REF] farther south along the Paposo fault.

DISCUSSION

Our study allows determining a mean slip rate along the Mejillones fault of 0.20 mm/yr during the past 60 ka. This is lower than the mean ~0.6 mm/yr slip rate estimated by [START_REF] Cortés | Paleoseismology of the Mejillones Fault, northern Chile: Insights from cosmogenic 10Be and optically stimulated luminescence determinations[END_REF] 16 for the time period comprised between 35 and 3 ka. However, it should be noted that [START_REF] Cortés | Paleoseismology of the Mejillones Fault, northern Chile: Insights from cosmogenic 10Be and optically stimulated luminescence determinations[END_REF] detailed estimates for the time intervals 35-14ka and 3ka-Present (respectively 0.46 ± 0.41mm/yr and 0.22 + 0.06 mm/yr, respectively) are not significantly different with our result. Therefore, we cannot confirm whether or not, the Mejillones fault has known an acceleration between 14 and 3 ka (i.e. 0.6 mm/yr) as proposed by [START_REF] Cortés | Paleoseismology of the Mejillones Fault, northern Chile: Insights from cosmogenic 10Be and optically stimulated luminescence determinations[END_REF]. On the other hand, the 0.2 mm/yr differential uplift estimated by [START_REF] Ortlieb | Late Quaternary Coastal Changes in northern Chile: Guidebook for a fieldtrip[END_REF] over the past 400 ka, between the Pampa Mejillones (uplift rate of 0.2 mm/yr) and the Morro Mejillones faulted block (uplift rate of 0.4 mm/yr) tends to confirm that the long term slip rate is closer to 0.2 mm/yr than 0.6 mm/yr.

If we consider that the faults within the AFS are moving contemporaneously with the strong subduction earthquakes, as it has been observed for instance during the 1995, Mw8.1 Antofagasta earthquake [START_REF] Delouis | Recent crustal deformation in the Antofagasta region (northern Chile) and the subduction process[END_REF], we may assume that the recurrence of surface ruptures along the AFS should be the same as the recurrence of subduction earthquakes. With a convergence rate between the Nazca plate and Chile of about ~8 cm/yr [START_REF] Demets | Effect of recent revisions to the geomagnetic reversal time-scale on estimates of current plate motions[END_REF], an average 150 ± 100 years recurrence interval for strong (Mw≥8) subduction earthquakes can be calculated, assuming a regular recurrence intervals and a coseismic displacement between 4 and 23 m (range of displacement estimated from Wells and Coppersmith (1994)'s general statistical function, linking Mw≥8 and the average coseismic displacement). Combined with our long term slip rate estimated along the Mejillones fault (0.2 mm/yr), a 150 ± 100 yrs return period would yield a mean vertical offset of 3±2 cm per event along the fault if this one moves during each subduction earthquake occurring bellow the peninsula. This is not consistent with the cumulative fault scarps observed along the fault, which are clearly the result of larger dislocations in surface, as also shown by the paleoseismic studies carried out by [START_REF] Cortés | Paleoseismology of the Mejillones Fault, northern Chile: Insights from cosmogenic 10Be and optically stimulated luminescence determinations[END_REF] suggesting that the cumulative fault scarps along the Mejillones fault are the results of surface rupturing events with vertical displacements comprised between 0.6 and 3.0 m. These last displacements are of the 17 same order of the co-seismic slip that we measured along the Paposo or the Salar del Carmen faults (Figure 8). This suggests that while the mean recurrence time of fault activation in the AFS region would be the same than those of the strong subduction earthquakes (ca 150 yrs), the recurrence time of the AFS faults, taken separately is longer. This also suggests that the different AFS faults can ruptured individually after several Mw8 earthquakes and/or that subduction earthquakes larger than the Mw 8.1 1995 Antofagasta event can occur in the area. In consequence, we can wonder whether the 2.5 m high coseismic fault scarp observed along the southern part of the Salar del Carmen fault, eastwards of Antofagasta, could not be related to the last major historical subduction earthquake in 1877 (M 8.6).

CONCLUSIONS

Our morphotectonic study along the Mejillones fault -mapping the topography and the geology of the offset alluvial fans, and dating them with in situ-produced cosmogenic 10 Beallows determining that this main normal fault of the Atacama fault system (AFS) has a 0.2 mm/yr longterm slip rate. This result suggests that the AFS is composed of slow slip rates faults despite the fast rate (8 cm/yr) of subduction. However, the potential co-seismic slip observed along the different faults of the AFS, combined with the fact they seem to be moving contemporaneously with subduction earthquakesas during the Mw 8.1, 1995, Antofagasta earthquakesuggest that a larger earthquake could occur episodically in the area.

In this line of reasoning, the 2007 Tocopilla (M 7.7) earthquake having released only 2.5 per cent of the moment deficit accumulated on the interface during the 130 years period following the 1877 historical event (Béjar-Pizzaro et al., 2010) could be regarded as a possible precursor of a larger subduction earthquake rupturing partially or completely the 500-km-long 1. 10 Be data: the 10 Be concentrations have been calculated relatively to the NIST Standard Reference Materiel 4325, using an assigned ratio of 10 Be to 9 Be of (2.79±1.5) x 10 -11 . 10 Be uncertainties (1σ) include a 3% contribution conservatively estimated from observed variations in the standard during the runs, a 1 statistical error in the number of 10 Be events counted, uncertainty in the blank correction. Associated 10 Be/ 9 Be blank ratio was (4.37±1.32) x 10 -15 .

Table 2. 10 Be concentrations extrapolated to surface assuming a neutron attenuation length of 150 g/cm 2 and samples depths and densities as presented in Table 1. Minimum exposure ages have been calculated using the CREp online calculator (see text). To combine and compare the different samples concentrations, the method proposed by [START_REF] Ward | Procedures for comparing and combining radiocarbon age determinations: a critique[END_REF] was applied (see text).

  FIGURES CAPTIONS

Figure 2 :

 2 Figure 2: Field views of the Mejillones normal fault scarp within its northern part (A) and its

Figure 3 .

 3 Figure 3. Aerial photograph (A) and corresponding geological interpretation (B) of the

Figure 4 .

 4 Figure 4. A: Field picture showing S2 and S1 alluvial surfaces displaced by the Mejillones fault

Figure 5 .

 5 Figure 5. Oblique view of the Digital Elevation Models of the three studied sites O1, O2 and O3,

Figure 6 :

 6 Figure 6: A: Picture showing surface S2 sampled (boulder A) at site O2 within the footwall. B:

Figure 7 .

 7 Figure 7. 10 Be concentrations from boulders sampled on S2 (boulders A to E) and S1 (boulders F

Figure 8 .

 8 Figure 8. A and B: Pictures showing the surface rupture observed along the Paposo Fault within

Figure 9 .

 9 Figure 9. A: Contours map of the scarp observed along the southern part of the Salar del Carmen

18

North Chile seismic gap. Keeping on quantifying morphotectonic and paleoseismic features along the AFS may help better asses this issue.

ACKNOWLEDGMENTS

This study was supported by the Institut National des Sciences de l'Univers (INSU), France, the ECOS-Sud C00U01 program and the Servicio Nacional de Geología y Minería (SERNAGEOMIN) in collaboration with the Institut de Recherche pour le Développement (IRD), France. We thank Nicolas Marinovic (SERNAGEOMIN) and Gérard Hérail (IRD) for their support and Maurice Arnold for his help at the Tandétron facility of CNRS-Gif (France).

10 Be analyses were supported by the CEREGE and by the Géosciences Laboratory at the University of Montpellier 2. We thank Manuel Suarez, James Jackson, Sylvain Bonvalot, Serge Lallemand, and the team "Risks" in Geosciences Montpellier for fruitful discussions. We also thank three anonymous reviewers and the Associate Editor for their comments on the manuscript.