M Adelinet

C Domínguez

J Fortin

S Violette

Seismic-refraction field experiments on Galapagos Islands: a quantitative tool for hydrogeology

Keywords: seismic-refraction, volcanic rocks, water table, acoustic velocities, porosity Highlights -Seismic-refraction acquisition on volcanic islands

Due to their complex structure and the difficulty of collecting data, the hydrogeology of basaltic islands remains misunderstood, and the Galapagos islands are not an exception. Geophysics allows the possibility to describe the subsurface of these islands and to quantify the hydrodynamical properties of its ground layers, which can be useful to build robust hydrogeological models. In this paper, we present seismic refraction data acquired on Santa Cruz and San Cristobal, the two main inhabited islands of Galapagos. We investigated sites with several hydrogeological contexts, located at different altitudes and at different distances to the coast. At each site, a 2D P-wave velocity profile is built, highlighting unsaturated and saturated volcanic layers. At the coastal sites, seawater intrusion is identified and basal aquifer is characterized in terms of variations in compressional sound wave velocities, according to saturation state. At highlands sites, the limits between soils and lava flows are identified. On San Cristobal Island, the 2D velocity profile obtained on a mid-slope site (altitude 150 m), indicates the presence of a near surface freshwater aquifer, which is in agreement with previous geophysical studies and the hydrogeological conceptual model developed for this island. The originality of our paper is the use of velocity data to compute field porosity based on poroelasticity theory and the Biot-Gassmann equations. Given that porosity is a key parameter in quantitative hydrogeological models, it is a step forward to a better understanding of shallow fluid flows within a complex structure, such as Galapagos volcanoes.

-Sensitivity approach of mechanical parameters on inverted porosities 1. Introduction The hydrogeology of volcanic island is complex and remains poorly understood. Their hydrodynamic functioning relies on a disrupted geological setting up. Alternating active volcanic phases which contribute to the building up of the main edifice, combining basaltic lava flows or pyroclasts events with quiet periods, where weathering processes are dominant and alter the fresh outcropped basalts. As a result, the internal structure of volcanic islands looks like a layer-cake; where massive or fractured basaltic lavas are interlayered with weathered material (ashes, clayed soils, weathered basalts or pyroclasts…). Building an efficient flow model at island scale requires a large set of data: climatological and hydrological monitoring to constrain the water cycle on surface but also hydrodynamic parameters of the multi-layered system. in this context, porosity is one of the key parameters for such flow modelling as it controls to the groundwater storage. On Galapagos Islands, a previous study provides porosity and permeability data of soils at different elevation points [START_REF] Adelinet | The relationship between hydrodynamic properties and weathering of soils derived from volcanic rocks, Galapagos Islands (Ecuador)[END_REF]. However, data is limited to the surface. The present study suggests a new approach, which allows the estimation of the field porosity of the substratum (at deeper levels) using its acoustic properties. In order to investigate acoustic and petrophysical properties of Galapagos subsurface, we chose to use seismic-refraction methods. The advantage of seismic methods compared to electrical ones for instance is the direct relationship between acoustic velocities and porosity whereas electrical methods are generally only useful to determine position of water table in such complex geological structures [START_REF] Revil | The hydroelectric problem of porous rocks: inversion of the position of the water table from self-potential data[END_REF][START_REF] Astm D | Standard Test Method for Laboratory Determination of Pulse Velocities and Ultrasonic Elastic Constants of Rock[END_REF].Indeed, accurate knowledge of seismic velocities helps to estimate the porosity of the groundwater formations. For instance in [START_REF] Garambois | On the use of combined geophysical methods to assess water content and water conductivity of near surface formations[END_REF] porosity is inferred by using both P and S waves compared with GPR results. We are aware that only combined geophysical approaches can provide well-constraint hydrogeological models. However, as many studies have been performed on Galapagos islands during the last decade, we chose to focus our study only on the acoustic and petrophysical behaviour of both studied islands, Santa Cruz and San Cristobal. Seismic refraction methods are generally used in oil industry for static corrections in processing of seismic reflection., and commonly in engineering applications (e.g. see [START_REF] Khalil | Engineering applications of seismic refraction method: A field example at Wadi Wardan, Northeast Gulf of Suez, Sinai, Egypt[END_REF]). In hydrogeology, it is less employed due to the lack of resolution in depth and also due to the cost of an intensive processing. Moreover, seismic-refraction has some limitations, such as, when low-seismicvelocity layers are overlain by high-seismic-velocity layers [Haeni, 1981]. Nevertheless, seismicrefraction surveys have been used to describe the velocity structure of some basaltic islands, such as Canary Island [START_REF] Bosshard | Crustal structure of the western Canary Islands from seismic refraction and gravity data[END_REF][START_REF] Banda | Features of crustal structure under the Canary Islands[END_REF] or Faeroe Islands [START_REF] Pálmason | Seismic refraction measurements of the basalt lavas of the Faeroe Islands[END_REF]. A very interesting study has been performed in Iceland combining seismic tomography field experiments with laboratory ultrasonic measurements to obtain a comprehensive picture of the velocity systematics according lithology [START_REF] Grab | Seismic velocity structure of a fossilized Icelandic geothermal system: A combined laboratory and field study[END_REF]. Generally, refraction survey carried out on 2D lines allows the mapping of velocities on a depth profile. Different limits could be seen as refractors for volumetric waves: interface between weathered and unweathered rocks, fracture areas and limits between dry and saturated layers [START_REF] Mari | Geophysics of reservoir and civil engineering[END_REF]. Indeed, water table is a very effective refractor with P-wave velocities in saturated rocks generally more than 1500 m/s [START_REF] Kearey | An introduction to geophysical exploration[END_REF]. In our case, seismic refraction is used to assess subsurface acoustic properties. In such case, targets are close to the surface, less than 50 m deep, and therefore issues on depth resolution and cost are overcome. Interpretation of layer velocity can be assured as we expect an increasing gradient of velocity with depth at this scale. Indeed, if we consider the largest scale of the entire island, this assumption could not be validate due to the structural alternating sequence of lava flow deposits and soil/weathered rock, which could potentially results in a low velocity layer underlying a high velocity layer. In this paper, rather than just present the velocity structure of the subsurface, we propose to go a step further by using velocity data to derive porosity. The present study is included into the Galapagos Islands Integrated Water Studies project (GIIWS), which begins in 2003 [d'Ozouville, 2007;d'Ozouville et al., 2008ad'Ozouville et al., , 2008b;;[START_REF] Auken | An integrated study of the hydrogeology of volcanic islands using helicopter borne transient electromagnetic: Application in the Galápagos Archipelago[END_REF]Pryet, 2011a[START_REF] Pryet | 3D resistivity gridding of large AEM datasets: a step toward enhanced geological interpretation[END_REF]Pryet et al., 2012a[START_REF] Pryet | Quantification of cloud water interception along the windward slope of Santa Cruz Island, Galapagos (Ecuador)[END_REF][START_REF] Violette | Hydrologic and thermal modelling of an active volcano: the "Piton de la Fournaise, La Réunion Island[END_REF][START_REF] Dominguez Gonzalez | Comparison of deep percolation rates below contrasting land covers with a joint canopy and soil model[END_REF]. Note that before this project, no baseline data existed for hydrogeological understanding on Galapagos. The keypoint of the GIIWS is to use several datasets in order to build robust conceptual hydrogeological models for Santa Cruz and San Cristobal islands. Dataset are collected using in-situ measurements and indirect data, such as geophysical measurements. The integration of a large datasets is very important to understand the hydrogeological working of such complex structure as described on a Micronesian example by [START_REF] Ayers | Hydrogeology of an atoll island: a conceptual model from detailed study of a Micronesian example[END_REF], on La Reunion [START_REF] Violette | Hydrologic and thermal modelling of an active volcano: the "Piton de la Fournaise, La Réunion Island[END_REF], on Mayotte island [START_REF] Vittecoq | Influence of successive phases of volcanic construction and erosion on Mayotte Island's hydrogeological functioning as determined from a helicopter-borne resistivity survey correlated with borehole geological and permeability data[END_REF] andon Martinique Island [START_REF] Vittecoq | Heterogeneity of hydrodynamic properties and groundwater circulation of a coastal andesitic volcanic aquifer controlled by tectonic induced faults and rock fracturing -Martinique Island (Lesser Antilles -FWI)[END_REF]. The present paper is structured as follow: First, the acquisition on sites and the layouts on both islands are extensively presented. Second, we figure out the different 2D velocity profiles deduced from raw seismic data. Third, we interpret our velocity data in terms of porosity data in the poroelasticity framework. Finally we discuss the results in a global hydrogeological framework.

Material and methods

Sites description

The archipelago is located on the Ecuador at about 1000 km west from the South America continent. Santa Cruz and San Cristobal Islands are respectively in the central and eastern part of the archipelago. Geologically they belong to the central sub province [Mc Birney and Williams, 1969;[START_REF] Bow | Geology and petrogeneses of the lavas of Floreana and Santa Cruz Islands: Galápagos Archipelago[END_REF][START_REF] Geist | Evolution of Galápagos magmas: Mantle and crustal fractionation without assimilation[END_REF]]. Both islands are built from similar basaltic rocks and are exposed to similar climate conditions. However there are some differences between islands. In terms of age, San Cristobal is older than Santa Cruz with the more recently erupted lavas dated of 2.4 M.y. [START_REF] Geist | Geology and petrogenesis of lavas from San Cristobal Island, Galapagos archipelago[END_REF] against 1.3/0.95-0.05 M.y. on Santa Cruz [START_REF] White | Petrology and geochemistry of the Galápagos Islands: Portrait of a pathological mantle plume[END_REF]. The weathered cover is also different: 1 m on Santa Cruz compared to 10 m on San Cristobal [START_REF] Geist | Geology and petrogenesis of lavas from San Cristobal Island, Galapagos archipelago[END_REF]. On both islands, physical properties of soils evolved according to the elevation and the rainfall regime [START_REF] Adelinet | The relationship between hydrodynamic properties and weathering of soils derived from volcanic rocks, Galapagos Islands (Ecuador)[END_REF][START_REF] Violette | Hydrologic and thermal modelling of an active volcano: the "Piton de la Fournaise, La Réunion Island[END_REF]. Soils are thicker, less porous and less permeable in altitude than near the coast.

In the frame of this study, two seismic campaigns have been carried out on Galapagos Islands. Three sites on Santa Cruz (SZ) Island were investigated in 2011 whereas six sites on San Cristobal (SC) Island were done in 2013. Experiments were carried out at the same time of the year, in the beginning of the cool season (July). Figure 1 presents maps of the sites location on both islands. Sites have been chosen according to their altitude and distance to sea. We defined three types of location: coastal, mid-slope and highlands (Table 1). Mid-slope and highland sites are always located on the windward side of the islands. Ground floors of sites present also different aspects (Figure 2): unweathered or weathered basalts, scoria cone, pyroclasts, soils, etc… Coastal sites have been chosen to estimate wave velocities of basaltic rocks where water table position is known, e.g. limit between dry basalts and basalts saturated with salt water due to ocean intrusion. Given the expecting layering in these sites, interpretation should be much easier: dry versus water saturated basalts or scoria. Windward highlands sites have been studied mainly to measure the thickness of soils, which is a key parameter to quantify storage in the water cycle and water fluw as a recharge to deep aquifers. Mid-slope sites have only been investigated on San Cristobal because the perched aquifers depth is compatible with the penetration depth of seismic-refraction method (several tens of meters). Indeed the emergence of springs indicates the existence of shallow groundwater on San Cristobal [Pryet et al., 2012a;[START_REF] Dominguez Gonzalez | Comparison of deep percolation rates below contrasting land covers with a joint canopy and soil model[END_REF]. On the contrary, Santa Cruz Island presents a buried perched aquifer, identified by helicopter borne geophysical method at more than 100 meters depth [START_REF] Adelinet | The relationship between hydrodynamic properties and weathering of soils derived from volcanic rocks, Galapagos Islands (Ecuador)[END_REF][START_REF] Auken | An integrated study of the hydrogeology of volcanic islands using helicopter borne transient electromagnetic: Application in the Galápagos Archipelago[END_REF]Pryet et al., 2012a], a depth that cannot be reached with sledge-hammer source.

Let now detailed the site location and the layout for both islands. SZ1 is very close to the seashore (Figure 1), on a 6.5 m high cliff. The seismic line was deployed perpendicular to the seashore. SZ2 is an urban site; the line follows the drawing of a future road that leads to a residential area. As observed in the picture of figure 2, the road is without asphalt. The line was also deployed perpendicularly to the seashore. SZ3 is located in a pasture with flourishing vegetation. This longest investigated line was deployed south-eastward along the steepest slope. During acquisition, weather was foggy and rainy on this last site. SC1 is a San Cristobal site, the nearest from the coast. The line was deployed parallel to the seashore. SC2 is an urban site on a road without asphalt, the line was deployed perpendicularly to the seashore. SC3 is close to the airport and to a scoria cone on which an open-pit mine is settled. Figure 3 presents a detailed map of windward mid-slope sites on San Cristobal (SC4 and SC5). SC4 is located downstream of a major water catchment of the island (Cerro Gato water catchment) and near a perennial stream of the same name. An outcrop of pyroclastic material has been described just above the seismic layout [START_REF] Izquierdo | Volcanic island evolution and their implications for the groundwater conceptual model development in San Cristóbal (Galápagos)[END_REF]. SC5 is located on a watershed beside the one of SC4 (Figure 3b andc). The seismic line has been deployed near a temporary river-bed. We will discuss further the implication of the mid-slope sites locations in terms of velocities in section 3.2. SC6 is on the same watershed than SC4 but in the top of it (altitude around 600 m).

Survey layout and processing for seismic-refraction

As far as waves are generated by a sound source and travel across a layered media, different processes occur: wave can be diffracted, reflected or refracted on interfaces. Seismic-refraction methods measure the shortest time of a compressional wave to travel down from the source, through the ground, and back up to sensors placed on the land surface. By measuring the travel times of the sound wave and applying the Snell-Descartes's law that governs the propagation of sound, the geometry and characteristics of the subsurface geology and/or hydrogeology can be inferred. Therefore, field data consist of measured distances and seismic travel times. From this timedistance information, velocity variations and depths to individual layers can be calculated and modelled.

The seismic-refraction survey carried out on Galapagos Islands involves a total of 9 lines. Figure 4 describes the layout used on Santa Cruz (SZ) and San Cristobal (SC) sites. Each line was shot in both forward and reverse order. The sound source is a 4 kg sledge-hammer striking a Teflon plate. The shot points (SP on Figure 4) are spread on the lines with minimum 4 end-off and 1 central shots. 24 geophones of natural frequency of 10 Hz are regularly spread on the lines and measure the wave acceleration in the vertical direction. The distance between two geophones (inter-trace IT on figure 4) varies between 1 and 5 m, depending on the studied site. Table 2 summarizes the distances and acquisition parameters used for each site. No frequency filters were applied during data acquisition.

Objectives were different according to the sites. For the coastal ones, the objective is to visualize the interface between dry and water saturated rocks having the visual control of the sea level. We know that P-wave velocities are not sensitive to varying salt content of the groundwater, however as soon as we are closed to the sea shore we can interpolate this interface as the top of the salty wedge. The objective in highland site (SZ3) is to estimate the soil thickness and the weathered basalt properties for highland sites, thus a longer line was used for SZ3.At San Cristobal sites, we had the experience of Santa Cruz experiments made 2 years before. In addition targets were shallower, thus long lines were not useful. As a consequence, we used the same layout with a line length of about 100 m and a 2 m inter-trace for all the sites on San Cristobal. The same layout allows us to compare results from one site to another in terms of depth penetration and P-wave velocities. 5. Trace normalized amplitudes are used for display purpose. The data show good quality where first arrival refracted waves are easily identified. General workflow of seismic-refraction processing is presented on Figure 6. Data were analyzed using the software programs Pickwin (ver 5.1.1.2) and Plotrefa (ver 3.0.0.6) from the SeisImager software package developed by Geometrics Inc. The first step consists in picking first break arrivals of each shot point (5 for Santa Cruz sites and 9 for San Cristobal sites). It corresponds to P-wave arrivals. The picking was performed manually (Figure 6.1) using Pickwin. Then arrival times are plotted against source-to-geophone distances, which results in time-distance or traveltimes curves (Figure 6.2). The tomographic method using Plotrefa [START_REF] Zhang | Non-linear refraction travel-time tomography[END_REF]] involves the creation of an initial velocity model. Four parameters are required: the depth to top of the lowest layer, the minimum and maximum velocities and the number of layers. According geophone spacing, depth to top of the lowest layer is equal to 10, 20 or 40 meters (respectively for 1, 2 and 5 meters geophone spacing). Minimum and maximum velocities are fixed from the traveltimes curves, adding 30 % for the maximum value. Generally, these values are set to 300 m/s (propagation of the sound in air) and 3000 m/s, respectively. The number of layers is always fixed to 20. Then, the tomographic method involves iteratively tracing rays through the model, comparing the calculated travel times to the measured travel times, modifying the model, and repeating the process until the difference between calculated and measured times is minimized. At the end of the tomography process, a global Root Mean Square (RMS) error is calculated, integrating all errors between calculated and measured traveltimes in a least square meaning. This parameter assesses the quality of the workflow. Table 3 presents the final RMS errors obtained after the tomography processing of the nine studied profiles. The final RMS errors lie between 0 and 3 ms for most of the profiles (SZ1, SZ2, SC1, SC2, SC3, SC4 and SC5). Profiles SZ3 and SC6 show increased RMS errors, which can be related to poorer data quality caused by worse acquisition conditions (bad coupling due to wet clayey soils for these highland sites) compared to the other profiles.Moreover, RMS errors calculated for Santa Cruz sites are higher than the ones for San Cristobal. This can be consequence of the number of shot points : 5 for Santa Cruz and 9 for San Cristobal leading to a better convergence. The final cell size of the tomography model depends on the geophone spacing for the horizontal component and on the resolution during raytracing processing for the vertical one. The vertical resolution is globally equal to the half of the horizontal one. Table 4 summarizes the cell size obtained for the studied profiles according to the geophone spacing and also the total number of cells. 6 provides the reliability of raytracing and the tomography processing of our study case. Observed and calculated traveltimes are very closed (Figure 6a) and raytracing is coherent (Figure 6b). Note that we compare the first arrivals only for the times that are inside the profiles. Indeed cells that are off sides are not covered by enough rays and then cannot be solved. Final RMS error is low, less than 1 ms, decreasing sharply just after one iteration (Figure 6c). Assessing the coherence of raytracing allows us to be confident with the penetration depth. Following results have been checked to be reliable in terms of depth.

P-wave tomography results

Before our results are presented, we want to note once again that seismic-refraction is not the more efficient tool to study the variation of salt content into water. Therefore, we would infer the sea intrusion beneath island only for the coastal sites, where we know that the salty wedge is at the same depth that the groundwater level. Indeed, near the coast, the freshwater layer is very thin with only a few centimetres thick. At other sites, the visible water saturated layers can only be described as groundwater levels.

3.1. On Santa Cruz Island Figure 8 presents the tomography profiles of the three sites investigated on Santa Cruz Island. Resolution in depths and in P-wave velocities is different due to the variable intertrace between geophones. All profiles are shown with a flat surface because we succeed to deploy the different profile without important topographic variations. SZ1 site is located on a cliff (6.5 m above the sea level) with a direct visual checking of the seawater level. Outcrops exhibit a massive basaltic lava flow with numerous metric fractures. The velocity profile precisely presents a strong refractor 6.5 m beneath the ground surface, which corresponds to the sea water level position. A first layer of very low velocity (around 500 m/s) is very thin (about 2 meters width) and corresponds to unconsolidated material. The second and third layers of velocity respectively around 1600 m/s and 2400 m/s, are separated by an interface located at the sea water level position. Thus, these layers correspond respectively to dry and seawater saturated layers. This site allows the possibility to check the application of seismic-refraction to visualize ocean intrusion. On SZ2 profile, the seawater level is also visible, 20 m deep beneath the ground surface, which corresponds to the altitude of the site. We have access to more information than for SZ1 concerning the velocity layering thanks to a higher penetration depth (due to longer seismic line). The thickness of the first layer of low velocity (500 m/s) is about 2 meters and it could correspond to embankment. Indeed, the seismic line of SZ2 has been deployed on a future road. The second very low velocity layer (around 1000 m/s) is several meters thick with a lateral variation from 2 to 6 meters thick. It could correspond to burden unconsolidated material. Beneath this layer we find again a low velocity layer with a magnitude similar to the one obtained on SZ1 (around 1600 m/s), which is considered as a dry layer. Finally, beneath the seawater level, we identify a water saturated layer with higher velocity (2400 m/s) which can be the top of the salty wedge. The profile obtained on SZ3 site is the longest one of the present study. Distance between geophones is 5 meters, and the total spread is about 200 m. Due to the poor signal to noise ratio, the identification of the first break is more difficult in this profile. The global RMS error is quite high (Table 3: less than 4 ms actually).However, results are coherent with the other ones and the penetration depth is about 50 m. A first velocity layer with P-wave velocity around 350 m/s is observed from the surface to about 15 m deep. This velocity value is closed to the one of P-waves into the air (340 m/s), that is why we assume that this first layer is an unconsolidated one.. Beneath this very low velocity layer, a medium of about 15 m deep with a velocity around 1500 m/s is visible. This layer covers another velocity layer with a positive gradient according the depth, from 1500 m/s to 2400 m/s at the deepest parts of the profile. Mostly In the western part of the profile, beneath about 30-40 m depth we find a high velocity layer (about 2400 m/s). Regarding results obtained on SZ1 and SZ2 we interpreted it as a water saturated layer. But the resolution is too low to be sure on this last interpretation. To conclude on Santa Cruz, we identified three layers, from top to bottom: unconsolidated material with very low velocities near the surface, a dry layer associated to low velocities (around 1600 m/s) and a water saturated one associated to high velocities (2400 m/s). The last one can be easily identified at SZ1 and SZ2. Note that, even for the high velocity layer, the absolute velocity values remain low. It is probably because of the high degree of fracturing in the basaltic rocks, which is in agreement with the outcrop observations made on the cliff beneath SZ1. Compressive acoustic waves are very sensitive to fractures, and their velocity decrease sharply with increasing fracture density.

3.2. On San Cristobal island Figures 9 and 10 present the P-waves velocity models obtained respectively at the three coastal sites and at the three mid-slope and highland sites on San Cristobal Island. Similarly to Santa Cruz, three velocity layers could be identified from the six profiles as follows:

-A very low velocity layer (around 500 m/s) interpreted as unconsolidated material -A low velocity layer (between 1400 and 1700 m/s) interpreted as dry rocks -A high velocity layer (between 2400 and 2700 m/s) interpreted as water saturated rocks, note that water can be salty or fresh depending of the location, coastal or highland respectively. San Cristobal coastal sites (SC1, SC2, SC3) have a flat topography. The velocity models of these sites present a thin unconsolidated layer near the surface (Figure 9). The thickness of this layer is less than 1 meter for SC2 and SC3 and about 2 meters for SC1. The ocean intrusion is visible for the three coastal sites at sea level. As SC1 is the nearest site from the seashore, the saturation interface is a strong refractor and could hide a limit of alteration located just above. Thus the first layer with low velocity is probably in fact a mixed layer composed with unconsolidated material and unsaturated rocks. This can explain its relative high thickness according to other coastal sites. Once again, at this site the sea water intrusion appears beneath the real sea level. Moreover, as SC1 is the nearest site from the seashore, it could be due to a tide effect. SC2 presents the typical layering made with a very low velocity layer (less than 500 m/s) of unconsolidated material, a low velocity one (around 1600 m/s) of unsaturated fractured basalts and a high velocity one (around 2700 m/s) associated to a water saturated layer. The saturation interface is shown at sea level as expected. SC3 presents the same layering than SC2 with a smoother transition between water saturated and unsaturated layers. Moreover, the velocity of the water saturated layer is lower for SC3 profile (around 2400 m/s) than for the SC1 and SC2 ones(around 2700 m/s). Results obtained for mid-slope and highland sites are contrasting (Figure 10). First, these profiles present a topography (more than 15 % slope for SC5 and SC6 for instance). The velocity profile deduced at the highest site (SC6) presents only one very low velocity layer associated to soil, with one sublimit inside, at about 30 m depth. Site SC5 is located at a mid-slope on the windward side as well. The first layer is thinner (about 3.5 meters). This soil layer is above a distinguishable low velocity layer (around 1300 m/s), which could be associated to dry material according to the results obtained at the coast. In the northern part of the profile, a layer with higher velocity (2000 m/s) is recognisable in depth. It could be the top of a partially saturated level. The fully saturated layer would be located below. Finally, tomography profile obtained on SC4 is the most interesting one. From a hydrological point of view, this site is located downstream of a major perennial spring of the island (Cerro Gato spring) and near a stream fed from this spring. The spatial configuration of seismic layout regarding the hydrological aspects is presented on figure 3. The first soil layer is about 3 m deep, comparable to SC5. Beneath, the second velocity layer of about 10 m width corresponds probably to an unsaturated layer. But on the contrary to SC5, there is clearly a third layer on the bottom part of the profile with high velocity (about 2400 m/s). This layer might correspond to a saturated material. The velocity value is indeed similar to the corresponding one measured on SC3.

Interpretation of tomography results

Geological and hydrogeological interpretations

From an acoustic velocity point of view, both islands present the same layering in the subsurface (i.e. in the first dozen of meters): an unconsolidated level near the surface with very low velocities associated to unconsolidated material, a low velocity layer associated to dry volcanic material and a high velocity level in depth, associated to water saturated volcanic material. The first layer located very close to the surface could be interpreted as debris at the coast (embankment or desegregated volcanic material) and clayed soil in the highlands. For instance, on SC6, It could be interpreted as a thick soil above highly weathered volcanic material. Indeed, at the end of the formation of San Cristobal, the island was partially covered by pyroclastic deposits according to [START_REF] Geist | Geology and petrogenesis of lavas from San Cristobal Island, Galapagos archipelago[END_REF]. The alteration of such material may have formed the soil we see on the seismic-refraction data. Coastal sites with visual controls on the geology and the saturation interface allow the quantification and calibration of velocity in basalts, while it is water saturated or not. This is the case of SZ1 and SC1, providing basaltic lava flows on outcrops. Variations in velocities according to these calibration data supply new information. For instance, the saturated layers of SC3 and SC4 have lower velocities than the saturated layers of SC1 (2400-2500 m/s against 2700 m/s). It could be explain by the nature of the basement, SC3 and SC4 are very close to a scoria deposit cone where pyroclastic deposits can be observed around the areas (Figure 3). Now, if we compare both islands, on Santa Cruz, the velocities respectively in dry and water saturated basalts are around 1400 and 2400 m/s whereas on San Cristobal they are 1400 and 2700 m/s, respectively. The velocity difference between the two islands could be explained by the elder of each one. Indeed, as mentioned in the introduction, Santa Cruz is younger than San Cristobal and probably their fractures in basalts are fresher and not filled with other material (able to reduce the porosity). On San Cristobal, cooling fractures should be filled by alteration material or oversaturated fluid precipitation.

From a hydrogeological point of view, SC4 shows important results. The existence of a high velocity layer related to a water saturated level is clear. This could be interpreted as an aquifer.. This result is in accordance with previous helicopter-borne geophysical study, which reveals that both islands have prominent low resistivity layers in the range of 30-130 ohm.m beneath the windward highlands [d'Ozouville et al., 2008b;[START_REF] Auken | An integrated study of the hydrogeology of volcanic islands using helicopter borne transient electromagnetic: Application in the Galápagos Archipelago[END_REF][START_REF] Pryet | 3D resistivity gridding of large AEM datasets: a step toward enhanced geological interpretation[END_REF]Pryet et al., , 2012a]]. These resistivity values are of particular interest because they are characteristic of basalt saturated with water on other islands [START_REF] Lienert | An electromagnetic study of Maui's last active volcano[END_REF][START_REF] Descloitres | Electrical structure beneath the eastern collapsed flank of Piton de la Fournaise volcano, Reunion Island: Implications for the quest for groundwater[END_REF][START_REF] Krivochieva | Integrating TDEM and MT methods for characterization and delineation of the Santa Catarina aquifer (Chalco Sub-Basin, Mexico)[END_REF][START_REF] Vittecoq | Influence of successive phases of volcanic construction and erosion on Mayotte Island's hydrogeological functioning as determined from a helicopter-borne resistivity survey correlated with borehole geological and permeability data[END_REF]2015]. The Cerro Gato stream is located 80 m away from the left limit of the profile (A). This layer could be a saturated layer that corresponds to a freshwater perched aquifer supplying the hydrographic network. The water table is quite horizontal with a very low slope (2° towards the SW). corresponding to the direction of the topography. The orientation suggests that the aquifer is fed by the stream. Compared to SC4, profile obtained on SC5 does not show any fully water saturated layer. Figure 3 shows that SC5 is near from a temporary river bed. As we investigated sites during the end of the hot season (in July), the river bed was dry. Consequently groundwater could be at deeper levels and not supplied by surface network at this period of the year.

Inferring porosity data from tomography P-wave velocities

A key point to build quantitative hydrogeological model is the assessment of porosity of the different geological layers, porosity being the storage capacity of an unconfined aquifer. Acoustic velocities associated with effective medium modelling can bridge the gap between the geological lithology description and the porosity assessment at field scale [START_REF] Adelinet | Deriving microstructure and fluid state within the Icelandic crust from the inversion of tomography data[END_REF][START_REF] Adelinet | Effective medium modeling: How to efficiently infer porosity from seismic data?[END_REF]. Indeed, P-waves (and also of course S-waves) are very sensitive to the rock heterogeneity, such as grain contacts, pores, cracks or fractures, and Galapagos volcanic islands are a highly fractured medium. However, according to the scale and to the geophysical method, we do not investigate the same size of heterogeneity, it depends on the wavelength. In our case, due to the acquisition layout (sledge-hammer source), we have a central frequency of about 100 Hz [START_REF] Keiswetter | A field investigation of source parameters for the sledgehammer[END_REF][START_REF] Feroci | Some consideration on shallow seismic reflection surveys[END_REF]. Given an average P-wave velocity of 2000 m/s, it results in an average wavelength of 20 m for our seismic experiment. According to the standard test method for the determination of velocities [ASTM D2845-08, 2008, the wavelength should be at least three times the average heterogeneity size. Reversely, that means that our method provides velocities affected by heterogeneities of several meters length (fractures specially). Moreover, considering that seismic resolution is defined by a rule of thumb as the quarter of the wavelength, it means that our acquisition set-up would not be accurate to image anything with a thickness smaller than 5 meters. Besides, as seismic data reveal the existence of saturated layer in depth, we have to discuss another frequency effect on elastic properties of volcanic material. Indeed, at low frequency, the fluid pressure is constant and unaffected by the seismic waves, it is the drained regime. At higher frequencies, the assumption of fluid pressure equilibrium becomes invalid, the fluid pressure is locally uniform but changes when the waves pass through ; the regime is called the undrained regime. Both drained and undrained regimes are relevant of the poroelasticity theory [START_REF] Gassmann | Elasticity of porous media[END_REF][START_REF] Biot | Theory of propagation of elastic waves in a fluid-saturated porous solid. II. Higher frequency range[END_REF]Murphy, 1986]. According [START_REF] Cleary | Elastic and dynamic response regimes of fluid-impregnated solids with diverse microstructures[END_REF], the cut-off frequency (f c ) between the drained and the undrained states depends on the rock intrinsic permeability k (m²), the drained bulk modulus K d (equal to the dry modulus, Pa), the fluid viscosity η (Pa.s) and a flow length L (m) as follows:

݂ = ସ×× ఎ× మ (1)
The poroelasticity framework and the associated equations allow calculating the drained saturated moduli from the drained dry ones. To apply such theory we need to know in which frequency range the field experiments are. We can evaluate the cutoff frequency for our investigated sites considering a flow length, L, equal to a half wavelength (around 10 meters) and a drained dry bulk modulus around 10 GPa. The touchy parameter is the rock intrinsic permeability. At the laboratory scale, andesite basaltic and scoria flows may have the same permeability around 10 -12 m² [START_REF] Saar | Permeability-porosity relationship in vesicular basalts[END_REF]. We expect that on the field and due to fractures, the permeability is higher. Based on a hydrogeological modelling, [START_REF] Dominguez Gonzalez | Comparison of deep percolation rates below contrasting land covers with a joint canopy and soil model[END_REF] provide range between 10 -9 and 10 -11 m² for San Cristobal sites, which are in agreement with literature on permeability in highly fractured basalts [START_REF] Bear | Dynamics of fluids in porous media[END_REF]. Furthermore, Ingrebitsen & Scholl [1993] provide near surface horizontal permeabilities around 10 -10 m² for the Hawaiian Kilauea volcano which is a good analogue for Galapagos islands [START_REF] Violette | Hydrologic and thermal modelling of an active volcano: the "Piton de la Fournaise, La Réunion Island[END_REF]. Now, if we take an average value of 10 -10 m² for our field sites, we obtain a cutoff frequency of 0.4 Hz between the drained regime and the undrained one. This value is lower than the sledge-hammer frequency of 100 Hz and then validates our poroelastic approach that we will describe now, i.e. we will assume that fluid pressure in the saturated basalts is affected by the seismic wave (undrained regime). Our goal is to relate P-wave velocities to porosity for each investigated site. In an isotropic framework, P-wave velocity (VP) can be written as a function of bulk and shear moduli, respectively noted K and G, as follow:

ܸܲ = ට ା ସ ଷ ൗ ீ ఘ ⟺ ܭ = ߩܸܲ ଶ -4 3 ൗ ܩ (2)
Where ߩ is the rock bulk density.

Given that our interest is the difference between dry (subscript dry ) and saturated (subscript sat ) state, we introduce ܭ∆ as the difference between the saturated and the dry bulk modulus:

ܭ∆ = ܭ ௦௧ -ܭ ௗ௬ = ߩ ௦௧ ܸܲ ௦௧ ଶ -4 3 ൗ ܩ ௦௧ -ߩ ௗ௬ ܸܲ ௗ௬ ଶ + 4 3 ൗ ܩ ௗ௬ (3)
Considering that we remain in the poroelasticity framework, we can use the first equation of Biot-Gassmann, which announces that the saturated shear moduli is equal to the dry one. So:

ܭ∆ = ߩ ௦௧ ܸܲ ௦௧ ଶ -ߩ ௗ௬ ܸܲ ௗ௬ ଶ (4)
So, the variation in bulk moduli is directly related to the variation of P-wave velocity according saturation state. Velocities in dry and saturated medium are available on sites where a water table is visible (sea water level or perched aquifer), i.e. sites SZ1, SZ2, SZ3, SC1, SC2, SC3 and SC4.

Moreover, density is related to the porosity ∅ according the saturation state by:

ߩ ௦௧ = ሺ1 -∅ሻߩ + ∅ߩ ௪ (5) ߩ ௗ௬ = ሺ1 -∅ሻߩ (6)
Where ߩ , ߩ ௗ௬ ܽ݊݀ ߩ ௦௧ are the densities of the matrix, the dry rock and the saturated rock respectively.

Then, ∆K can be expressed as a function of porosity and variation in velocities (noted ܭ∆ due to the use of P-wave velocities):

ܭ∆ ሺ∅ሻ = ሺ1 -∅ሻߩ ൫ܸܲ ௦௧ ² -ܸܲ ௗ௬ ²൯ + ∅ߩ ௪ ܸܲ ௦௧ ² (7)
Besides, the second Biot-Gassmann equation relates the saturated and the dry bulk moduli as:

ܭ ௦௧ = ܭ ௗ௬ + ఉ² ∅ାሺఉି∅ሻ ಼ ಼ బ (8)
Where β is the Biot coefficient, K f and K 0 the water and matrix bulk modulus respectively. Besides, the Biot coefficient is equal to:

β = 1 - బ (9)
Then the variation in bulk modulus can be expressed in a different way as a function of porosity (noted ܭ∆ ீ due to the use of Biot-Gassmann theory):

ܭ∆ ீ ሺ∅ሻ = ఉ² ∅ାሺఉି∅ሻ ಼ ಼ బ (10)
No assumption is made on nature of porosity inclusions in the Biot-Gassmann equations: it is a macroscopic approach using the total pore space. We keep in mind that for water catchment the effective porosity is need. Our approach could be coupled with a structural geological field work to analyse precisely the fracture network. The next step would be to calculate the porosity proportion identified as crack porosity and introduce it in the effective medium modelling.

Finally from ( 7) and ( 10) we have a relationship between the difference of P-wave velocities, between saturated and dry states, and the porosity. The leading idea consists in minimizing a given objective function to determine the porosity. The objective function J is define as:

ܬሺ∅ሻ = ቀ1 - ∆ ೇು ሺ∅ሻ² ∆ ಳಸ ሺ∅ሻ² ቁ ଶ (11)
It is a one-term function, which quantifies the data mismatch in a least-squares sense. The minimization process is iterative; the process is repeated until the objective function is small enough. The optimization parameter is the porosity. Other parameters remain constant (Table 5), especially the Biot coefficient which is fixed to 1 in a first approach, i.e. we assume that K dry is much smaller than K 0 . Matrix modulus is taken equal to 25.6 GPa from experimental data performed in laboratory on Santa Cruz basalts as both islands present the same general lithology [START_REF] Loaiza Ambuludi | Étude des propriétés physiques et du comportement mécanique des Basaltes: étude cinétique, mécanistique et toxicologique de l'oxydation dégradante de l'Ibuprofène par l'oxydation avancée électrochimique et photochimique[END_REF]. Table6 summarizes input data (dry and saturated P-wave velocities) and output data (porosities) for the different sites on which we have both velocities. Note that SC3 and SC4 are respectively close to scoria cone and pyroclastic debris. Galapagos lava flows are described as riches in vacuoles due to gas trapping inside magma [START_REF] Loaiza Ambuludi | Étude des propriétés physiques et du comportement mécanique des Basaltes: étude cinétique, mécanistique et toxicologique de l'oxydation dégradante de l'Ibuprofène par l'oxydation avancée électrochimique et photochimique[END_REF]. Inverted porosities are in range of acceptable values for such effusive volcanic material [START_REF] Adelinet | Du terrain au laboratoire, étude des propriétés élastiques du basalte[END_REF][START_REF] Schaefer | Geomechanical rock properties of a basaltic volcano[END_REF][START_REF] Siratovich | Physical property relationships of the Rotokawa Andesite, a significant geothermal reservoir rock in the Taupo Volcanic Zone, New Zealand[END_REF]. Laboratory porosity measurements made on fresh Santa Cruz basalts provide values between 10 and 12 % [START_REF] Loaiza Ambuludi | Étude des propriétés physiques et du comportement mécanique des Basaltes: étude cinétique, mécanistique et toxicologique de l'oxydation dégradante de l'Ibuprofène par l'oxydation avancée électrochimique et photochimique[END_REF]. We can classify results of table 4 according three porosity facies: low porosity on basaltic San Cristobal sites (SC1 and SC2), intermediate porosity for scoria San Cristobal coastal site (SC3), and high porosity of Santa Cruz sites and mid-slope San Cristobal pyroclastic site(SZ1, SZ2, SZ3 and SC4). According to geological interpretation made before, the lowest porosity is associated with coastal basaltic sites (SC1 and SC2).

In order to discuss much precisely the absolute values of porosity we need to have a sensitivity approach on the parameters introduced in the objective function. We identified two major sinks of errors on porosity computations: the values of Biot coefficient arbitrary fixed to 1 and the matrix bulk modulus, which has been chosen from experimental laboratory data. Thus we perform simulations with different couples of Biot coefficient and matrix bulk modulus whereas other parameters remaining constant. We choose velocity values of site SC1 as representatives of good quality data.

Figure 11 presents the sensitivity results. The range of variation for tested parameters is 0.6 -1 for Biot coefficient and 10 -50 GPa for matrix bulk modulus in order to see the influence of both parameters on porosity. The computed porosity range is between 0 (computed for very low matrix bulk modulus, not realistic with volcanic rocks) and 0.17. The dependence of porosity with bulk modulus is the highest for Biot coefficient equal to 1. For realistic values of Biot coefficient and matrix bulk modulus (respectively 0.8-1 range and 20-40 GPa range, white square on Figure 9), the porosity is more dependent on the variation in Biot coefficient than on the matrix bulk modulus. Note that porosity increases with higher values of bulk modulus, which could be incoherent at the first glad. However, in this sensitivity approach the macroscopic field velocities remain constant (data from SC1 profile). Then, if the matrix is stiffer (high bulk modulus), the effective medium modelling needs to increase the porosity in order to fit with velocities which remains the same all over the tested values of bulk modulus.

In order to reduce this range of variability we introduce a variable Biot coefficient into our inverse modelling. For that we need to express the dry bulk modulus (K dry of eq. 9) as a function of known parameters. As we have only measure P-wave velocities on field, we assume a constant Poisson's ratio (noted ν) which is a more stable coefficient. In this case the dry bulk modulus (K dry ) is expressed from dry P-wave velocity and ν:

ܭ ௗ௬ = ଵା ଷሺଵିሻ × ρ ୢ୰୷ × VP ୢ୰୷ ² (12)
Injecting the dependency of ρ ୢ୰୷ with porosity we obtain:

ܭ ௗ௬ = ଵା ଷሺଵିሻ × ሺ1 -∅ሻ × ρ × VP ୢ୰୷ ² (13)
We can now reformulate the objective function by introducing the new K dry expression within the Biot coefficient (eq. 9). Table 7 presents the new porosities calculated from minimization of the new J function and with a Poisson's ratio of 0.25 accordingly with literature data [START_REF] Schultz | Limits on strength and deformation properties of jointed basaltic rock masses[END_REF][START_REF] Adelinet | Du terrain au laboratoire, étude des propriétés élastiques du basalte[END_REF][START_REF] Loaiza Ambuludi | Étude des propriétés physiques et du comportement mécanique des Basaltes: étude cinétique, mécanistique et toxicologique de l'oxydation dégradante de l'Ibuprofène par l'oxydation avancée électrochimique et photochimique[END_REF]. This table presents also the post-calculation of the Biot coefficient and the K dry values for each site. Inverted porosities are smaller than the previous ones computed with a Biot coefficient equal to 1, in agreement with the sensitivity approach given in Figure 9. Moreover, the post-processed Biot coefficient is relatively constant between sites, with a mean value of 0.90 and a standard deviation of 0.007.

Finally we plot the porosity in the map dry P-wave velocities versus saturated P-wave velocities (Figure 12) using the last optimization process with variable Biot coefficient. We show three domains. The first one is the low porosity area (porosity less than 8 %) in which basalts of coastal San Cristobal sites are present (SC1 and SC2). This area is interpreted as the domain of basaltic rocks weakly weathered. The second porosity area concerns intermediate values (between 8 and 13 %) in which SZ1, SZ3, SC3 and SC4 are present. Due to the presence of SC3 and SC4, we interpret this area as the domain of scoria and pyroclasts. The limit between the two areas could also be a limit of weathering: a weathered basaltic lava flow could have the same porosity. A third domain is also recognizable with higher porosity (more than 13 %). It could be attributed to very weathered or fractured rocks. The lowest computed porosity is the one of SC1 and SC2 with 5 %. It is the site closest the seashore with evidence of basaltic lava flows as underground. The most similar site on Santa Cruz is SZ1 with a porosity of 11 %. The difference between porosities should rely once again on the age of the islands. Lava flows on San Cristobal are older than on Santa Cruz. Porosity had time to be filled with secondary material, especially into the pores and fractures. On Figure 12, we notice that all the three sites of Santa Cruz are close and the highest value of porosity is obtained for SZ2 (14 %) which is the closest to the coast. Then we can assume that lower velocities recorded on Santa Cruz interpreted as high porosity values could be due to the high degree of fracturing of the young basaltic lava flows forming this island.

Conclusion

Seismic-refraction has been successfully used on Santa Cruz and San Cristobal Islands. It was a challenge to obtain acoustic information in the subsurface in such complex structural area. Moreover, penetration depths are enough to image different hydrogeological structures, salted wedge and a freshwater perched aquifer especially. Even if obtained P-wave velocity values obtained remains low we are able to distinguish different velocity layers and interpret them in terms of water saturation state. One step further should be the acquisition or the specific processing [START_REF] Foti | Experiments of joint acquisition of seismic refraction and surface wave data[END_REF][START_REF] Williams | Comparison of P-and S-wave velocity profiles obtained from surface seismic refraction/reflection and downhole data[END_REF][START_REF] Grelle | Seismic refraction methodology for groundwater level determination:"Water seismic index[END_REF][START_REF] Pasquet | 2D characterization of near-surface V P/V S: surface-wave dispersion inversion versus refraction tomography[END_REF][START_REF] Uhlemann | Landslide characterization using P-and S-wave seismic refraction tomography -The importance of elastic moduli[END_REF] to have access to shear-wave tomography, specially using the surface wave processing (Multi-Channel Analysis of Surface Waves, [Park et al., 1999]). Nevertheless, our study provides the mapping of low and high velocity layers in different geological context (basaltic lava flow or pyroclastic deposits) and at different altitudes. We went a step further by interpreting differences between dry and watersaturated P-wave velocities in terms of porosity thanks to the poroelasticity theory. From this work, absolute porosity values could be attributed to Galapagos subsurface material according elevation and geological facies: unweathered and weathered basaltic lava flows, scoria and pyroclast materials. This data could be very helpful in the building of flow models for the Galapagos Islands.

FiguresFigure 1 -

 1 Figures

Figure 2 -Figure 3 -

 23 Figure 2 -Pictures of investigated sites according altitude. 753

Figure 5 -

 5 Figure 5 -Shot gathers recorded for profile SC3. Trace amplitudes are normalized for display purpose.

Figure 6 -

 6 Figure 6 -General methodology of seismic-refraction processing made on the Galapagos data. SeisImager Geometrics softwares were used for first break picking (Pickwin module) and tomography iterative process (PlotRefra module).

Figure 7 -

 7 Figure 7 -Example of raypath processing to ckeck the tomography inversion quality. Example of SC2 profile.

Figure 8 -

 8 Figure 8 -Results obtained on Santa Cruz Island. Contrast between dry and water-saturated layers is observed. The interface is interpreted as the salted wedge near from the coast. For the highland site SZ3, the transition is smoother and can be interpreted more as a progressive water saturation profile in weathered material.

Figure 9 -

 9 Figure 9 -Velocity profiles obtained for coastal sites on San Cristobal Island. As for Santa Cruz profiles, a sharp contrast between low and high-velocity layer is observed corresponding to the contrast between dry and saturated layers in depth. The scoria site (SC3) presents velocities lower than basaltic sites (SC1 and SC2).

Figure 10 -

 10 Figure 10 -Velocity profiles obtained for land sites on San Cristobal Island. SC4 presents both dry and saturated velocity layer. The interface is interpreted here as the top of a freshwater saturated aquifer. Whereas SC5 is located also in the mid-slope of the island, it presents only dry velocities in depth. Limits of refraction resolution are achieved for highland site SC6 without any refractor in depth. Only unconsolidated material has been investigated (probably soil and highly weathered volcanic material).

Figure 11 -

 11 Figure 11 -Sensitivity approach on the effect of Biot coefficient and matrix bulk modulus on inverted porosities using SC1 input velocity data.

Figure 12 -

 12 Figure 12 -Modelling porosity from field P-wave velocities. Different San Cristobal lithologies are identified according inverted porosities: basaltic bedrocks (SC1 and SC2) or scoria / pyroclast cones (SC3 and SC4). The three sites of Santa Cruz are close, probably due to high fracturation in the bedrock.

Table 1 -Description of sites investigated by seismic-refraction on Santa Cruz and San Cristobal islands ISLAND SITES ALTITUDE m a.s.l. DISTANCE TO SEA m LOCATION TYPE DESCRIPTION

 1

	Santa Cruz	SZ1 SZ2 SZ3	6.5 20 393	20 1495 7950	Coastal Coastal Highlands	Perpendicular to seashore, on a 6.5 m high cliff, basaltic lava flow without soil Compact soil with scoria deposits Thick clayed soil layer with abundant vegetation, humid condition for acquisition
		SC1	1	135	Coastal	Parallel to seashore, compact soil with scoria debris (probably substratum is scoria)
	San Cristobal	SC2 SC3 SC4	8 10 160	350 915 1295	Coastal Coastal Mid-slope	Compact soil with rock debris (the substratum is probably basaltic rock) Close to a scoria cone (red mine), compact soil with debris (the substratum is probably scoria deposit) Soil and weathered basalt, near a pyroclastic outcrop
		SC5	230	1790	Mid-slope Silty to clay soil
		SC6	590	3870	Highlands Clayed soil (more than 4 m of thickness)

Table 2 -Layout description of seismic-refraction survey on Santa Cruz and San Cristobal islands

 2

									Acquisition
	ISLAND SITES	geophone spacing (m)	Seismic spread (m)	Geophones	Strike bases	Acquisition system	time / sample interval
									(ms)
	Santa	Cruz	SZ1 SZ2 SZ3	1 2 5	33 66 165	1C -10 Hz 1C -10 Hz 1C -10 Hz	5 5 5	DAQLINK III DAQLINK III DAQLINK III	500 / 0.5 500 / 0.5 500 / 0.5
	San Cristobal	SC1 SC2 SC3	2 2 2	96 96 96	1C -10 Hz 1C -10 Hz 1C -10 Hz	9 9 9	GEODE GEOMETRICS GEODE GEOMETRICS GEODE GEOMETRICS	500 / 0.5 500 / 0.5 500 / 0.5

Table 3 -Mean RMS errors obtained for the nine profiles after the tomography process. These values show that the models fit our travel time data well.

 3

	Site	SZ1	SZ2	SZ3	SC1	SC2	SC3	SC4	SC5	SC6
	RMS error (ms)	1.024	2.365	3.757	0.864	0.262	0.966	0.836	1.104	3.846

Table 4 -Horizontal (X size) and vertical (Z size) components of the tomography cell for each profile. The total number of cells is also mentioned.

 4

	Site	SZ1	SZ2	SZ3	SC1	SC2	SC3	SC4	SC5	SC6
	X size (m)	1	2	5	2	2	2	2	2	2
	Z size (m)	0.5	1	2.5	1	1	1	1	1	1
	Number of cells	486	630	517	432	432	432	432	432	432
	Using as example the SC2 profile, figure

Table 35 -Constant parameters

 35

	ρ_matrix	ρ _water	K_water
	(kg/L)	(kg/L)	(GPa)
	2.7	1	2.25

Table 46 -Input P-wave velocities and calculated porosity with constant Biot coefficient

 46

	sites	SZ1	SZ2	SZ3	SC1	SC2	SC3	SC4
	VP_sat (km/s) 2.39	2.33	2.39	2.69	2.72	2.53	2.43
	VP_dry (km/s) 1.30	1.36	1.35	1.34	1.38	1.33	1.38
	J evaluations	10	11	9	13	13	11	10
	J minima	1.44E-13 3.02E-12 7.12E-11 2.12E-11 4.31E-11 1.25E-10 7.62E-12
	inverted porosity	0.15	0.18	0.16	0.08	0.08	0.11	0.15

Table 57 -Porosity calculations using variable Biot coefficient and constant Poisson's ratio (0,25)

 57

	sites	SZ1	SZ2	SZ3	SC1	SC2	SC3	SC4
	VP_dry (km/s)	1.30	1.36	1.35	1.34	1.38	1.33	1.38
	new_porosity	0.11	0.14	0.12	0.05	0.05	0.08	0.11
	K_dry (GPa)	2.25	2.40	2.42	2.55	2.71	2.44	2.55
	Biot coefficient	0.91	0.91	0.91	0.90	0.89	0.90	0.90

Acknowledgements

This work is a part of the GIIWS project funding by the ANR-blanc 2010 GIIWS Ref. 601-01. The authors are grateful to the Galápagos National Park Service, the Charles Darwin Foundation, the Municipalities of Santa Cruz and San Cristóbal, and the Consejo de Gobierno del Régimen Especial de Galápagos (CGREG) for local collaborations and logistical support. Christian Domínguez participation was supported by a PhD Scholarship from the Ecuadorian Government through the Secretaría Nacional de Educación Superior, Ciencia y Tecnología (SENESCYT). The helpful comments and corrections from reviewers are gratefully acknowledged.


Hollow Beta zeolite single crystals for the design of selective catalysts

Ana Rita Morgado Prates a , Céline Pagis a,c , Frederic C. Meunier a , Laurence Burel a , Thierry Epicier b , Lucian Roiban b , Siddhardha Koneti b , Nicolas Bats c , David Farrusseng* a , Alain Tuel

Supporting Information Placeholder

We present a Pt@Beta catalyst with a unique, well-controlled location of metal nanoparticles that allows more efficient use of this rare, expensive metal in catalysis. The zeolite crystal has an inner cavity, leaving a thin zeolite shell where the metal nanoparticles are encapsulated, thereby ensuring a relatively small diffusional path length and high selectivity. Hollow Beta is obtained by a controlled dissolutionrecrystallization method. The location of the particles was not only revealed by electron tomography 3D reconstruction but was further confirmed by a model hydrogenation reaction of aromatics.

Zeolites are microporous materials with molecular sized pore dimensions 1 . Thanks to chemical and physical properties of these aluminosilicate crystals, zeolites have been widely used in petrochemical and oil refining industry as catalysts, adsorbents and molecular sieves. ZSM-5, Y and Beta zeolites with MFI, FAU and *BEA framework topologies, respectively, are amongst the most studied ones 2 . In particular, *BEA is a large pore zeolite, that can be easily synthesized with different Al contents, which makes this zeolite particularly attractive as catalyst in petrochemical industry. 3 The diffusion of molecules through the microporous network of the zeolite can be extremely slow, which would tend to reduce reaction rates and the useful fraction of crystals for catalysis, i.e., the volume of zeolite that is actually used 4 . In the case of zeolites containing precious group metal (PGM) within the zeolite framework, the result is that a considerable fraction of expensive metal is wasted 5 .

New classes of zeolite crystal morphologies have been developed in order to decrease transport limitations: zeolites with a shorter critical length such as nanosized zeolites 6,7 , or hierarchical micro-mesoporous zeolite, in which the mesopores work as "molecular highways" 8 . For both approaches, the diffusional path length is smaller, which is as expected. Conversely, the external surface increases, usually in association with the mesoporous specific surface area. Hollow zeolite single crystals are a new class of zeolite morphology consisting in a bulk zeolite containing a very large inner cavity 9 . The morphology of hollow single crystal zeolites departs from nanosized and hierarchical porous zeolites because the critical length (which corresponds to the shell thickness) is reduced while the external surface remains the same (Figure 1e). Note that state-of-the-art hollow zeolites are generally made of an aggregation of small seed crystals (Figure 1d) resulting in a polycrystalline micro-mesoporous shell. 9,10 Figure 1 -a) bulk zeolite, b) nanozeolites, c) hierarchical porous zeolite, d) hollow polycrystalline zeolite, e) hollow zeolite single crystal, and f) metal nanoparticles encapsulated within a hollow zeolite single crystal.

There are several protocols for encapsulating metal nanoparticles within zeolites while leaving no nanoparticles on the external surface of the crystal 5,[11][12][13] . This encapsulation approach has been shown to be particularly efficient against sintering by coalescence, while improving reaction selectivity 14 . In cases with fast reaction kinetics, we can assume that the nanoparticles located in the center of the crystal do not participate in the catalytic conversion because of internal diffusion limitations. Therefore, the effectiveness factor is impaired 15 . A catalyst that would not contain metal nanoparticles too far away from surface or in the center of the crystal would be an appropriate design to address the waste of expensive metals. Nonetheless, to the best of our knowledge, the synthesis of encapsulated nanoparticles with such a spatial gradient has not been achieved. Encapsulating metal nanoparticles within a hollow zeolite single crystal (Figure 1f) provides an alternative design in which all nanoparticles are located close to the surface as compared to a bulk crystal with nanoparticles encapsulated throughout the entire volume of the crystal. So far, only metal nanoparticles encapsulated in hollow MFI single crystals have been reported 5,16 . As expected, these catalysts have shown remarkable size-and diffusion-driven selectivity as well as sintering resistance. 9,[15][16][17] In the case of the hydrogenation of toluene on Pt nanoparticles encapsulated in an extremely thin zeolite shell, the effectiveness factor was demonstrated to have reached 100% versus 10% for the corresponding classic bulk impregnated zeolite 15 .

The formation of the hollow zeolite single crystal requires a bulk crystal with a spatial concentration gradient, for example Al zoning, in which the center of the crystal is richer in Si and can be preferentially removed by desilication, leaving a hollow zeolite structure. 18 Such a spatial gradient with an Al-rich surface is generally observed in ZSM-5 crystals prepared with TPA + cations as the template, but as far as we are concerned, it is not usually observed for other zeolites, such as *BEA. A typical alkaline desilication treatment over Beta zeolite creates a disorganized mesoporous pore system, but not an internal cavity 4 .

For a ZSM-5 with such an Al-rich surface, the treatment with an alkaline reagent such as Na2CO3 solution generates an inner cavity thanks to the selective desilication of the center of the crystal. Alternatively, desilication can also be performed with alkaline TPAOH solution at higher temperatures, which enables not only the desilication of the center but also the recrystallization of the dissolved species on the outer surface of the crystal (dissolution-recrystallization method), thereby improving the regularity of the zeolite walls. 15 The latter method is very appropriate for the encapsulation of nanoparticles inside the inner cavity 5 .

To the best of our knowledge, hollow Beta zeolite single crystals with encapsulated nanoparticles have never been reported. Indeed, in the case of Beta zeolite, the absence of Al zoning (or of any other property with a spatial gradient) prohibits the use of the dissolution-recrystallization strategy employed in the synthesis of the aforementioned hollow zeolite single crystals with encapsulated nanoparticles.

Herein, we present the synthesis of hollow Beta zeolite single crystals (Pt@Beta) in which Pt nanoparticles are encapsulated in the cavity and inside the zeolite shell. The Pt nanoparticles are demonstrated to be accessible only through the *BEA zeolite by means of a model catalytic hydrogenation test. The novel synthesis consists indeed in a dissolution-recrystallization process, but over CIT-6 (zincosilicate with *BEA structure) containing Pt nanoparticles, instead of using the same bulk crystal type, as in the methods presented before for MFI zeolite. In this final step, Beta zeolite is formed during the recrystallization stage around the CIT-6 surface, which is subsequently dissolved (Figure 2). The synthesis of Pt@Beta is inspired by a seed-assisted synthesis approach described by the literature 19 , which used CIT-6 crystals as Beta seeds (Figure 2). CIT-6 is a zincosilicate microporous material with the same *BEA framework topology as Beta. The preparation of the CIT-6 crystals is performed according to a procedure described in the literature 20 . A synthesis mixture is prepared with the following composition: SiO2 : 0.05 LiOH : 0.65 TEA : 0.03 Zn(CH3CO2)2 : 30 H2O (check Supplementary information for more details). After preparing an aqueous solution of LiOH, TEAOH and Zn(CH3CO2)2, Ludox HS40 colloidal silica is added dropwise. The resulting mixture is stirred for 1h at 40°C, to yield a clear solution. Crystallization is performed in a 48 ml Teflon-lined stainless-steel autoclave, under static conditions at 145°C for 4 days. The X-ray powder diffractogrammes of the as-synthesized samples show a typical *BEA pattern as expected, with neither impurities nor other crystalline phases (see Figure S1). SEM images of bulk CIT-6 crystals show a homogeneous distribution of crystals with tetragonal pyramidal morphology and sizes between 0.5 and 1.1 µm (Figures 3a andS2.

a) b)

Figure 3 -SEM images of bulk CIT-6 crystals (on the left) and Pt@Beta crystal (on the right). Because these samples were not coated with a metal layer, as in classic SEM sample preparation, one can obtain information on the inside of the crystals. The white/grey color depends on the density of atoms perceived by the electon beam, and therefore the number of electrons detected. This contrast directly informs on the hollow nature of the crystals.

As a second step of the process, these as-made bulk zincosilicate CIT-6 crystals are impregnated with a Pt precursor solution and then reduced by a NaBH4 solution. This sample is referred to as Pt/CIT-6. XRD results show that after impregnation, Pt/CIT-6 is highly crystalline. The two small peaks at 40° and 46° for Pt/CIT-6 are attributable to platinum nanoparticles; ICP analysis reveals a final Pt content of 2.6 wt.%. The location and size of platinum nanoparticles impregnated on the surface of CIT-6 were investigated by TEM (Figure S3), which shows that the platinum nanoparticles (particle size between ~10-20 nm) are successfully anchored on the external surface of the crystals. Nitrogen adsorption/desorption (77 K) performed on the calcined Pt/CIT-6 solid presents a type I isotherm typical of a purely microporous crystal (Figure S8).

As described before, in the case of hollow MFI crystals with encapsulated nanoparticles, the impregnation is carried out on the calcined bulk MFI crystals, followed by an H2 reduction, and finally the dissolution-recrystallization method to create the cavity. This process, however, cannot be applied to the CIT-6 crystals. Actually, if CIT-6 crystals are calcined before the hollow Beta synthesis, the CIT-6 is dissolved during the hydrothermal treatment before any Beta crystallization can occur on its surface. As-made CIT-6 crystals containing the organic template present the advantage of being more stable in alkaline medium, thereby allowing the crystallization of Beta zeolite on the surface before its complete dissolution.

Pt@Beta is obtained by hydrothermal treatment of the as-made Pt/CIT-6 in the presence of an organic-free aluminosilicate gel (SiO2 : 0.01 Al2O3 : 0.3 NaOH: 20 H2O) in which CIT-6 represents 10 wt.% of the total amount of silica source. Firstly, an aqueous solution of NaOH and NaAlO2 is prepared. Then, CIT-6 crystals are added to the solution, which is gently stirred until the dispersion becomes homogeneous. Finally, Aerosil 200 is slowly added with vigorous stirring. The hydrothermal treatment is performed at 150 °C for 22 h. During this time, Beta zeolite starts to recrystallize on the CIT-6 surface, which afterward is dissolved, leading to hollow single crystals.

As revealed by XRD and electron microscopy images (Figures S1,S2, S4, and Figure 3b and Figure 4), after hydrothermal treatment a highly crystalline hollow Beta crystal is obtained, in which the Pt nanoparticles are encapsulated inside the cavity and/or in the zeolite shell. The Pt@Beta presents a crystal size between 0.6 and 1.2 µm. The size of the cavity corresponds approximately to that of the dissolved CIT-6 crystal, leaving an average shell thickness of 90-190 nm. The observed crystal shape is a truncated tetragonal pyramidal morphology, pine tree shaped on each sides, which is typical for Beta zeolite crystal morphology 4 .

HRTEM images and STEM images over sliced hollow Beta crystals also clearly reveal the presence of mesopores in the internal part of the microporous wall (see Figure S6) with an average size between ~5 and 15 nm. Nonetheless, it must be noted that the external part of the hollow zeolite is entirely microporous, because the mesopores are not connected to the external surface. The homogeneity of the synthesis can be validated by SEM observations (Figure 3b and Figure S2). Specific sample preparation was used to highlight the hollow nature of the crystals (SI). When comparing SEM images of bulk CIT-6 and hollow *BEA under identical conditions, it is obvious that all *BEA crystals possess a large internal cavity, as evidenced by the contrast between the center and extremities of the crystal.

The N2 physisorption isotherm of Pt@Beta at 77 K (Figure S7) presents a type I adsorption branch at low pressure, which is typical of a microporous material, and a type IV branch with a hysteresis that closes at p/p° =0.53, which can be associated with condensation phenomena in bottle-neck pores [21][22][23] . Indeed, as mentioned above, HRTEM and STEM images over sliced hollow Beta crystals clearly show the presence of mesopores in the internal part of the microporous wall (see Figure S6). These results show that this material consists of a hollow single crystal in contrast to a polycrystalline aggregation such as published elsewhere 24 .

Chemical analysis of the Pt@Beta crystals gives Si/Al and Si/Zn molar ratios of 8 and 237, respectively. The remaining Zn concentration is very low, clearly indicating that the dissolution of CIT-6 during the hydrothermal treatment is almost complete. The remaining fraction of Zn could result from Zn species reincorporated either in the hollow Beta framework or in the pores as extra-framework species to compensate for the charge of the aluminum atoms.

The Pt loading on the Pt@Beta is 1 wt.%, and its dispersion is around 3%, as calculated from the particle size distribution determined from TEM images. TEM observations over a large population of Pt@Beta crystals, consistently show the presence of Pt NPs inside each crystal, and the absence of visible Pt NPs on the external border of the crystals, (see Figure 4 and Figure S4 for example). This is an evidence of the internal location of the Pt NPs. Transmission electron tomography with corresponding 3D reconstruction was carried out on an isolated hollow crystal, confirming the absence of visible Pt nanoparticles on the external surface (Figure 5 and Figure S7). On the other hand, the nanoparticles are clearly observed on the inner surface of the cavity, and inside the zeolite wall, in the mesoporous region, with an average diameter of 28 nm (Figure 5-b). Pt@Beta crystal show to be resistant to calcination at 750°C for 10 hours. The particle size distribution of the Pt NPs remains very similar regarding the as-made Pt@Beta (Figure S4 andS5).

The location of the nanoparticles inside the wall and cavity was also confirmed by a model hydrogenation reaction using reactants with different diffusion rates. The catalytic hydrogenations of toluene and mesitylene were chosen as model reactions to demonstrate that there are no particles at the zeolite surface. Toluene and mesitylene hydrogenations were carried out using a tubular quartz plug flow reactor placed in a tubular furnace at atmospheric pressure. A mass of 10 ± 0.1 mg of catalyst, held between quartz wool plugs, was used. Each substrate was fed separately using a saturator maintained at 0 °C for toluene and 17 °C for mesitylene, leading to partial pressures of 0.91 kPa and 0.20 kPa, respectively. The reactor effluent was analyzed in a 10-cm path length gas cell fitted in a Bruker Tensor 27 FT-IR spectrophotometer. The activity of Pt@Beta in the hydrogenation of these two substrates was compared with that obtained over a conventional Pt catalyst supported on silica (Pt/SiO2), for which all Pt nanoparticles are accessible without any diffusion constraints (Figure 6). The Pt/SiO2 data were reported in an earlier paper 16 and are reused here for the sake of comparison.

Figure 6 -Activity (turnover frequency) of Pt/SiO2 catalyst (black symbols) and Pt@Beta catalyst (red symbols), in the hydrogenation of toluene (circles) and mesitylene (triangles). Toluene pressure = 0.9 kPa, balance H2. Mesitylene pressure = 0.06 kPa and 0.2 kPa, for Pt/SiO2 and Pt@Beta respectively, balance H2.

The Pt/SiO2 catalyst was active for the hydrogenation of both toluene and mesitylene, showing very similar turnover frequencies (Figure 6). The turnover frequency and apparent activation energy values (e.g., TOFPt/SiO2 = 6 x 10 -2 s -1 at 80 °C and 55±3 kJ/mol for toluene, Figure S11) are consistent with the literature data for both reactions 16 .

Catalytic activity of Pt@Beta for toluene hydrogenation shows similar TOF values to those observed in the case of Pt/SiO2 (e.g., at 80 °C, TOFPt@Beta = 4 x 10 -2 s -1 ). Similarly, the value for the apparent reaction activation energy (Ea = 47 kJ/mol) is consistent with those reported in the literature for Pt/SiO2 16,25 (Figure S10). On the other hand, mesitylene conversion was negligible (conversion in the error range of ca. 0.01%), and the corresponding TOF is therefore essentially nil over Pt@Beta, in contrast to Pt/SiO2. Even though the mesitylene pressure used was slightly different for the two catalysts, this should not affect the TOF values since the reaction order with respect to the aromatic on Pt-based catalysts is typically zero or close to zero 25,[START_REF] Bond | Metal-Catalysed Reactions of Hydrocarbons[END_REF] .

The TOF values of Pt/SiO2 and Pt@Beta are similar. Based on our latest studies 15,16 , we can rationalize that no significant diffusional limitation occurs in the toluene hydrogenation (i.e., effectiveness is high) thanks to the small shell thickness. In contrast, we can conclude that severe diffusional limitation occurs for the bulkier mesitylene for Pt@Beta. Although mesitylene is able to enter into the *BEA pore system, its diffusion in the microporous framework is significantly slower than that of toluene due to its larger kinetic diameter, which is 0.86 nm compared to 0.61 nm for toluene (Figures S12-S13). The fact that toluene (but not mesitylene) is readily hydrogenated clearly demonstrates the absence of nanoparticles on the external surface. It also evidences the single crystal nature of this material, as well as the absence of a mesopore network connecting the external surface and the inner cavity. A mesopore system would allow the diffusion of mesitylene through the zeolite shell, and as consequence its hydrogenation.

Figure 2 -

 2 Figure 2 -Multistep synthesis of Pt@Beta starting from CIT-6

Figure 4 -

 4 Figure 4 -TEM images of Pt/CIT-6 (left) and Pt@Beta zeolite with Pt nanoparticles encapsulated inside the wall and inside the cavity (right).

Figure 5 -

 5 Figure 5 -3D models of a Pt@Beta zeolite crystal reconstructed from electron tomography. The zeolite is represented in grey and Pt nanoparticles are represented in yellow. External view (left), transparent view (right).

  a aIRCELYON, 2 avenue Albert Einstein, 69626 Villeurbanne Cedex bUniv. Lyon, INSA-Lyon, Université Claude Bernard Lyon 1, MATEIS, UMR 5510, CNRS, 7 avenue Jean Capelle, 69621 Villeurbanne Cedex, France cIFP Energies Nouvelles, Rond-point de l'échangeur de Solaize, BP3, F-69360 Solaize, France

Herein we have presented a Pt@Beta catalyst with a unique, well-controlled location of metal nanoparticles in the zeolite shell. This novel catalyst design makes it possible to use a greater portion of metals involved in the catalytic reaction, reducing the waste of expensive metals and improving the effectiveness factor.

Acknowledgements

The authors would like to thank Cecile Daniel for the scientific services provided for this study. The authors would like to thank CLYM (Centre Lyonnais de Microscopie) for granting access to the FEI ETEM Titan 300-80 keV and Jeol 2100F.

ASSOCIATED CONTENT Supporting Information

The Supporting Information is available free of charge on the ACS Publications website. Information on synthesis, XRD, SEM, TEM, N2 isotherms, toluene and mesitylene adsorption isotherms, catalytic tests, FT-IR transmission spectra, Arrhenius-type plot of catalytic results.

AUTHOR INFORMATION

Corresponding Author

david.farrusseng@ircelyon.univ-lyon1.fr

Notes

The authors declare no competing financial interests.

For Table of Contents Use Only Hollow Beta zeolite single crystals for the design of selective catalysts

Ana Rita Morgado Prates a , Céline Pagis a,c , Frederic C. Meunier a , Laurence Burel a , Thierry Epicier b , Lucian Roiban b , Siddhardha Koneti b , Nicolas Bats c , David Farrusseng* a , Alain Tuel a

We present a Pt@Beta catalyst with a hollow morphology and encapsulated metal nanoparticles (NPs) inside the zeolite shell, ensuring a small diffusional path length and high selectivity. Materials were obtained by a controlled dissolution-recrystallization method. The location of the NPs was not only revealed by electron tomography 3D reconstruction but was further confirmed by a model hydrogenation reaction of aromatics.