Petrophysical characterisation of volcanic ejecta to constrain subsurface lithological heterogeneities: implications for ediﬁce stability at basaltic volcanoes

Piton de la Fournaise (PdF) is an active and dominantly e ﬀ usive basaltic volcano. The geologic record preserves evidence of rare explosive eruptions, associated with the seaward sliding of the steep east ﬂank and the collapse of large calderas emplacing lithic-rich breccias. Breccias of PdF o ﬀ er the opportunity to sample a wide range of crustal lithologies. Unaltered blocks are representative of plutonic (gabbros), sub-volcanic (dolerites emplaced in sills and dykes), and volcanic (lavas) units. We document a large variability in density, porosity, P-wave velocity (dry and wet), and uniaxial compressive strength (UCS). The large variation in P-wave velocity and UCS results from the wide ranges in texture and lithology. Notably, some of the dense gabbroic units having experienced several cycles of natural reheating are comparatively weak. We infer that volcano instability should not be interpreted solely in terms of altered rock units. Rather, the large petrophysical heterogeneity of crustal rocks at PdF, and likely at many other basaltic volcanoes, must be considered when interpreting monitoring data and assessing volcano stability. Le Piton de la Fournaise (PdF) est un volcan basaltique actif et à dominante e ﬀ usive. Les archives géologiques conservent des preuves d’éruptions explosives rares, associées au glissement vers la mer du ﬂanc est escarpé et à l’e ﬀ ondrement de grandes calderas produisant des brèches riches en lithiques. Les brèches de PdF o ﬀ rent l’oppor-tunité d’échantillonner une large gamme de lithologies crustales. Les blocs non modiﬁés sont représentatifs des unités plutoniques (gabbros), subvolcaniques (dolérites placées dans sills et les dykes) et volcaniques (laves). Nous documentons une grande variabilité de la densité, de la porosité, de la vitesse de l’onde P (sèche et humide) et de la résistance à la compression uniaxiale (UCS). La grande variation de la vitesse de l’onde P et du UCS résulte des larges gammes de texture et de lithologie. Notamment, certaines des unités gabbroïques denses ayant subi plusieurs cycles de réchau ﬀ ement naturel sont relativement faibles. Nous en déduisons que l’instabilité du volcan ne doit pas être interprétée uniquement en termes d’unités de roches altérées. Au contraire, la grande hétérogénéité pétrophy-sique des roches crustales à PdF, et probablement dans de nombreux autres volcans basaltiques, doit être prise en compte lors de l’interprétation des données de surveillance et de l’évaluation de la stabilité du volcan.


Introduction
Caldera collapses and flank failures punctuate the history of volcanoes worldwide and represent major highly hazardous events in their evolution [Poland et al. 2017]. Moreover, vertical caldera collapse and flank deformation or collapse can be coeval and thus mechanistically linked [Hunt et al. 2018;Martí 2019]. Never-Corresponding author: dimuro@ipgp.fr theless, their link to magma transfer and storage in the plumbing system, together with the nature of weakness zones responsible for volcano collapses still needs to be fully elucidated.
For primarily effusive volcanoes, caldera and pit collapse events are usually related to magma withdrawal from the shallow part of the plumbing system during flank eruptions or lateral intrusion. The 2007 caldera-forming eruption of Piton de la Fournaise, the 2018 Lower East Rift Zone eruption of Kīlauea vol-cano (Hawai'i, USA) and the recent (2003)(2004)(2005)(2006)(2007)(2008)(2009)(2010)(2011)(2012)(2013) effusive events of Stromboli, all associated with draining of the summit lava lake/conduit and eventual collapse of the summit area are excellent examples Neal et al. 2019;Di Traglia et al. 2014]. The surface expression is expected to depend significantly on shape, volume, and depth of the magma reservoir. However, the nature of magma storage beneath most calderas is subject to significant speculation as storage zones can be complex and span a large range of depths , and references therein]. Moreover, volcano stability and shape evolution are strongly dependent on the rheology of the deforming crust. Volcanoes can undergo slow deformation and spreading over extended time periods, punctuated by events of rapid slumping and flank failure [De Vries and Francis 1997]. Instability factors include slope steepening and overloading, peripheral erosion, hydrothermal alteration, and rainfall infiltration and, next to magma emplacement, they can be triggered by hydrothermal pressurisation and large earthquakes (for reviews of processes and trigger mechanisms, see McGuire [1996] and Borgia et al. [2000]).
A combination of geophysical and geochemical approaches can often shed light on the present state of active plumbing systems beneath volcanoes. The results of such approaches are often interpreted in terms of deep zones of magma storage revealed by zones of low seismic velocity, overlain by shallower and often smaller low-velocity zones often linked to ground deformation. However, the detailed material state of these shallow zones, where most deformation and seismic sources are located, is still unclear. They may include small melt and fluid volumes (sills and crystal mushes), cooling dykes, faults, altered rocks, and circulating hydrothermal fluids.
Geophysical techniques serve to reveal clues about the internal structure of volcanic edifices and allow us to monitor possible changes. The measurements are commonly obtained using sensors at Earth's surface. The inversion of geophysical datasets is, however, dependent on the heterogeneity of the underlying crust built up from layers of variable thickness, lateral extent, and petrophysical properties. Thus for robust data interpretation, the relationship between the physical properties of subsurface rocks and their link to geophysical measurements must be constrained by accurate laboratory measurements. Seismic tomography utilises variations in seismic velocities, which in turn depend on a variety of factors such as spatially varying lithology, fluid and melt fraction, temperature, and anisotropy [Lees 2007]. Volcano inflation and a decrease in bulk edifice strength, eventually leading to eruptions or summit/flank failure, have been demonstrated to result in changes to the seismic velocity of the edifice [Brenguier et al. 2008;Olivier et al. 2019]. The vertical and lateral distribution of rock heterogeneities can influence fluid transfer and rock deformation and eventually the mechanism and spatial distribution of crustal seismicity. Notably, weak crust can deform almost aseismically. Magnetotelluric techniques rely on electrical conductivity, which is a function of rock porosity and connectivity, fluid and melt fraction, and composition [Chave and Jones 2012]. Geodesy measures the strain evolution of the Earth's surface, which is also typically related to the evolution of pressure sources (fluids, melts), the influence of which on the surrounding rocks is highly dependent on their physical properties [Dzurisin et al. 2019].
Many if not all basaltic volcanoes have experienced several collapse events (lateral flank and vertical caldera collapses). Piton de la Fournaise (PdF) is one of the most active basaltic volcanoes on Earth, characterised by predominantly effusive activity from its summit or ephemeral vents on its flanks that are structurally affected by rift zones [Roult et al. 2012;Michon et al. 2015]. Several nested U-shaped calderas that open seawards prove episodes of major instability that have repeatedly punctuated the more typical styles of eruptive activity [Gillot et al. 1994]. Each event of this kind can be linked to proximal breccias and medial fall deposits of ash and lapilli [Ort et al. 2016].
Since the formation of the youngest caldera (Enclos Fouqué), a new intracaldera cone has grown, the summit of which has been repeatedly affected by destructive events such as small caldera and pit crater collapses [Michon et al. 2013]. The last summit collapse event occurred in 2007 (1 km wide Dolomieu caldera) and was associated with a large seaward sliding of the eastern flank (1.4 m;[Froger et al. 2015]), but did not produce a breccia deposit.
These breccias offer a great opportunity to extensively sample what is an otherwise largely inaccessible lithological variability. Here, we have studied the breccias of the Plaine des Sables and Enclos Fouqué caldera collapse events of PdF in order to characterise the lithological and petrophysical variability of the volcano and obtain new constraints on the link between magma storage, transfer, and edifice collapse. and 0.029 Ma ago, Piton des Neiges entered a more explosive phase emitting differentiated magmas with an alkaline affinity [Salvany et al. 2012, and references therein]. Intense erosion in a tropical setting has deeply incised the volcano with the formation of deep valleys and large depressions: the Cilaos, Mafate, and Salazie "cirques".
PdF represents a younger (<0.535 Ma) edifice built on the southeast slopes of Piton des Neiges [Gillot and Nativel 1989] by dominantly effusive activity [Albarède et al. 1997;Michon et al. 2013]. The edifice is built on an older and poorly exposed volcanic edifice named Les Alizés ]. Lavas cropping out at the bottom of deeply incised valleys on the flanks of the PdF edifice (Rivière des Remparts; Rivière de l'Est) are plagioclase-rich basalts, locally called "Pintades". The Pintades are significantly distinct from the subsequently emitted plagioclase-poor lavas of PdF, where olivine is the main type of phenocryst [Albarède et al. 1997]. PdF is one of the most active intraplate basaltic shield volcanoes in the world, with an average of one eruption every 10 months since 1985 [Roult et al. 2012]. Northeast and southeast rift zones delimit the steep and mobile eastern flank (Figure 1). High cliffs and deep valleys cut into the edifice and testify frequent mass-wasting events with vertical and lateral edifice failure and caldera collapses that have punctuated the eruptive activity [Gillot et al. 1994]. The oldest identified caldera collapse events are those of the Rivière de Rempart (~250 ka) and Morne Langevin calderas (~150 ka, Michon et al. [2016]). More recently, the edifice of PdF was shaped by the Plaine des Sables (65-40 ka) and Enclos Fouqué calderas (~5-2.5 ka) [Ort et al. 2016]; ). Pyroclastic ash, lapilli, and breccia layers mantle the top of the caldera cliffs and correspond to the products of highly explosive activity produced during the volcano collapses [Ort et al. 2016].
The large-scale collapses of PdF ( Figure 1) have been attributed to structural weakness derived from gravitational instability of a growing volcanic edifice (landslides with intrinsic faults) and mechanical instability (caldera collapse) following magma withdrawal during eruptions or a combination of the two [Bachèlery 1981;Duffield et al. 1982;Gillot et al. 1994;Oehler et al. 2004;Oehler et al. 2008;Merle and Lénat 2003;Merle et al. 2010]). Di Muro et al. [2014] attributed smaller summit collapses and lateral flank sliding to lateral magma transfer and magma draining from a stacked sequence of shallow sills and a crystal-rich mush pressurised by a new input of deep volatile-rich magma.

Outcrop description
Collapse-related breccias have the potential to contain lithic blocks from an extended portion of the stratig-raphy that may otherwise not be exposed and therefore provide an excellent opportunity to sample the lithological variability of a volcanic edifice. For this study, we collected unaltered blocks >10 cm in diameter with the aim of documenting the largest possible range of macroscopically distinguishable lithologies erupted during the most recent major explosive events at PdF. We sampled pyroclastic breccia deposits emplaced on the margins of the Plaine des Sables and Enclos Fouqué calderas ( Figure 2). The samples are named after their stratigraphic units and following the stratigraphic study of Ort et al. [2016]. "A" samples derive from the breccia at the Morne Langevin cliff and are related to the formation of the older Plaine des Sables caldera ( Figure 1). "B" samples are blocks from two outcrops of the younger Bellecombe breccia, covering the lower (at Piton de Bert) and upper part of the stratigraphy (outcrop 17) ( Figure 1). The two locations vary substantially in the abundance in gabbroic clasts (less abundant at Piton de Bert).
On the rims of Plaine des Sables and Enclos Fouqué calderas, the sampled units consist of beds of matrixsupported, commonly poorly sorted breccias with variable amounts of ash matrix and angular to sub-rounded clasts of all grain sizes. These layers are usually <1 m thick and drape the pre-existing topography formed by older lavas. They are covered by thin (less than a few decimeters) lapilli beds and soils. In the breccias, most blocks of vesicular lavas and dense dolerites are smaller than 40 cm in diameter, while most gabbroic and dunite clasts are smaller than 10 cm. The rarely observed dunite clasts were too small to be characterised in terms of their physical properties. Many blocks show cooling joints, cracks, layering, and variable degrees of alteration. Such clasts were discarded in order to select physically homogeneous and unaltered rock samples for laboratory testing. Component analysis of lapilli clasts (<10 mm) and larger blocks show that unaltered and mildly altered clasts represent the majority (60-70 %). Three samples served as reference samples for the selected lithologies ( Figure 2): Sample "C" is an olivine-rich porphyritic lava (locally called "oceanite") emplaced shortly before the first Bellecombe eruptive event at the site of "Petite Carriere" [Ort et al. 2016]. Sample "D" from the central part of a meter-thick, subvertical dyke cropping out in the Cirque de Cilaos (Piton des Neiges edifice) ( Figure 2). Sample "E" is an aphyric lava from the July 2015 eruption [Thivet et al. 2020].
All studied samples were described macroscopically and thin sections were prepared for most samples.

Macroscopic and microscopic sample description
Fourteen blocks were selected in order to obtain homogeneous samples that cover the largest possible range in crystal content and type, vesicularity, pore number density, and pore size ( Figure 3; Table A1). Analysis of rock cylinders drilled from the central parts of each clast allowed us to select the less altered blocks and achieve the "original" petrophysical parameters (density, porosity, P-wave velocity, and uniaxial compressive strength) ( Figure 3). Thin section observations confirm that the samples are fine-grained and have a bulk basaltic composition typical of PdF products, with variable proportions of olivine, plagioclase, clinopyroxene, and oxides. On the basis of their phenocryst assemblage, the samples fall into four different groups: an aphyric group, an olivine-rich group, a plagioclase-rich group, and a cotectic group (olivine + clinopyroxene + plagioclase) (Table A1). Rock samples are dominantly homogeneous in terms of crystal size distribution, with the exception of the sample A04, which shows some evidence of mingling between portions with variably finegrained crystal sizes. At the scale of a laboratory sample, all analyzed rocks can be considered homogeneous.
In the Plaine des Sables sample set ("A" samples), two samples are highly vesicular (A01 and A02), with pores up to 2 cm in diameter. The longest dimension of phenocrysts is small, usually <1 mm, and they are represented dominantly by plagioclase in sample A01 and olivine in sample A02. Two other samples (A03 and A04) are vesicle-poor, porphyritic, and contain a large amount of plagioclase phenocrysts up to 1 mm long, while olivine (up to 0.6 mm) and clinopyroxene (up to 0.2 mm) phenocrysts are much less abundant and smaller. Sample A03 contains small scattered mmsized "gabbroic" aggregates of plagioclase, olivine, and clinopyroxene.
In the Bellecombe sample set ("B" samples), three samples are vesicular (samples B01, B02, and B03). Sample B01 is rich in plagioclase phenocrysts (up to 1 mm), sample B02 is rich in olivine phenocrysts (up to 2 mm) and smaller plagioclases, and sample B03 is aphyric and exhibits a fine-grained matrix of randomly oriented plagioclase (<100 µm), oxides, and glass. Sample B04 shows mm-sized plagioclase and smaller olivine phenocrysts in a finer-grained matrix of olivine, plagioclase, and oxides. Sample B05 is a dolerite that exhibits spectacular 1-2 mm large variolites with a well-organised radial array of olivine and plagioclase. These two samples are classified as finegrained troctolite or micro-troctolite (olivine + plagioclase). Samples B06 and B07 are fine-grained pyroxene gabbros (most crystals are <1-0.5 mm). At the scale of the thin section, sample B06 contains thin layers rich in fine-grained (<300 µm) plagioclase that alternate with thick medium-grained bands exhibiting an ophitic texture (plagioclase laths included in mm-sized pyroxene, partly replaced by oxides). Sample B07 has a more homogeneous sugary texture. Samples B04, B05, B06, and B07 are fully representative of intrusive rocks within our sample set (Figure 4), and they likely correspond to small intrusive bodies (dykes or sills) or to marginal, rapidly cooled parts of larger intrusive bodies. The micro-troctolite and micro-gabbroic samples correspond to high and intermediate cooling rates, respectively.
The occurrence of samples rich in plagioclase phenocrysts in both the "A" and "B" sample sets was un- expected as these lithologies are almost absent in the recent eruptive products of PdF [Boivin and Bachèlery 2009] and are only known from old units cropping out at the base of Rivière des Remparts [Albarède et al. 1997;Valer et al. 2017]. This is the first time that the occurrence of these lithologies have been identified on the Enclos Fouqué rim, suggesting that they extend over a large area of the PdF edifice. Following the local nomenclature, these lavas can be called "micro-Pintades". The pre-Bellecombe oceanite sample (sample "C") contains a large amount of large olivine phenocrysts (up to 5 mm) and minor pyroxene and plagioclase. As we sampled the inner part of the lava flow, its matrix is fully crystallised and fine-grained (<100 µm long crystals). The dyke sample (sample "D") from Piton des Neiges (olivine-phyric basalt) shows sparse mm-sized euhedral olivine crystals in a coarse-grained matrix of plagioclase, olivine, pyroxene, and oxides. Its dolerite texture is very similar to that of the microtroctolites. With respect to the lava samples, the matrix of the dykes contains very little amounts of vesicles and glass and their crystal size distribution is coarser (plagioclase laths are 100-500 µm long). Our sample set thus spans the full range of basalt paragenesis described at PdF, from magmatic products typically erupted by recent eruptive activity, covering (1) aphyric, cotectic, and olivine-bearing rocks, (2) a set of relatively rare, plagioclase-rich samples corresponding to or predating the early stages of activity of PdF volcano, and (3) a set of plagioclase-rich intrusive rocks. Thin sections confirm the unaltered nature of the selected rocks. Only minor modification is observed, such as zeolite infilling some of the pores within two plagioclase-rich samples (samples A01 and B01) and secondary minerals (chlorite at the expense of olivine and possibly former interstitial tachylitic glass) in samples B04 and B05. The amount of basaltic glass is mod-erate in the 2015 reference lava sample (sample "E") ( Figure 4), but almost absent in most samples from the explosive breccias, whose matrices are often fully crystallised.

Experimental methods
Cylindrical samples, 20 mm in diameter, were cut and precision-ground such that their length was 40 mm. Before analysis, the samples were washed in distilled water and dried at 80°C for 24 h. The bulk sample density was calculated from the bulk volume and dry mass of each cylinder. The connected porosity of each sample was calculated using the bulk sample volume (calculated using the sample dimensions) and the skeletal volume provided by a pycnometer (Quantachrome Pycnometer 1200 E) at the Department for Earth and Environmental Sciences (Munich, Germany). The total porosity of each sample was calculated using the bulk sample density (calculated using the dry mass and sample dimensions) and the density of a powdered aliquot of the sample (calculated using the dry mass of the powdered aliquot and its volume, measured by the pycnometer). The isolated porosity was determined by subtracting the connected from the total porosity.The P-wave velocity of the samples was measured along the axis of each cylinder at ambient pressure and temperature, using a custom-built rig with ultrasonic wave transducers at both ends, an oscilloscope, and a pulser/receiver. P-wave velocities were measured on oven-dry samples and samples vacuum-saturated with water. P-wave velocities were measured at the Department for Earth and Environmental Sciences (Munich, Germany).
Finally, the samples were deformed in uniaxial compression at a strain rate of 10 −5 s −1 until macroscopic

Petrophysical characterisation of volcanic ejecta
Di  failure at the Institut de Physique de Globe de Strasbourg (Strasbourg, France). Axial displacement and force were measured using a linear variable differential transducer and a load cell, respectively. These measurements were converted to axial strain and stress using the sample dimensions. The displacement accumulated within the loading chain was subtracted from the measured displacement. All uniaxial deformation experiments were performed on oven-dry samples at ambient pressure and temperature. The uniaxial compressive strength (UCS) is taken as the maximum axial stress the sample withstood before succumbing to macroscopic failure. Static Young's modulus was calculated using the stress-strain data within the elastic portion of the stress-strain curve Heap et al. [e.g. 2020a]. We adopt the convention that compressive stresses and strains are positive.

Results
Aphyric lava of the 2015 eruption (sample E01) has an average total porosity of 0.28, while the average porosity is only 0.13 in the pre-Bellecombe oceanite lava (sample C01). The total porosity of our reference lava samples falls in the lower range of the variability (0.15-0.50) documented by [Di Muro et al. 2014] for 2007 lavas (aphyric to very olivine-rich lavas). As expected, porosity is very low (only 0.04) in the dyke sample from Piton des Neiges. Bulk dry density increases from 2.2 g cm −3 in modern aphyric lava, to 2.9 g cm −3 in olivine rich oceanite, to 3.0 g cm −3 in the dense olivine basalt dyke of Piton des Neiges. In the samples from the explosive breccias, total porosity ranges from 0.25 to 0.47 in the vesicular A01 and A02 samples and from 0.21 to 0.44 in the vesicular B01, B02, and B03 samples. The density of these samples ranges from 1.9 to 2.2 g cm −3 in the vesicular "A" (Plaine des Sables) samples and from 2.1 to 2.6 g cm −3 in the vesicular "B" (Bellecombe) samples. These vesicular rocks fall in or are slightly above the range of porosity and density documented in dense modern subaerial lavas of PdF [Thivet et al. 2020]. Total porosity is only as low as 0.09-0.13 in the dense (2.9 g cm −3 ) plagioclase-rich blocks of the Plaine des Sables unit (samples A03 and A04) and is 0.10-0.17 in the intrusive micro-gabbro (density = 2.7-2.6 g cm −3 ) and the olivine-rich micro-troctolites (density = 2.9-3.0 g cm −3 ) of the Bellecombe unit.
In our dataset (Table A1), the A02 lava sample represents the most porous and least dense sample, and the Piton des Neiges dyke and the B04 dolerite microtroctolite sample are the least porous and most dense samples ( Figure 5). The micro-gabbros have a very moderate total porosity of 0.10-0.17. Thin section observations show that these samples do not contain interstitial vesicular glass and are fully crystallised. Their low porosity is thus mostly related to cracks, crystal cleavages, and intergranular crystal contacts that allow

Figure :
Connected porosity as a function of total porosity for the rocks from La Réunion (data available in Table A ). The : line (i.e. where total porosity equals connected porosity) is shown on the graph. "Dykes" correspond here to D and micro-troctolite samples.
for the infiltration of fluids and there is no network of interconnected pores. The large amount of plagioclase in these rocks together with the large number of interfaces/cracks confers a slightly lower density with respect to olivine-rich dense sheeted bodies (dykes).
Connected porosity is plotted as a function of total porosity in Figure 5 (data provided in Table A1). The connected porosity of the samples we collected ranges from 0.03 to 0.40. These data show that there is little isolated porosity in vesicular samples with a total porosity above 0.25. This is consistent with results obtained on recent lavas and scoria of modern activity that contain coalesced pores [Thivet et al. 2020]. Dense samples with a total porosity below 0.25 (e.g. blocks A02, A03, B02, B04, B05, and D01) contain significant isolated porosity (average = 0.03--0.08). Samples from micro-troctolite block B05 are characterised by the highest values of isolated porosity (up to 0.10). In this block, thin sections show the occurrence of secondary mineral replacing former interstitial glass and possibly filling some pores. Some sample heterogeneity in terms of crystal content and phase proportions together with variable amounts of cracks and small isolated pores produce minor data scattering for these dense and crystal-rich samples.
Dry (black symbols) and wet (blue symbols) P-wave velocity is plotted as a function of connected porosity in Figure 6a and dry P-wave velocity is plotted as a function of density in Figure 6B (data provided in Table A1). The dry P-wave velocity of the samples ranges from~1 to~6 km s −1 . In general, P-wave velocity decreases as a function of increasing porosity ( Figure 6A) and increases as a function of increasing bulk sample density ( Figure 6B). The highest P-wave velocity values are thus measured in micro-troctolites and Piton des Volcanica 4(1): 41 -66. doi: 1 .3 9 9/vol. 4. 1.4166 Neiges dyke and the lowest in the vesicular lavas. However, the link between porosity/density and P-wave velocity is complex ( Figure 6A and Figure 6B). Among the dense or low-porosity samples, we note that the microgabbro samples have a very low P-wave velocity of only 1.3-2.1 km s −1 , much lower than other samples of the similar porosity. Scattering exists also in the vesicular rocks, with aphyric sample E01, plagioclase-phyric sample B01, and olivine basalt sample A02 showing a low dry P-wave velocity and B03 aphyric sample and A01 plagioclase-phyric sample having higher values. The dense olivine-rich oceanite samples plot on the best-fit linear trend of the data (the dashed line in Figure 6A). The plot of wet-to-dry P-wave velocity ratio as a function of connected porosity ( Figure 6C) shows that saturation with water does not significantly influence the P-wave velocity of most dense samples with a porosity below 0.25 and of some vesicular samples (A01 and B03). Wet velocity could not be measured in the dense micro-gabbros, but saturation does not affect the P-wave velocity of the dense micro-troctolite or the Piton des Neiges dyke. In the vesicular lavas, P-wave velocity in most samples with a porosity above 0.25 is higher when saturated with water ( Figure 6C). The Pwave velocity almost doubles (average wet/dry ratio: 1.91) in the highly vesicular A02 lava sample. However, the very small change observed in some vesicular samples (A01 and B03) suggest that besides bulk vesicularity, bubble size and spatial distribution together with connectivity and partial infilling by secondary minerals (A01) control the rock permeability and water saturation.
Representative stress-strain curves are shown in Figure 7, and are typical of those for volcanic rock deforming in compression (described in, for example, [Heap et al. 2014a]). Axial stress is first a non-linearly increasing function of axial strain, often attributed to the closure of pre-existing microcracks. This is followed by a quasilinear elastic stage. Axial stress is then non-linearly increasing function of axial strain, attributed to the formation and growth of microcracks (i.e. inelastic deformation). A peak stress can then be observed, which is followed by a stress drop. It is during this stress drop that microcracks coalesce to form a macroscopic fracture. Figure 8 and Figure 9 show plots of uniaxial compressive stress (UCS) and Young's modulus, respectively, as a function of connected porosity. The UCS and Young's modulus of the samples collected ranges from~15 to~390 MPa and~2 to~70 GPa, respectively.
The average strength and Young's modulus are highest in the dense Piton des Neiges dyke and PdF microtroctolite samples, and are lowest in the vesicular lavas ( Figure 8 and Figure 9; Table A1). However, and despite their low porosity, the strength and Young's modulus is also very low in the granular micro-gabbros. These data show that UCS and Young's modulus decrease non-linearly as a function of increasing poros- Dry (black symbols) and wet (blue symbols) P-wave velocity as a function of connected porosity for the basaltic rocks from La Réunion collected for this study (data available in Table A ). The dashed line is the best-fit line of the dry P-wave velocity data.
[B] Dry P-wave velocity as a function of dry bulk sample density.
[C] The ratio of wet to dry P-wave velocity as a function of connected porosity.
ity: decreases are large when porosity is increased in dense rocks below a porosity of 0.25, and small in vesic-

Petrophysical characterisation of volcanic ejecta
Di  ular rocks (lavas) when porosity is increased above a porosity of 0.25. In the group of vesicular lavas, the lowest values (<10 GPa) are measured in slightly altered samples with some zeolite infilling in the vesicles (B01 and one measurement from A01 sample). Among the dense samples, the slightly altered chlorite bearing samples (B04 and B05) do not differ from the unaltered ones. The micro-gabbros represent a notable exception because, in spite of their high density and low porosity, they have low values of UCS and Young's modulus when compared to other samples of a similar porosity (Figure 8 and Figure 9; Table A1).   Table A ).

Figure :
Young's modulus as a function of total porosity for the basaltic rocks from La Réunion collected for this study (data available in Table A ).

Volcano heterogeneity constrained by petrophysical data
Our data show that the P-wave velocity ( Figure 6; Table A1), UCS ( Figure 8; Table A1), and Young's modulus ( Figure 9; Table A1) of basaltic rocks from PdF decrease as a function of increasing porosity, and that Pwave velocity can increase when saturated with water ( Figure 6;  (Figure 8). Thin section observations show that the gabbro samples are unaltered, and so we conclude that their granular texture and the occurrence of a fine network of microcracks and cleavage planes are the most likely causes of their relative weakness and low P-wave velocity. The Young's modulus of volcanic rocks is reviewed in Heap et al. [2020b]. These authors show that the Young's modulus of a wide variety of volcanic rocks (andesite, dacite, basalt, pyroclas- Mt Etna (Zhu et al., 2016) Pacaya (Schaefer et al., 2015) Kīlauea (Bubeck et al., 2017) Harat Almadinah (Al-Harthi et al., 1999) Stromboli (Apuani et al., 2005)  ocean basalts (Christensen et al., 1980) Deccan basalts (Vedanti et al., 2015) Harat Almadinah (Al-Harthi et al., 1999) Paraná-Etendeka (  tic rocks) decreases as a function of increasing porosity. The number density of microstructural defects (microcracks, pores) increases as porosity increases, which serves to reduce material stiffness (Young's modulus). The scatter in the UCS (Figure 8) and Young's modulus ( Figure 9) data reflect variations in the proportions of pores and microcracks (porosity is a simply a scalar and does not account for such geometric differences) and variations in crystal content, both of which are known to influence the strength and stiffness of volcanic materials [e.g. Heap et al. 2016;Zorn et al. 2018;Heap et al. 2020b]. However, we highlight here that the type of phenocrysts does not seem to play a major role, as the Young's modulus of a very olivine-rich oceanite lava is the same as the plagioclase-bearing or aphyric lavas (Table A1). Slightly altered vesicular lavas records the lowest Young's modulus (<10 GPa) in the high porosity rocks. We compare our new P-wave velocity, UCS, and Young's modulus data for PdF with previously published data for basalts in Figure 10. The published values for basalt are in general agreement with our new data for basaltic rocks from PdF ( Figure 10). Of note, the micro-gabbro samples have a P-wave velocity much lower than other samples of the same porosity and do not follow the general trend of decreasing P-wave velocity as a function of increasing porosity defined by the compiled data ( Figure 10). A low P-wave velocity for low-porosity rocks is likely the consequence of a high microcrack number density, and therefore low pore number density. Microcracks have a low volume compared to pores, due to their high aspect ratio, and therefore do not contribute significantly to the porosity of a given sample. However, microcracks can greatly decrease elastic wave velocities [e.g. Vinciguerra et al. 2005;Stanchits et al. 2006]. P-wave velocities in gabbros can be enhanced upon exposure to high-temperature, as experimentally shown by [Keshavarz et al. 2010]. Partial replacement of pyroxene by oxides in our samples could record re-heating of the micro-gabbros after their cooling and solidification. We also note that the UCS ( Figure 10B) and Young's modulus ( Figure 10C) of the low-porosity basalts from PdF are, in general, higher than those previously measured for basalts. Values of UCS and Young's modulus for basaltic rocks from PdF can reach~400 MPa and 70 GPa, respectively ( Figure 10B and C; Table A1). However, this likely reflects a sampling bias towards porous samples and/or samples that have not been deeply buried, rather than that basalts from PdF are stronger and stiffer than the rocks that comprise other basaltic volcanoes, such as Pacaya volcano (Guatemala), Mt Etna (Italy), and Kīlauea.
We highlight that our mechanical data (Table A1 and  Table 1) were collected at room pressure and temperature, and at a strain rate of 10 −5 s −1 . First, it is well established that the compressive strength [e.g. Kennedy and Russell 2012;Heap et al. 2015;Zhu et al. 2016;Ryan et al. 2020] and Young's modulus [e.g. Heap et al. 2020b, and references therein] of volcanic rock will increase as a function of increasing effective pressure (i.e. depth). Second, even in the absence of glass, an increase in temperature can influence the strength and Young's modulus of volcanic rock [e.g. Benson et al. 2012;Bakker et al. 2016;Coats et al. 2018;Heap et al. 2018;Heap et al. 2020b]. For example, Heap et al. [2018] show that the UCS and Young's modulus of andesite from Volcán de Colima (Mexico) both increased as temperature increased from 25 to 700°C. Third, although the influence of strain rate on the mechanical behavior of rocks in the brittle regime is small compared to its influence on viscously deforming rocks and magmas, an increase in strain rate from 10 −6 to 10 −1 s −1 increased the UCS and Young's modulus of andesite from Volcán de Colima from 36.6 to 60.5 MPa and from 12.0 to 19.0 GPa, respectively [Lavallée et al. 2019]. We also note that an increase in effective pressure (i.e. depth) will likely increase the P-wave velocity of the studied volcanic rocks [e.g. Vinciguerra et al. 2005;Benson et al. 2006;Stanchits et al. 2006;Nara et al. 2011]. Nevertheless, the methods used herein are a standard way to assess these physical and mechanical properties and, as shown in Figure 10, we benefit from the ease of comparison with pre-existing datasets. However, and also due to the scale-dependence of these physical properties, we urge caution when directly comparing these data with those collected using geophysical field methods.

Implications for volcano stability
The morphology of PdF volcano results from a combination of construction, deformation, and destruction processes that occurred during its lifetime (~0.5 Ma). Submarine flanks have gentle slopes ranging between 2 and 5°on average and are covered with products related to the erosion and collapse of the subaerial part of the edifice , and references therein]. Slope angles increase near the coast to 8-15°a nd peak higher on the volcano's flanks (20-35°) (Figure 11) . The summit area, made up primarily by the plateau of the Enclos Fouqué caldera, exhibits a low slope inclination of 2-8°. This morphology is not observed on the western flank where the edifice is buttressed by the older and larger Piton des Neiges edifice. The Enclos Fouqué caldera hosts a steep central cone frequently affected by summit collapses of calderas and pit craters. Michon et al. [2016] review the models proposed to explain the formation of the U-shaped Enclos Fouqué caldera and the potential link between central calderas and the collapse of the steep eastern flank. They highlight the possible role of dense intrusive bodies in forming gliding horizons [Merle and Lénat 2003]. Di Muro et al. [2014] linked the fluidisation and lateral transfer of crystal-rich dense cumulates to the seawards sliding of the eastern flank during large eruptions. Volcano monitoring shows that the steep eastern flank of PdF slides at a constant low rate (<3 cm yr −1 ), which can accelerate during the phases of deep magma ascent that precede eruptions [Peltier et al. 2016]. Meter-scale sliding events of limited swaths of the eastern flank have been frequently recorded during lateral eruptions, such as the 2007 caldera-forming distal event [Froger et al. 2015], and have been linked with lateral magma injection below the steep eastern flank [Got et al. 2013;Di Muro et al. 2014;Albert et al. 2019] (Figure 11). Large vertical caldera collapses, often associated with highly explosive events and the lateral sliding of the eastern flank, have punctuated the activity of this dominantly effusive volcano [Ort et al. 2016]. During the last large caldera and flank collapse events (Enclos Fouqué-Bellecombe events), lithic clasts of variably evolved melts and gabbroic intrusive bodies-representing the whole transcrustal plumbing system-were erupted [Upton et al. 2000;Ort et al. 2016]. In these breccias, gabbros (plagioclase + clinopyroxene + olivine) are the main type of intrusive rock while dense dunite (olivine) and wherlite (olivine + clinopyroxene) recording lithostatic pressures corresponding to depths >5 km are rare. These intrusive rocks were likely stripped from the old hypovolcanic intrusive complex identified by geophysical studies below the western part of the Enclos Fouqué caldeira and the Plaine des Sables area (PSc in Figure 11; from [Gailler et al. 2009]). [Gillot et al. 1994] suggested that La Réunion Island has been periodically affected by repeated par-  Gillot et al. [1994] speculated that the repeated intrusion of magma bodies and dykes inside the volcanic edifices might have contributed to increase their disequilibrium. In the model of Gillot et al. [1994], the successive dyke intrusions of magma create the arcuate NE-SE rifting of PdF and contribute to the weakening of the steep volcano flank on the eastern side of the volcano. The part of the magma ascending at crustal level that is not erupted develops shallow crustal (micro)gabbroic intrusions and intrusive sheets, possibly at the interface between the older volcanic structures and the new volcanic pile/cone built inside the preceding caldera. Old mafic-ultramafic intrusive sheets have been drilled on the western volcano flank of PdF and they underlie beneath a 1 km-thick sequence of mostly subaerial lavas Rançon et al. 1989] (Alc in Figure 11. Hyaloclastite beds have been identified in the range 400-800 m b.s.l. in the eastern flank drilling and the deepest bed separate the vesicular lava pile from the underlying intrusive body  . Recent petrological studies confirm that magmas injected in the crustal section of the volcano plumbing system experience extensive gabbroic fractionation, while mafic dunite and wherlite fractionation trends are dominant at deeper levels (mantle/crust underplating and deeper) [Boudoire et al. 2019].
In conclusion, several studies suggest that, both at Piton des Neiges and PdF, the top of gabbroic units underlying the thick pile of vesicular lavas, facilitates detachment, by localising magma injection in sills, hydrothermal circulation, and shear deformation.[Famin and Michon 2010] describe a detachment fault in outcrops of the northern part of Piton des Neiges, focusing both ductile and brittle deformation with a seawards orientation. Petrological and structural analysis of the PdN outcrops show that the detachment is located between the top of an old body of layered gabbro and a 200 m-thick cover of zeolitised debris-avalanche breccia. The upper 10 m of the gabbro are characterised by hydrothermal alteration in the greenschist facies and an intense deformation. A~50 m-thick pile of tens of 0.2-1.5 m-wide olivine-rich sills, with variable degrees of greenschist facies alteration and deformation, occurs between the breccia and the underlying gabbro. A large number of basaltic dykes crosscut the layering of the gabbro, while only a few dykes are observed in the overlying units. The very low abundance of basaltic dykes in the breccia and overlying units suggest that frequent magma intrusions occurred before one or more flank collapses had occurred between 2 and 1 Ma ago. The PdN gabbro, the sheeted sills, and the breccia are today covered by~2.3 kmthick sequence of young (<0.19 Ma) differentiated lavas and tephra from the late phase of activity at Piton des Neiges volcano. All these results suggest that repeated injection of basaltic magma in sills can drive fast, cointrusive slip, while hydrothermal alteration is only localised and linked to multiple events of sill emplacement and cooling and can promote slow inter-intrusion creep [Famin et al. 2016]. This model is consistent with the relatively low amount of altered rocks (<30-40 %) we found in the studied pyroclastic breccias. Famin and Michon [2010] suggest that low-angle intrusion of olivine-rich magmas along a plane of weakness can explain the frequently observed, syn-eruptive, seaward deformation of the unbuttressed eastern flank of PdF. [Chaput et al. 2017] conclude that, in contrast to Hawaiian volcanoes, La Réunion volcanoes are brittle edifices built on a weak substratum and that their long-term deformation cannot be attributed to the deformation of ductile bodies within the edifices. Nevertheless, the several hyaloclastite units occurring between the lava pile and the underlying intrusive body can represent a potential ductile layer facilitating the seawards sliding of the PdF volcano flank, as suggested on other basaltic volcanoes [Schiffman et al. 2006] (Figure 11). Piton de la Fournaise is a relatively thin edifice build on the older Les Alizé volcano and soil or breccia units might occur between them and represent weak units [Del Potro and Hürlimann 2008].
The models of [Gillot et al. 1994] and [Famin and Michon 2010] both propose that multiple dyke intrusions play an important role in weakening the volcano edifice. Famin and Michon [2010] stress the role Petrophysical characterisation of volcanic ejecta Di  of olivine-rich magmas in destabilising the steep eastern flank. Despite their differences, all these models require the occurrence of ductile layers (e.g. sediments; magma injected in sills; hydrothermal alteration) to enhance volcano deformation and destabilisation. Our data show that the petrological composition of lavas (olivine-or plagioclase-bearing, phenocrystrich or phenocryst-poor) affects rock density, especially when the samples contain abundant large phenocrysts, but that it does not play an important role in controlling their stiffness, for which porosity is the most important parameter. As shown for recent eruptions at PdF, the degree of pre-eruptive degassing [Di Muro et al. 2014;Gurioli et al. 2018] and the magma ascent rate [Thivet et al. 2020] modulate the vesicularity and the density of the eruptive products. Long-lasting, shallow pre-eruptive storage and cooling, and slow ascent rate results in open system degassing and densification of the magma [Di . Plastic deformation of voluminous dense olivine cumulates at the base of the volcanic edifice has been invoked to explain flank deformation at Hawaiian volcanoes [Clague and Denlinger 1994]. At PdF, large-scale ductile deformation has also been invoked to explain the deformation of the eastern flank but it has been attributed to hydrothermally altered rock units and/or pore fluid pressurisation [Merle and Lénat 2003]; [Michon and Saint-Ange 2008]. Indeed, a ductile failure mode in volcanic rocks can be promoted by hydrothermal alteration [Mordensky et al. 2019]. However, hydrothermal alteration can also decrease permeability and increase pore fluid pressure ]. An increase in pore fluid pressure can also serve to embrittle rocks by decreasing the effective pressure acting on the rock [Farquharson et al. 2016].
At PdF, the olivine-rich oceanite lavas formed by the remobilisation of olivine-rich cumulates can be treated as a physical analogue of the dunitic bodies at Kīlauea. However, drilling and geophysical data suggest that their volume at shallow level is likely to be small. Petrological studies of both PdN and PdF show that large intrusive bodies are often gabbroic and not dunitic and that ultramafic (dunite-wherlite) units occur in the lower part of the hypovolcanic bodies (Figure 11). This is a consequence of the evolved average composition of the magmas injected at shallow level in the crustal section of La Réunion, which is thinner than at Hawaiian volcanoes [Boudoire et al. 2019]. Our experimental dataset shows that gabbros are comparatively weak because of their high density of microfractures and cleavage of the primary minerals. Episodes of thermal stressing can further weaken these granular rocks. We found the largest difference in rock stiffness between the dense fine-grained plagioclase (micro-pintade), olivine basalts, and olivineplagioclase (troctolite) on one side, and the vesicular lavas, on the other side. The dense rocks can be emplaced in subvertical dykes or in subhorizontal sills and have experienced high-intermediate rates of cooling and extensive degassing. As postulated by Famin and Michon [2010], the interface between a thick lava pile and a sequence of sills and dykes represents a major contrast in rock properties (stiffness, density, and strength), which can favour further magma arrest, the circulation of hydrothermal fluids, and rock alteration as well as the brittle destabilisation of the lava rock pile. We add that vesicularity is a major controlling parameter and that the (micro)gabbroic intrusive bodies represent weak bodies and can favour the formation of detachment surfaces. Most importantly, repeated reheating of the gabbros by multiple sill injections (Figure 11) can substantially weaken these rocks [e.g. Keshavarz et al. 2010], irrespective of occurrence of hydrothermal alteration.
Our study focuses on unaltered intrusive and effusive rocks in order to highlight the primary variability of rock properties. A study on the effect of hydrothermal alteration on rock physical properties would require a constraint on the spatial scale of this process. Pervasive hydrothermal alteration can potentially permeate both cooling magma bodies and host rocks and might even reduce the contrast in rock properties. However, as stressed by Famin et al. [2016], fluid mobility and rock alteration are focused near the fault/detachment planes and do not produce large-scale rock property modifications. As the sill volumes are relatively small and therefore the amount of related fluids and heat are very limited, we have to conclude that the faults identified at Piton des Neiges or hypothesised at PdF have been active for a long time and act as semi-permanent drainages of fluids and heat. In this scenario, meteoric fluids represent a major source of fluids with respect to magmatic fluids. Our results therefore allow us to constrain the large-scale heterogeneity of the volcanic rocks composing the volcanic edifice at La Réunion, with the exception of localised modifications related to hydrothermal alteration along long-lasting active faults. A corollary of this reasoning is that the very existence of the faults can stem from primary rock heterogeneities forced by the load of the huge volcanic pile. All vesicular lavas we studied are characterised by a UCS of 30-120 MPa (average 57 ± 31 MPa), a stress range that overlaps with the measured entrapment pressure of most melt and fluid inclusions (average 30 MPa; ) hosted in olivine crystals of PdF magmas. As soon as PdF magma pressure overcomes this (lithostatic) threshold range, fractures can easily propagate into the lava pile, promoting dyke injection and magma ascent. This pressure threshold is also consistent with the shallow location of the plumbing system collapsed to form the 2007 summit caldera, following a substantial drainage of magma during the eruption [e.g. Peltier et al. 2009;Di Muro et al. 2014]. Deeper reservoirs, located within the oceanic crust and deeper [Boudoire et al. 2019], are likely overlain by a thick sequence of dense stiff sills, as ob-Volcanica 4(1): 41 -66. doi: 1 .3 9 9/vol. 4. 1.4166 served at PdN. These deep reservoirs can elastically deform during the periodical refilling and draining without collapsing, as large lithostatic pressures (Table A1) are needed to break these stiff finely-crystallised rocks. The occurrence of very frequent and small summit collapses (pit craters, <1 km summit calderas) at PdF [Michon et al. 2013] is consistent with the drainage of small volumes of evolved magmas stored at shallow depth, i.e. inside the weak and easily breakable vesicular lava pile (Figure 11). Gallart et al. [1999] obtained a large-scale P-wave velocity model along a SW-NE transect crossing the La Réunion edifice. In this model, velocities in excess of 5 km s −1 are associated with the oceanic (basaltic crust) basement and the underlying underplating of dense olivine-clinopyroxene-rich (dunite-wherlite) rocks. Pwave velocities in the dense olivine-rich underplating layer are in the range 7.4-7.9 km s −1 , well above the range we document for the breccias of Bellecombe and Plaine des Sables (5.9-1.3 km s −1 ; Table A1). In the volcanic edifice, Gallart et al. [1999] identify clear heterogeneities with high-velocity (4.2-6.3 km s −1 ) material attributed to an intrusive body surrounded by lowvelocity (<4 km s −1 ) lavas and pyroclasts that form most of the volume of the edifice. Our P-wave velocity data for the vesicular PdF lavas (Table A1 and Table 1) are in agreement with seismic velocities inferred for the volcanic pile forming La Réunion island. Moreover, Pwave velocities we measure in dense plagioclase-rich rocks (micro-gabbros, dyke, micro-troctolite; Table A1 and Table 1) are also in the range of the dense core of the volcanic edifice identified by geophysical studies ( Figure 11).

Implications for inversion of geophysical data
Geological and geophysical studies have provided evidence for a large gabbroic intrusive complex at the core of PdN volcano, the main edifice that forms the subaerial part of La Réunion island. Gailler and Lénat [2012] suggested that this gabbroic core can extend to the base of the La Réunion edifice at~4 km below sea level (b.s.l.). A 3 km-deep drillhole near the subaerial base of the eastern flank of PdF provided evidence for a 1 km-thick lava pile, underlain by 1.5 km of fine-medium grained gabbro (average crystal size 1-3 mm), 0.3 km of alternating wherlite and dunite, and 0.2 km of alternating dunite and micro-gabbro ]. This dense and dominantly gabbroic core has been attributed to an old volcano (Les Alizès volcano; Alc in ( Figure 11) built on the SE flank of Piton des Neiges before the birth of PdF [Gailler et al. 2009;Gailler and Lénat 2012]. Gravimetric [Gailler et al. 2009] and seismic tomographies ( [Prôno et al. 2009]; [Nercessian et al. 1996]) have provided evidence for the occurrence of small dense intrusive complexes below PdF summit cone. Prôno et al. [2009] showed that Pwave velocity increases with depth from~3.2 km s −1 at the summit of PdF to~5.3 km s −1 at 6 km b.s.l. This range is in agreement with our experimental P-wave velocity data on lavas, dykes, sills, and gabbroic intrusive blocks scattered in the explosive breccias (Table A1). Gailler et al. [2009], Prôno et al. [2009], and Nercessian et al. [1996] identified a small high-velocity volume (4 km s −1 ) below the summit craters, between the topographical surface and 1 km above sea level (Figure 11). We propose that this high-velocity volume consists of a dense plug of dykes and sills, surrounded by lavas with a lower velocity (3.2 km s −1 ). These authors also identified low-velocity zones inside the highvelocity volume, located above and below sea level and attributed them to possible zones of magma/fluid storage ( Figure 11).
Our experimental data fully support the interpretation of vesicular lithologies surrounding a denser core. The dense core can be formed by gabbroic/doleritic sills and dykes with variable amounts of olivine. Our data for gabbroic rocks show that these dense rocks can have low P-wave velocity ( Figure 6A; Table A1), likely because of their high density of microcracks, cleavages and their specific intergranular locking. The zones with low seismic velocity could thus represent areas with a gabbroic composition, linked to persistent magma storage. However, the thermal cycles experienced by these holocrystalline rocks will dramatically affect their physical properties. This interpretation can thus solve the apparent inconsistency between the large magma storage volumes (>1 km 3 ; [Prôno et al. 2009]) inferred on the basis of seismic tomography and the widespread seismicity recorded in the same volumes by [Massin et al. 2011]. We thus propose that these low-velocity zones correspond to dominantly gabbroic bodies with variable amounts of other lithologies (dykes) and small volumes of fluids (melts; hydrothermal fluids). Older and more dense high-velocity bodies identified on the western volcano flank correspond to the old hypovolcanic complex related to the pre-Plaine des Sables volcanic edifice (Psc in Figure 11). The occurrence of "gabbroic" cores (i.e. plagioclaserich) at La Réunion clearly contrasts with results obtained on the island of Hawai'i, where the zones of high seismic velocity that coincide with high gravity anomalies are interpreted as dense dunite (i.e. olivinerich) bodies within or below magma storage volumes [Clague and Denlinger 1994]. Ductile deformation of the dunite bodies has been invoked to help transfer stress to the distal flank and mobilise the seawardssliding flank of Kīlauea. The difference between the two systems is consistent with the average evolved composition of magmas at La Réunion [Albarède et al. 1997;Boudoire et al. 2019] with respect to Hawai'i and the smaller size of the volcanic system, resulting from a lower magma output rate. Identification of specific lithologies on the basis of geophysical data is a complex task. Our data show that, in a P-wave velocitydensity space, a poorly-vesicular olivine-rich lava body Petrophysical characterisation of volcanic ejecta Di  Figure : Synthetic structure of PdF volcano. : deep magma ascent path feeding the shallow magma reservoirs below the summit Dolomieu caldeira. : dyke injection provoking the compression ( ) and sliding of the vesicular lava pile forming the eastern volcano flank [Got et al. ]; and : possible sill injection paths below the eastern volcano flank feeding distal eruptions like that of [Di Muro et al. ]. cannot be discriminated from a rapidly cooled dense microcrystalline intrusive body (emplaced in a subhorizontal sill or sub-vertical dyke) and an olivine-rich troctolite dyke/sill cannot be discriminated from a finegrained plagioclase-rich dyke/sill. In vesicular lavas, their average cooling rate (which decreases with increasing lava flow thickness) and burial pressure can result in large scatter in P-wave velocity as a function of either density or porosity. Moreover, micro-fracturing plays an important role as nominally dense rock such as gabbro can have a low P-wave velocity ( Figure 6A; Table A1).
A detailed interpretation of spatial distribution of lithologies must also take into account the role of fluid saturation on seismic velocities. We show that P-wave velocity of dense dykes and fine-grained intrusives (troctolites) is not affected by fluid saturation (Figure 6C; Table A1). On the contrary, a strong velocity increase is observed in vesicular lavas ( Figure 6C; Table A1). Routine volcano monitoring of volcanic activity on the basis of the temporal change in seismic velocity [e.g. Brenguier et al. 2008] is partly hindered by the effect of fluid circulation (hydrothermal or surface/water table fluids) in cracks and fractures [Rivet et al. 2015;Savage et al. 2015]. Our data help to quantify the effect of an increase in the height of the water table due to volcano pressurisation and inflation or a seasonal change of the water table depth. This new petrophysical dataset, summarised in Table 1, can thus be used to filter these secondary effects in order to identify magmatically and tectonically driven processes.

Conclusion
Basaltic shield volcanoes are large edifices with gentle slopes, often <10°and a relatively homogeneous bulk chemical composition. Our new petrophysical data show that, in spite of the little changing bulk composition, such volcanoes can be characterised by rocks with highly variable rock physical and mechanical properties, a factor that could play a major role in their destabilisation. This is the first time that a wide variety of samples from PdF collected at different localities have been characterised in terms of their physical and mechanical properties. The laboratory results are considered representative of the natural variability observed at PdF and can be applied to the basaltic edifices. The data show that basaltic rocks can vary considerably in terms of their density, porosity, P-wave velocity, and uniaxial compressive strength. Because individual rock types are not characterised by specific physical property values, and that properties are highly dependent on thermal history and microcracking, it is challenging to use physical property measurements to identify lithology. In addition, it can also be concluded from our data that variable physical properties do not necessarily require hydrothermal alteration, and that thermal cycles can produce a spectacular change in the properties of dense granular materials (i.e. the gabbros). The results show that even unaltered rock can be weak enough to produce instabilities, which can then affect the morphological evolution of basaltic edifices. Therefore, the large petrophysical heterogeneity detailed here for rocks comprising the edifice at PdF should be considered when interpreting/modeling basaltic edifices with geophysical methods. Our new petrophysical database for intrusive and surface volcanic rock lithologies of PdF should assist the interpretation of geophysical measurements collected for volcano monitoring and eruption forecast.

Author contributions
ADM and UK conceptualised the study. ADM and FS performed field work and petrological analysis. FS, UK, and MH performed laboratory analysis. All authors discussed the data and revised the manuscript.

Acknowledgements
This project was funded by ANR national French project "SlideVOLC" (https://anr.fr/ Project-ANR-16-CE 4-4, 18.07.2019). ADM was supported by a CAS (Center for Advanced Studies) fellowship for his work at LMU laboratory. MJH thanks Thierry Reuschlé. DBD acknowledges the support of ERC 2018 ADV Grant 834225 (EAVESDROP). We thank L. M. Rossetti and F. Di Traglia whose insightful reviews significantly improved the clarity of the manuscript and J. Browning for handling the review process.

Data availability
Rock samples, rock cores and thin sections are stored at Piton de la Fournaise volcanological observatory (OVPF/IPGP) and available upon request to the corresponding author; Raw data of petrophysical measurements summarised in appendix are available upon request to the corresponding author. The raw data from the uniaxial compressive strength (UCS) experiments performed at Université de Strasbourg (France) can be downloaded at: https://doi.org/1 .6 84/m9. figshare.14363276.v1.

Copyright notice
© The Author(s) 2021. This article is distributed under the terms of the Creative Commons Attribution 4.0 International License, which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made.  A Appendix 1