Abstract

Slow-moving clayey earthslides frequently exhibit seasonal activity suggesting that deformation processes do not only depend on slope and intrinsic geomechanical parameters. On the contrary, seasonal motion patterns are frequently observed with acceleration during the wet season and deceleration during the dry season. Within landslides, it is mainly the phreatic water table which is monitored. However, in the case of deep-seated landslides made of heterogeneous lithological units and with several slip surfaces, the characterization of the phreatic water table does not allow to relate satisfactorily the activity of the landslide with environmental parameters such as rainfall and subsequent water infiltration at depth. This paper presents a seasonal analysis of water infiltration within a slow-moving clayey landslide. Results of an extensive geotechnical and geophysical prospect are first exposed. Then, rainfall and water table level time series are analysed for two water tables using the cross-correlation technique: the phreatic water table located a few metres deep and a water table located above a shear surface located 12 m deep. Results show that water infiltrates faster down to the deepest water table. Then, time series were split between “dry” and “wet” seasons and the effective rainfall was computed from the original rainfall time series. Cross-correlation results show that the phreatic water table responds identically to rainfall in both seasons. On the contrary, the water table located above the shear surface has a very contrasting behaviour between summer (mainly drainage) and winter (behaviour similar to the phreatic water table with storage of water during a few weeks). This difference in behaviour is in agreement with the landslide kinematics.

1. Introduction

Slow-moving clayey earthslides frequently exhibit seasonal activity suggesting that deformation processes do not only depend on slope and intrinsic geomechanical parameters. On the contrary, seasonal motion patterns are frequently observed with acceleration during the wet season and deceleration during the dry season [17]. On the one hand, acceleration phases are classically related to an increase in pore water pressure following water infiltration and an increase in water table level during wet seasons. On the other hand, deceleration phases are related to water drainage, lowering of the water table level, and subsequent decrease in pore water pressure during dry seasons. This might also be related to potential dilation of fine-grained material in the vicinity of the shear zone [8, 9], or to liquefaction/fluidization of water-saturated sediments close to the rupture zone [1013]. However, this seasonal dichotomy is not systematic and slow-moving landslides can be reactivated and fail without warning during dry seasons [14]. Water flux is mainly considered vertical from the surface down to depth with infiltration occurring through the matrix [15]. As such, 1D models for water infiltration are widely used to evaluate and/or predict shallow landslide response to rainfall. This led to an establishment of empirical relationships relating the total amount of rainfall over a period (thus defining a so-called intensity-duration index) to shallow landside activity, at a local, regional, and/or global scale (e.g., [1621]). Based on rainfall thresholds, these relationships were further applied to the build-up of early-warning systems [22, 23] and/or to landslide forecasting [24, 25]. However, these approaches are devoted to shallow landslides, with disrupted mass a very few metres thick and with a relatively homogeneous motion for water infiltration. They do not take into account possible and rapid infiltration through fissures. Furthermore, they cannot account for complex deep-seated landslides with possibly several shear surfaces located down to several tens of metres. The evolution of ground mechanical characteristics (and the physical processes at stake) close to rupture surfaces which control landslide acceleration and/or fluidization is still an open question, and several works were devoted to this specific task in the past decade [1013, 2630]. However, and even if rainfall and subsequent water infiltration are known to be a major driving mechanism for landslide acceleration (e.g., [31]), the experimental observation of water infiltration deeper than a few metres remains poorly reported.

In the Trièves area (French Western Alps), numerous slow-moving clayey earthslides are observed and have been studied for more than 30 years [2, 3, 32, 33]. In this area, conceptual models based on field observations and tests proposed that water infiltrates preferentially through subvertical surface cracks which are further connected to shear surfaces at depth [5, 34]. This model was confirmed both numerically [35] and experimentally [36]. From a numerical point of view, it remains difficult and challenging to integrate dense fissure networks to hydrological models [37]. From an experimental point of view, multidisciplinary approaches are often required to evidence infiltration phenomena. Recently, geophysical (resistivity and seismic) methods along with hydrological and geological prospecting were applied on the Avignonet landslide in the Trièves area [36]. Results showed that fissures served as preferential paths for rapid water infiltration from the surface down to the phreatic water table located a very few metres below the surface. Similar results were obtained by Travelletti et al. [38] during a controlled rain infiltration experiment monitored by electrical resistivity within black shales (Super-Sauze landslide, SE France). However, because of superficial water tables and limited resolution at depth, these observations were limited to very shallow infiltration. Using temperature logging in a borehole naturally filled with water, Bièvre et al. [36] noticed a small positive temperature anomaly at 10 m depth in summer which they related to a fast infiltration of water along a fissure/shear surface. Finally, most of water table datasets within deep-seated landslides of a thickness of several tens of metres concern phreatic (i.e., superficial) water tables (among many others, [6, 39]). All these results indicate that it remains difficult to evidence experimentally in situ fast motions of water infiltration deeper than a very few metres. To overcome this issue, an original experiment was designed in this work. Two piezometers (equipped with pore pressure sensors) were installed at different depths: the first one to monitor the phreatic water table located at a few metres depth and the second one to monitor a water table located a few metres deeper, immediately above a shear surface located at a depth of 12 m. Water tables and rainfall were monitored for 4 years from late 2012 to late 2016. The first objective of this work is to present the results of a multidisciplinary study to document a slow-moving clayey landslide. The second objective is to expose original results which highlight experimentally the differential water infiltration modality for the two water tables and, more specifically, to evidence the role of fissures as preferential and fast water infiltration paths at depth. The third objective of this work is to evaluate the use of effective rainfall (instead of total rainfall) to establish correlations with water table time series and to test the possible advantage of this dataset corrected from evapotranspiration. Finally, the final objective of this work is to establish seasonal correlations between time series and to evidence a possible variability of water infiltration from the surface and further underground water flow as a function of the period of the year.

2. Study Site

The study site is located in the western Alps in SE France (Figure 1(a)), near the city of La Mure. It is located on the eastern flank of the NNE-orientated valley of La Bonne river, the basement of which is made of marly limestones of Aalenian age. This flank of the valley is the western border of the Beaumont plateau, a more or less flat structure made of quaternary sediments and located at an average elevation of around 850 m above sea level (asl). The physiography of the study area is mainly controlled by its quaternary history. During the Last Glacial Maximum (LGM), the valley was dammed twice by the Bonne glacier coming from the east [40]. These successive phases of advance and retreat led to the deposit and erosion of lithologically contrasted layers, including moraines, glaciolacustrine, and deltaic sediments. Over the Jurassic bedrock and from the base to the top, the stratigraphy is as follows [40]. At the end of the interglacial Würm I sequence, alluvial sediments with a thickness of around 50 m were deposited by the paleo-Bonne river. These sediments are mainly made of sands and gravels but also contain pebbles. They are locally cemented. At the beginning of the next glacial period (Würm II), the valley was dammed by the Bonne glacier and a lake was created. Fine-grained sediments made of alternating layers of clays and silts (so-called laminated clays) settled in the lake over a thickness of around 20 m. These levels are associated to the Trièves clays in terms of age and of lithological content. During the peak of this glacial episode (LGM), the progradation of the glacier towards the west was accompanied by the deposition of 40 m-thick ground moraines with a high clay content. Following an episode of warming, the glacier melted and the lake was established once again. Around 60 m of laminated clays were first deposited. They are overlain by 40 m of lacustrine sands and then deltaic sands and gravels. This sedimentary sequence marks the progressive filling of the lake in the area. These deposits are overlain by a few metres of moraines containing big blocks and which are indicative of a further glacial episode (Würm III), with somewhat a more limited extension than the previous one. By analogy with the Trièves area, the definitive melting of the glacier allowed the Bonne river to erode into these soft sediments and to initiate landsliding [33]. Numerous landslides were reported and studied along the Bonne valley since more than 35 years [25, 26, 32, 41, 42]. Synthetically, it was shown that these landslides were settled within laminated clays and revealed several shear surfaces at depths ranging from a few metres to up to 30 m to 40 m. Moreover, it was evidenced that the activity of these landslides was seasonal, with higher deformation occurring during autumn and winter and which originates from snowmelt and heavy rain falls.

Over the study site, the thickness of sediments over the Jurassic basement is estimated to be around 100 m. The vertical lithological heterogeneity induces heterogeneity in terms of hydrogeological and mechanical properties. The average slope of the valley is 18°. Two main landslide zones (labelled 1 and 2 in Figure 1(b)) were observed and monitored since the late 1970s. Several shear surfaces were evidenced by geotechnical prospecting, the two most active being located within the first 10–15 m below the surface. The first one is at a depth of around 5 m and is located at the interface between moraines and laminated clays. The second one is at a depth of around 12 m and is located within the laminated clays. These depths are similar to those of the most active shear surfaces affecting landslides in the Trièves area [27, 33, 36]. This landslide corresponds to a clay/silt compound landslide according to Hungr et al. [28].

On the study site, these shear surfaces affect a road with a consequent traffic of 1750 vehicles per day, 7% of which are trucks (http://www.enroute.mediterranee.developpement-durable.gouv.fr). Recent remediation works conducted in 2009 allowed improving traffic conditions and adapting the road to less intense deformation in zone 2. In this area, present-day motions are lower than 2 cm/yr. However, remediation works did not allow to decrease the motion of zone 1 which is still highly active. Remediation works consisted of retaining permeable structure along the uphill flank of the road. The road was also shifted a few metres uphill (towards the east) with the objective of moving it out of the main crack in zone 1 and in zone 2. Drainage structures were dug in the upper part of the slope to reduce and maintain the phreatic water table at a low level for increasing the security factor of the site (by decreasing pore pressure). Geotechnical design of those structures was performed by modelling and calculation (see Section 3.5). The purpose was to get a safety factor close to 1.3. This is lower than the safety factor of 1.5 usually required according to French and/or European standards [43]. However, the purpose was not to reach a full security factor. This would have led to a designing of huge retaining structures, hard to build on the site with the consequence of a significant increase in the cost of works.

3. Material and Methods

The site was investigated using various ranges of in situ techniques to determine the geological, geophysical, geotechnical, and hydrogeological characteristics of the site. Measurements started in September 2010 and ended in late June 2017. The different techniques used in this work are detailed in the following section, and a synthesis of the experiments is provided in Table 1.

3.1. Geological, Geophysical, and Hydrogeological Prospecting
3.1.1. Coring

A 25.2 m-long borehole, labelled s1 in Figure 1(c), was cored to study the lithological characteristics of the site. The average recovery was 72%. Core samples were visually analysed. This allowed identifying the distinct lithological units. It was also possible to visually identify disrupted levels corresponding to shear surfaces. Finally, samples were analysed in terms of grain size distribution and geotechnical parameters (identification and tests) for further geotechnical modelling.

3.1.2. 2D Electrical Resistivity Tomography

2 electrical resistivity tomography (ERT) profiles were conducted on the site, namely, pe01 and pe02 (Figure 1). This method is a classical technique for landslide investigation (for a review, see [44, 45]). Resistivity measurements were conducted with a multichannel ABEM Terrameter LS. pe01 was conducted along a NW-SE direction, with 62 electrodes with a constant spacing of 5 m. The objective of this profile was to gain information on the geological setting. 1105 measurements were acquired in a multigradient configuration [46]. pe02 was conducted at the top of zone 1 along a NNE-SSW direction. It consisted of 64 electrodes with a constant spacing of 2 m. One further roll-along (12 electrodes) allowed gathering a 150 m-long profile. The objective of this profile was to get a high-resolution image of the subsurface at the top of zone 1. Measurements were acquired in a multigradient configuration and also in a dipole-dipole configuration (see [47]). The 2964 experimental measurements were inverted jointly.

Apparent resistivity data were inverted using the boundless electrical resistivity tomography (BERT) algorithm developed by Günther et al. [48]. This code is based on finite element computations with irregularly shaped triangles [49]. The regularization strength (lambda) was manually tuned to get a chi2 value close to 1, indicating that the data were fitted within a default error of 3%. Inversion provided satisfactory results with relative root mean square (RRMS) errors below 5% after 3 and 2 iterations for pe01 and pe02, respectively.

3.1.3. Electromagnetic Profiling

Electromagnetic (EM) profiles were acquired using an EM31 conductivity metre (Geonics Ltd., Mississauga, Canada). The fundamentals of the technique are provided in reference books (e.g., [50], and specifics about the device can be found in McNeill [51]). In clayey landslide investigation, this technique is used to map the spatial distribution of apparent resistivity over the first metres below ground, noticeably to detect the presence of electrically conductive materials such as quick clays [52] or to highlight the geology below a weathered and disrupted cover [53]. Measurements were acquired with horizontal coplanar loops. This configuration mainly provides information on ground properties in the very first metres but with sensitivity to EM properties down to a maximum depth of around 5 m to 6 m. The EM device was coupled to a global navigation satellite system (GNSS) to georeference the measurements. A total area of 0.2 km2 was covered. Experimental data, with an average coverage of 0.5 pt/m2, were interpolated using a kriging algorithm with a mesh size of 5 m × 5 m. The result is a map depicting the spatial distribution of the apparent resistivity over the first metres below the surface.

3.1.4. Water Infiltration Tests

Finally, two infiltration tests were conducted at the top of landslide zone 1 with the objective of evaluating the steady-state infiltration rate (SIR) of two contrasted levels. The tests are based on the Beerkan infiltration method developed by Braud et al. [54] and improved by Lassabatère et al. [55]. It consisted in infiltrating known volumes of water with a variable head through a single annular ring of a diameter of 12.5 cm. The ring was set at a depth of 1 to 2 cm and further coated with clay to avoid lateral loss. First, 20 cl of water was poured and the time needed for the complete infiltration was measured. A second volume of water was immediately poured, and, once again, the time needed for the complete infiltration was measured. On the study site, three to five successive infiltrations were needed to achieve the experiments. These were conducted after removing the first 0.3 m of the ground to avoid root effects. For the first volumes, the water is infiltrated through the unsaturated zone in transient regime while, for longer times, the water is infiltrated through the saturated zone in steady-state regime. SIR is derived by computing the slope of the cumulated water height as a function of time, in steady-state regime (last experimental points). This method was successfully applied by Bièvre et al. [56] on the Avignonet landslide, in similar lithological units, to characterize the hydraulic properties of contrasted geological layers, and results by these authors will be compared to those obtained in the present work.

3.2. Landslide Deformation and Kinematics
3.2.1. Surface Kinematics

Monitoring of surface kinematics is a key point to landslide studies, noticeably to establish relationships between surface motion and the evolution of environmental parameters (such as rainfall and water table [6, 33, 5759]). 50 to 60 geodetic stations were monitored using a real-time kinematic (RTK) differential GNSS (D-GNSS) which provided measurements with an error of a few mm relative to the base station. For each campaign, the location of the base station, which acquired measurements for at least two hours, was postprocessed using the Trimble Business Center software (geospatial.http://trimble.com). It provided maximum errors of 2.5 cm (highest error on the vertical component). It was difficult to set permanent geodetic stations in the south-eastern part of the site because of farming activity (cattle). In this area, metallic plates connected to a drainage system were used as permanent stations.

The monitoring started with around 50 stations in late November 2014. The number of stations progressively increased, and 60 stations were fully operational starting from mid-March 2016 until late June 2017. Annual average total displacement rates were computed between mid-March 2016 and late June 2017 using the three directional components. Data were interpolated using a kriging algorithm with a mesh size of 2.5 m × 2.5 m to produce a spatially continuous map depicting the average annual displacement rates. This small mesh size was chosen to evidence potential lateral variations where numerous stations were close by. Cumulated displacements were also analysed for three stations (r01 to r03) located in the south-eastern part of the slide, in the vicinity of the upper limit of landslide zone 1. Finally, horizontal displacement rate was also computed at station r01 which is located close to borehole i3. It was compared with the motion deduced from inclinometer data (see below). D-GNSS measurements were also used to build a digital elevation model (DEM) of the study site in late July 2015. Topographic cross sections were extracted from this DEM. Finally, the mapping of some structures, such as cracks, was also achieved using this technique.

3.2.2. Deformation and Kinematics at Depth

Two boreholes were drilled to monitor the evolution of the deformation at depth as a function of time using inclinometer probes. Details about this kind of instrumentation can be found in reference books (e.g., [60]), and recent applications to landslides can be found in Uhlemann et al. [59] and in Gullà et al. [61]. Measurements and processing were conducted periodically (typically between 1 and 3 months) following French standard NFP94-156 [62]. Succinctly, a grooved tube is introduced into a borehole and is further sealed to the ground. The measuring probe is inserted into the tube, and sprung wheels are used to guide the probe along the grooves. The inclination with respect to the vertical axis is measured by accelerometers. Measurements were conducted each 0.5 m, from the base to the top of the tube. Using trigonometric relationships, the inclination at each depth step is turned into a horizontal displacement with reference to the first measurement in time. Displacements are then expressed as the cumulation of individual displacements at each depth step from the base to the top of the tube. The first measurement in time serves as a reference measurement.

3.3. Meteorological and Hydrogeological Time Series

The effective rainfall or groundwater recharge is estimated with a soil-water balance based on characterization of evapotranspiration () and parameters characterising the recharge area (soil available water capacity (SAWC), runoff, and vegetation coefficient ). The effective rainfall is estimated according to the workflow calculation proposed by Vallet et al. [63]. This workflow requires the records of precipitation, air temperature, relative humidity, solar radiation, and wind speed within or close to the study area. Daily rainfall and average wind velocity data were obtained from a meteorological station operated by the French Office of Meteorology (Météo-France) and which is located around 2.5 km NW of the study site. Daily solar radiation data were measured by a station operated by Météo-France and located 30 km south of the study site. Temperature and atmospheric pressure data were obtained from a station located in the Avignonet landslide in the Trièves area, 12 km WNW of the study site, and at an elevation of 680 m asl (i.e., comparable to the study site). This station is operated by the French National Observatory on Landslides OMIV [64]. Hourly temperatures were turned into average daily data. Daily minimum and maximum temperatures were extracted from the hourly time series.

The estimation of the effective rainfall with the soil-water balance requires the determination of three parameters which are , , and SAWC. and are related and are expressed as follows: where is the reference evapotranspiration and is determined using the Penman-Monteith method [65]. depends on the stage of development of the vegetation (cattle grass on the study site). It varies from a minimum value during winter (0.85) to a maximum value during summer (1) according to Allen et al. [65]. An average value of 0.94 was used in this study. SAWC is deduced from the nature of the soil. From field observation and laboratory tests, it was estimated to be 30 mm.

Two pore pressure sensors were installed in late September 2012 in two distinct boreholes (p1 and p2 in Figure 1(c)) and with two distinct recording units. The sensor in p1 was set at a depth of 8.5 m, 0.5 m below the base of the morainic cover, to monitor the phreatic water table. The sensor in p2 was set at a depth of 12.5 m to monitor a water table located above a shear surface at 12 m depth. Both pressure sensors were installed in closed boreholes (see [60]). The monitoring systems were also ventilated to remove the effect of atmospheric pressure variation [66].

3.4. Cross-Correlation of Rainfall and Water Table Time Series

The linear relationships between rainfall and water table were further studied using cross-correlation. This mathematical technique is classically applied in hydrogeology, notably to study the linear dependency between rainfall and aquifer recharge [6770]. In landslide studies, this technique has however been little applied. Helmstetter and Garambois [71] used it on the Séchilienne rockslide (French Alps) to correlate seismic catalogues and environmental parameters (rainfall). They evidenced very short reaction times (<1 day) between rainfall and the initiation of rockfalls and of micro-earthquakes within the landslide.

Daily values were first standardized (subtraction of the mean and normalization by the standard deviation of the dataset) to get time series with a mean centred on zero and with a variance of 1. The 95% confidence interval , in this case, depends only on the length of the dataset [72, 73] and is classically of the form where is the length of the dataset. It was modified to take into account the decrease in the length of the dataset as the time series to be correlated is moved apart from the centre of the reference time series. For positive lags, it is expressed by where is the 95% confidence interval value at position along the time series. For negative lags, it is symmetric of for positive lags.

3.5. Stability Calculations

Stability calculations were conducted at the limit equilibria using the global perturbation method of the normal stress along a failure curve [74] with a Mohr-Coulomb criterion. This method was developed in France and is applicable to any shape of slip surface. A classical 2D model has been implemented to illustrate field observations in a simplified way. The model was built from geological observations (lithology, geometry), and geomechanical parameters were derived from laboratory tests. The topographic cross section (red arrow in Figure 1(b)) was derived from a DEM built from D-GNSS measurements. The geometry of the geomechanical layers was derived from boreholes located close to the cross section (i.e., i1, i2, s1, p1, p2, and other boreholes not presented in this work).

The rock riprap retaining structure was not modelled. This was chosen to compute the actual stability of this particular area of the site from a general point of view and a general stability site study consideration. The riprap structure is made of rock blocks which are not connected to each other. This structure was designed to act as a drainage structure but also as a protection for the road and its users. It has no stabilizing role.

Computations were conducted with the software Geostab 2013 [75]. An initial calculation was conducted with fixed mechanical parameters, especially the cohesion and the internal friction angle of the soil. These parameters were obtained from residual shear tests which provided residual cohesion and friction angle of the upper part of the laminated clays. The failure curve was designed based on inclinometer data and on the location of cracks observed in the field (see Figure 1). Then, two water table levels were taken into account for the calculations. Time series from pressure sensors p1 and p2 indicated a low level corresponding to the lowest measured level (in summer) and a high table level, corresponding to the highest level measured (in winter).

4. Results

4.1. Geological, Geophysical, and Hydrogeological Characterization

The detailed analysis of the samples from coring s1 allowed establishing the lithological characteristics of the site. From the surface down to a depth of 8 m, the lithology is made of a mixture of sands, gravels, and pebbles within a clayey matrix. Endogenous and limestone pebbles suggest that this level is partially made of moraines. These moraines were most probably reworked by surface processes and by the activity of the landslide. From an 8 m depth down to 22.3 m, alternating laminae of clays and silts, with a unitary thickness of the order of the cm, were found. These levels correspond to glaciolacustrine deposits. Distorted strata were observed at depth ranges of 10.7–11.4, 15–16 m, and 18.9–20.4 m. One missing level (no core recovery) between 7.2 m and 7.7 m deep probably corresponds to disrupted soils. Below 22.3 m and down to 25.2 m (end of the coring), another heterogeneous level was observed. It is made of a dark clayey matrix containing sands and gravels. Several pebbles with diameters ranging between 2 cm and 8 cm were observed. These pebbles were made of dark limestone and also endogenous rocks and were marked by scours, which suggests that these levels correspond to morainic deposits.

Figure 2 presents ERT profiles pe01 and pe02. The 315 m-long profile pe01 (Figure 2(a)) shows the presence of layers with a resistivity higher than 100 Ω·m (and up to more than 500 Ω·m) for elevations above 720 m asl. Below this elevation, resistivity is lower than 50 Ω·m, which corresponds to the laminated glaciolacustrine clays identified in s1 (Figure 2(a)). The interface between these low-resistivity levels and the overlying formations appears subhorizontal. The base of the clays seems to be located at an elevation of around 670 m to 680 m, and the underlying formations exhibit higher resistivity (around 150 Ω·m to 200 Ω·m). The interface between these two units also appears subhorizontal with undulations. However, it must be stressed out that these levels at the bottom of the geophysical images are very poorly investigated and that, consequently, it is difficult to interpret them confidently. Using the log established by Monjuvent [40] on an outcrop located around 700 m NE of the landslide, it was possible to define the stratigraphy of the study site. The upper resistive layers could correspond to a composite unit made of, from the base to the top, lacustrine laminated fine sands, deltaic sands and gravels of Würm II age and, finally, moraines of Würm III age. The underlying laminated clays correspond to glaciolacustrine deposits of Würm II age. These deposits are similar to the prone-to-landsliding Trièves clays, in terms of age and geomechanical properties. The badly resolved more resistive layers located below 670 m to 680 m asl could correspond to the Lower Moraines unit. This interface is located at least 20 m below the morainic levels detected at the base of s1 (Figure 2). However, this interface appears continuous on ERT profile pe01 and, below the base of s1, low resistivity levels are clearly identified. This apparent discrepancy remains difficult to interpret, and further prospecting is required to solve this issue.

ERT profile pe02 was conducted at the top of landslide zone 1 (location in Figure 1(b)), perpendicular to the slide motion. It shows the presence of two main units. The top unit is made of around 2 to 3 m of material with a resistivity of around 150 Ω·m to 200 Ω·m. However, between 70 m and 105 m along the profile, the resistivity is higher (200 Ω·m to up to 500 Ω·m) with a thickness increasing to around 6 m at a distance of 80 m to 90 m along the profile. Below are found lower resistivity levels (30 Ω·m to 100 Ω·m) which could correspond to the glaciolacustrine laminated clays.

The map in Figure 3 presents the spatial distribution of apparent resistivity over the first metres below the surface. Apparent resistivity up to 75 Ω·m is observed in the Northern and Western parts of the site (i.e., downslope). This location and elevation correspond to clayey units as observed on nearby outcrops [40]. Resistivity then increases from around 75 Ω·m to up to 200 Ω·m southward and eastward (upslope). This location corresponds to the morainic cover associated with lacustrine fine sands and deltaic sands and gravels observed on ERT profiles (Figure 2). Apparent resistivity values are slightly higher than those obtained with ERT for both moraines and clays. This probably originates from superficial effects (surface processes, farming activity, landslide activity, etc.) which slightly increase the porosity of the ground and, consequently, its resistivity. More interestingly, the relatively high resistivity zone in the south-eastern and uphill parts of the figure appears to be connected to the upper part of landslide zone 1. This observation is in good agreement with ERT profile pe02 which shows, at the same location, a thickening and an increase in resistivity of the morainic unit at the surface (Figure 2(b)). These higher resistivity levels are interpreted as originating from a higher proportion of coarse elements (and, correlatively, a lower clay content). Consequently, it suggests that a channel of higher permeability could preferentially drive the water flux from the upper moraines, deltaic gravels, and lacustrine sands located above around 720 m asl (Figure 2(a)) down to the top of landslide zone 1. This is materialized by a blue arrow in Figure 3.

To confirm the hypothesis of a high permeability channel, water infiltration tests w1 and w2 (location in Figure 3) were conducted within high and low resistivity levels (moraines and clays, resp.). Results are presented in Figure 4. The curves depict the cumulated infiltrated water as a function of time. The steady-state regime was established relatively shortly after the experiments started but the last three values were used to compute the saturated infiltration rate (SIR). The obtained SIR values are 0.57 × 10−5 m/s and 1 × 10−5 m/s for w1 in the moraines and for w2 in the laminated clays, respectively. This factor of two confirms the difference in permeability and, hence, in grain size, between the two units. However, it must be stressed out that measurements were conducted at a depth of 0.3 m and that, at this depth, surface processes at stake (weathering, farming, landslide activity, etc.) prevent from considering these values as representative of the subsurface at depth. Furthermore, measurements were conducted using a ring with a diameter of 0.125 m. As such, measurements are very local. Despite these limitations, the factor of two evidenced from these tests give an order of magnitude of the difference between the two units.

4.2. Surface and Depth Kinematics

Figure 5(a) presents the map of surface displacement rates measured at the 60 geodetic stations between mid-March 2016 and late June 2017. Velocities range between less than 2.5 cm/yr (i.e., within the error bar) in the Western and Northern parts of the slide and up to 15 cm/yr in the south-eastern part of the slide.

Figure 5(b) presents the topographic cross section along with displacement rates measured at experimental geodetic stations. The displacement rate is almost 12 cm/yr at the south-east, close to the landslide limit. Then, the rate increases to 15 cm/yr downhill and further decreases to 6.2 cm/yr on top of the rock riprap. On the road close to the base of the riprap, the displacement rate decreases to 2.5 cm/yr. These results suggest that the high rates observed in the SW part of the slide are presently strongly attenuated by the riprap structure. It also suggests that this structure, made of individual rock blocks, might undergo mid- to long-term irreversible deformation (block falls).

Finally, Figure 5(c) details displacements measured at geodetic stations r01 to r03 (location in Figure 5(a)). r01 and r02 are located in the immediate vicinity of the upper part of landslide zone 1. They exhibit relatively linear trends with displacement rates of 11.3 cm/yr (r01) and 9 cm/yr (r02). Considering the limited number of measurements and the width of the error bars (+/− 2.5 cm), it is not possible to determine if the motions are affected by seasonal changes. Unlike r01 and r02, r03 is located away from landslide zone 1. Except between mid-March and mid-June 2016 where it shows a motion of almost 16 cm/yr, it exhibits no movement since mid-June 2016.

Results from inclinometer monitoring are presented in Figure 6. The first inclinometer (i1, location in Figure 5(c)) was set to a depth of 25 m. It was definitely sheared 25.5 months after it was installed (from mid-October 2010 to late November 2012). Results show the presence of 2 shear surfaces at depths of around 6 m and 12 m. From 13 m to 17 m deep, the curves suggest the presence of one further shear surface. Finally, from 17 m down to the bottom of the inclinometer, the curves are not stable, which suggest that the inclinometer was not set in a stable area with no deformation. Maximum measured displacements are around 15 cm, half of which are accommodated by the first shear surface at 4 m depth.

Inclinometer i2 (location in Figure 5(c)) was set to a depth of 40 m. It was sheared 17 months after its installation (reference measurement in mid-October 2013 and last measurement mid-March 2015). The average horizontal motion close to the surface (at 0.5 m depth) is 7.4 cm/yr. At geodetic station r01, located 10 m W of i2 (Figure 1(c)), the average horizontal motion was 5.6 cm/yr. Considering the error on measurements with GNSS systems (around 2 cm along the horizontal plane) and also that the monitored periods do not match exactly, these motions derived from variable instruments are in good agreement. In the same way as i1, the same deformation mechanisms are observed from the surface down to a depth of around 15 to 17 m. Below, a shear surface is observed at a depth of 29–30 m. Finally, below a depth of 29–30 m, no noticeable movements are observed, which suggests that the base of the landslide is located at this depth. From a quantitative point of view, 52% of the total deformation is accommodated by the shear surface at around 5 m depth, 25% at around 12 m depth, 15% at around 15 m depth, and finally, 8% at 30 m depth.

For both i1 and i2, the maximum deformation was observed at a depth of around 4 m to 5 m. Cumulated displacements at this depth are presented in Figures 6(c) and 6(d) for i1 and i2, respectively, as a function of time. Average motions for i1 and i2 at 4 m depth are 7.2 cm/yr and 7.8 cm/yr, respectively.

The curves for both i1 and i2 exhibit a nonlinear increase of the cumulated displacements. However, separating annual seasonal cycles between, on the one hand, a winter season from November to March (blue bars in Figures 6(c) and 6(d)) and, on the other hand, a summer season from April to October, leads to the discrimination of relatively linear, and distinct, displacement behaviours for both inclinometers. Average displacement rates were computed for each season and appear to be consistent for both inclinometers: from 7 cm/yr to 16 cm/yr during winter and 0 cm/yr to 2.5 cm/yr during summer.

4.3. Stability Calculations

The geotechnical model derived from geological and deformation results was then used to compute the slope stability and to determine a safety factor along the studied failure curve. Three geotechnical strata were taken into account (Table 2): upper moraines and lacustrine sands from the surface down to 3–5 m (Moraines), upper glaciolacustrine clays with residual mechanical parameters (upper laminated clays) down to a depth of 8–12 m, and lower glaciolacustrine clays (lower laminated clays) with classical mechanical parameters. Geotechnical parameters (Table 2) were derived from direct shear tests on samples originating from coring s1. These parameters were compared to previously published results in the area and showed a good agreement [24]. Residual parameters were obtained with residual shear tests. The difference between upper and lower laminated clays originates from results on coring s1. Actually, this depth of 12 m corresponds to the location of a shear surface (Figure 6).

The geotechnical model is presented in Figure 7. Several failure curves were computed by the software. The selected failure curve was chosen as the one matching field observations (location of shear surface at depth in inclinometers and location of cracks in the field) but also as the one providing the lowest safety factor. The selected shear surface corresponds to the most superficial share surface observed in inclinometer i1 at a depth of 5 m (Figure 6(a)) and was considered because it is the most active (50% of the total deformation).

Then, two calculations were performed which differ only by the location of the phreatic water table (lower in summer than in winter). This water table is located immediately above the interface between the moraines and the glaciolacustrine sediments. Results indicate a safety factor of 1.23 during summer (low water table level) and of 1.15 during winter (high water table level). From a global point of view, this means that stability is not critical but is however below the 1.5 value required for long-term stability [43]. Furthermore, it indicates that an increase in the water table level decreases the global stability of the site. This decrease of the safety factor in winter is in agreement with inclinometer and D-GNSS data (Figure 6) which show an increase in kinematics and in deformation during the winter season. This could also explain the opening of cracks observed during the winter season in the field.

However, it must be stressed out that these computations were conducted for the most superficial shear surface located at 5 m depth and that at least three supplementary shear surfaces, encompassing around 50% of the total deformation, were observed on inclinometer results from i2 down to a 30 m depth (Figure 6(b)). Furthermore, the difference between those two results remains small and highlights the progressive degradation of the stability conditions in this area. This might originate from the progressive change of the geometry of the upper part of the slope (towards SE) which becomes more unfavourable. At this stage, however, it is only possible to conclude that, despite remediation works, this area of the landslide is substable but that the global stability remains tenuous. Moreover, the model and the computation have to be improved. Supplementary geotechnical data (such as mechanical parameters at different depths and additional inclinometer data located in a profile along the slope) are required. They would noticeably allow taking into account the strong anisotropy of the laminated clays [42, 76] and improving the geometry of the failure curve. In agreement with Van Asch et al. [26], there is also a more or less strong discrepancy between parameters defined at the scale of the laboratory and at the scale of the field. Finally, the global stability of the area could be computed using other calculation methods which would be able to encompass the different shear surfaces and water tables observed on inclinometer data (Figure 6) and piezometers (Figure 8), respectively.

4.4. Cross-Correlation Results

Figure 8 presents the data used for the quantitative hydrogeological study. The total amount of rainfall is depicted in Figure 8(a) under the form of daily rainfall and of the cumulated rainfall over the studied time range. Rainfall has been further separated into summer and winter rainfall, on the basis of the kinematics analysis from inclinometer data (Figure 6). The time series shows that rainfall is regularly spread along the year with average amounts of 1 m per year (0.9 m in 2016 to 1.2 m in 2013). In the same way as for total rainfall, effective rainfall is presented in Figure 8(b). The average amount of effective rainfall over the studied period is 0.19 m per year (0.12 m in 2016 to 0.26 m in 2013). From the total rainfall time series (Figure 8), it is observed that rainfall is more or less regularly spread along years. However, the effective rainfall in Figure 8(b) shows a slight discrepancy between winter and summer. In winter, effective rainfall is more important especially because of reduced evapotranspiration (no vegetation, lower temperatures) with an average amount of 0.1 m between November and the end of April. In summer, evapotranspiration is more important and induces a reduced effective rain (average amount of 0.09 m over the studied period). Furthermore, run-off is favoured in summer whereas, on the contrary, infiltration is predominant in autumn and winter. This suggests that the variation in kinematics could be related to seasonal effects related to water infiltration within the ground.

Figure 8(c) depicts the variation with time of the two studied water tables. The phreatic water table monitored by piezometer p1 is located at an elevation of around 716 m to 716.5 m asl, 2 m above the interface between the morainic cover and the laminated clays (which is also the base of the aquifer). It shows maximum annual amplitudes of around 0.3 m to 0.4 m. The water table monitored by p2 is located at an elevation of around 713.3 m asl, around 2 m above the shear surface located at 12 m depth (corresponding elevation of 711.5 m asl). It shows a lower amplitude of variation (0.15 m to 0.3 m) compared to the phreatic water table. Both water tables exhibit the same trends. The recharge occurs between mid-to-late November and early January and is relatively rapid (2 to 3 weeks). The following decrease is slow and may last up to several months. Both water tables show similar variations (among others, labels 1 to 5 in Figure 8(c)), which suggests that they react almost simultaneously to the same meteorological events.

To evidence more quantitatively linear relationships between rainfall and water table levels, these time series were then cross-correlated and the results are presented in Figure 9. Computations were conducted for both total and effective rainfall time series. Figures 9(a) and 9(c) depict the cross-correlation over the whole time lags for p1 and p2, respectively. They reveal only one peak with a confidence greater than 95%. Both cross-correlations exhibit a positive peak (i.e., a rainfall event induces an increase of the water tables) but with different lags. Figures 9(b) and 9(d) present close-ups around these peaks. For p1, the first peak has a cross-correlation coefficient (CC) value of 0.1 for a positive lag of 5 days (Figure 9(b)), whereas, for p2, it has a CC of almost 0.13 for a positive lag of 2 days (Figure 9(d)). This first observation implies that p2 and, hence, the water table located above the shear surface at 12 m depth, reacts to a rain event 3 days before the phreatic water table. To control this observation, p1 was cross-correlated to p2 and the results are presented in Figure 9(e) (whole cross-correlation curve) and in Figure 9(f) (close up around 0 lag). The results show several significant peaks with a one-year-long periodicity, suggesting that the two water tables are both controlled by the same environmental parameters (rainfall). The first positive peak is observed at a positive lag of 3 days with a strong CC of 0.64. This lag confirms previous observations, and the relatively strong CC once again suggests that the two water tables behave similarly. It must be stressed out that the three time series were conducted with independent devices (different recording units). The shape of the curves after the first peak also brings insight into the hydrogeological behaviour of these two reservoirs. Following the first peak at five days for the phreatic water table p1 (Figure 9(b)), CC remains stable at a value of around 0.1 until a lag of 21 days and then it decreases rapidly to a nonsignificant value for a lag of almost 30 days. This indicates that rainfall (and subsequent water infiltration) has a constant effect on the phreatic water table up to three weeks after a rain event (i.e., the water table stays at a constant level). Furthermore, the decrease in the water table after a lag of 21 days suggests that a drainage of this aquifer is initiated. Differently, the water table monitored by p2 shows a regular decrease in CC just after the first peak, and nonsignificant CC values are reached for a lag of around 42 days (Figure 9(d)). These results suggest that contrarily to the phreatic water table which behaves like a short-term reservoir, the water is drained directly after reaching the water table located above the shear surface at 12 m depth. From a general point of view, there is no discrepancy between cross-correlations conducted with total and effective rainfall time series. This suggests that, when correlating water tables with complete rainfall time series, the use of effective rainfall is of very limited added value.

To further evaluate a possible influence of the seasonal variability of water table level to rainfall, the time series were further distinguished between spring and summer (April to October, further labelled summer) and between autumn and winter (November to March, further labelled winter) seasons. The distinct time series were then cross-correlated, and the results are presented in Figure 10. Results for winter and summer using effective rainfall are presented in Figures 10(a) and 10(b), respectively. In winter, p1 (resp., p2) shows a peak with cc around 0.145 (0.14) 5 (3) days after a rain event (Figure 10(a)). The correlation coefficient for p1 (resp., p2) remains stable up to 22 (19) days after rain events and then decrease to nonsignificant levels around 33 (51) days after rainfall. These globally similar results suggest that the two water tables have a comparable behaviour in winter. On the contrary, they exhibit very contrasted evolutions in summer (Figure 10(b)). p1 shows a first peak () 6 days after a rain event. The influence of rainfall on the water table level lasts around 17 days, and then the cross-correlation coefficient decreases to a nonsignificant level 27 days after rain events. p2 reveals a first peak () 3 days after rainfall and a subsequent fast decrease of the correlation coefficient to nonsignificant levels 6 days after rain events.

p1 and p2 in winter and p1 in summer exhibit comparable behaviours, despite that p2 reacts always before p1 (3 days compared to 5 to 6 days). The first peak is followed by a plateau-like shape, indicating that the infiltrated rain maintains the water table at a more or less constant level for around 17 days (p1 in summer) to 19 days and 22 days (p2 in winter and in summer, resp.). This limited but effective reservoir effect exists for both water tables in winter and for the phreatic water table (p1) in summer. On the contrary, the deepest water table in p2 reacts quickly in summer (3 days after rainfall) and the water table decreases directly after the peak with CC back to nonsignificant values 3 days after the peak. It indicates that, in summer, the infiltrated water is not stored in the water table, even briefly. Instead, it is rapidly drained from the water table located above the shear surface at 12 m depth.

Seasonal cross-correlations were also conducted using total rainfall and water table time series. Results are presented in Figures 10(c) and 10(d) for winter and summer, respectively. The shapes of the correlation curves between effective and total rainfall exhibit similar trends, especially in winter (Figures 10(a) and 10(c)). They mainly differ by their correlation coefficient, especially for summer time series. Using total rainfall, the cross-correlation in summer does not show a significant peak and a further plateau-like shape for the phreatic water table (Figure 10(d)) whereas it does when using effective rainfall (Figure 10(b)). In an opposite way, the correlation curve of p2 with total rainfall in summer reveals a short-duration peak followed by a 40-day-long plateau-shaped positive correlation with a weak but significant value of 0.08 (Figure 10(d)). The consequent interpretation deeply differs since, in the case of correlations using total rainfall time series, it could be suggested that rainfall has a long-lasting influence on the water table located above the shear surface at 12 m depth. However, from the processing of these seasonal data with effective rainfall, it is evidenced that, once infiltrated, the water is rapidly drained out of this water table. These observations confirm the interest of using effective rainfall to evaluate their influence on water tables. As a matter of fact, the interpretation of the behaviour of the water tables in summer would have been opposite using total or effective rainfall.

Finally, seasonal water table time series were cross-correlated and results are presented in Figures 10(e) and 10(f) for winter and summer, respectively. In winter, a positive and strong correlation is observed for a lag of 2 days (), indicating that the two water tables have similar evolutions with time. On the contrary, in winter, the peak is badly defined and is observed for a lag between 0 days and 2 days with weak correlation coefficients ( and , resp.). This observation suggests that, in summer, there is no linear dependency between the two time series and that the two water tables have different behaviours.

To summarize, it has been evidenced that water infiltrates faster down to the water table located above the shear surface at 12 m depth than down to the phreatic water table located at around 7 m depth. This fast infiltration is favoured by a preferential infiltration path corresponding to a shear surface. Splitting time series between dry and wet seasons reveals a comparable behaviour of the phreatic water table (both seasons) and of the deeper water table in winter. Infiltrated rainwater reaches both water tables within a few days and is stored during a few tens of days (typically around 20 days) before the water is drained. On the contrary, in summer, the water reaches the deep water table within a limited time (typically 3 days) and is immediately drained out of the water table within a comparable time delay.

4.5. Conceptual Hydrogeological Model

The synthetic hydrogeological and geotechnical conceptual model of the study site is exposed in Figure 11. The cross section is located in Figures 1(b), 3 and 5. It can be considered as representative of landslide zone 1. The geometry of the lithological units was derived from geological and geophysical prospecting, and the location of the shear surfaces was deduced from inclinometer data. Time lags deduced from cross-correlation using effective rainfall (Figures 10(c) and 10(d)) were used to compute average infiltration rates (or, equivalently, velocities of pore pressure transmission from the surface down to the top of the water table). From the surface down to the phreatic water table, assuming a vertical distance of 7 m and a lag of 5 days in winter and 6 days in summer, the average infiltration rate is 1.48 × 10−5 m/s (1.62 × 10−5 m/s in winter and 1.35 × 10−5 m/s in summer). From the surface down to the water table above the shear surface at 12 m depth, assuming a 20 m-long path along the shear surface (Figure 11) associated to a lag of 3 days (both in winter and in summer) drives to an infiltration rate of 7.7 × 10−5 m/s.

Results from the infiltration rates and from the water infiltration tests in Figure 4 were compared to similar tests conducted on the Avignonet landslide [36, 56]. Results are synthesized in Table 3. They show that infiltration rates in the moraine and along fractures are in good agreement over the two sites. In the same way, SIR values in moraines and in the laminated clays are identical. This comparison tends to confirm the results obtained in this work and suggests that the established conceptual model for water infiltration is realistic.

Results from both sites show that the water infiltration rate along fissures and shear surfaces is faster than through the morainic matrix. This rate is faster by a factor of at least two and up to five at Avignonet and Charlaix landslides, respectively.

5. Synthesis and Conclusions

An extensive geophysical, geotechnical, and hydrogeological study was performed between 2010 and 2017 on the slow-moving Charlaix clayey landslide in the French western Alps. This landslide showed two main sliding zones with motions up to 10 to 15 cm/yr. Following remediation works achieved in 2009, one specific zone continued to exhibit motions of the same order of magnitude. The objective of this study was to establish a geological and geotechnical model following numerous geotechnical and geophysical prospecting conducted on the study site and to understand the seasonal variations of displacement rates observed in situ.

Geological and geophysical results showed that the site is complex with a heterogeneous vertical and lateral arrangement of lithological layers. A thick (20 m to 30 m at least) permeable reservoir uphill, furthermore connected to the upper part of landslide zone 1, was evidenced (Figures 24). Geotechnical prospecting showed the presence of 4 main shear surfaces at depths of 4–5 m, 12 m, 15 m, and 30 m with more than 75% of the deformation being accommodated within the first 12 m below ground surface (Figure 6). Surface and underground kinematics (Figures 5 and 6) revealed a strong seasonality of the deformation with almost no displacements during dry seasons (April to October) and motions up to 16 cm/yr during wet seasons (November to March). Simplified stability computations taking into account the shallowest shear surface only and variable water table levels (highest in wet season) showed a decrease of the safety factor during the wet season, with a reduced global stability (Figure 7). This result is in agreement with surface and underground deformation observations which revealed higher velocity during these periods.

From geophysical and geological prospecting, it was hypothesized that a semipermeable reservoir exists uphill (Figure 2(a)) and which could feed landslide zone 1 with a significant amount of water. The drainage of water from this reservoir to the upper part of landslide zone 1 would be favoured by the presence of a permeable channel (Figures 2(b)4). However, no increase related to a water flux from this reservoir could be evidenced from the analysis of the time series (Figure 8).

Further hydrogeological monitoring was conducted during four years between late 2012 and late 2016 to characterize the evolution of two water tables (Figure 8). The first is the phreatic water table located at around 7 m depth, and the second is deeper (around 12 m depth) and located immediately above a shear surface. The cross-correlation of piezometer time series with effective rainfall showed that the deepest water table reacts twice faster to effective rainfall than the phreatic water table does (Figure 10). This confirms the role of fissures and shear surfaces as preferential infiltration paths for rainfall from the surface down to 12 m at least (in this work). Furthermore, the seasonal cross-correlation of the same time series reveals a globally similar behaviour of the phreatic water table in summer and in winter. Infiltrated rainwater reaches the water table within a few days (5 to 6 days) and is stored during a few tens of days (around 20 days) before the water is drained out of the aquifer. On the contrary, the second water table showed a very contrasted seasonal behaviour. In winter, it behaves similarly to the phreatic water table but, in summer, it acts as a drainage structure with no significant effect of rainfall on the water table six days after a rain event. Nevertheless, the hydrogeological analysis of the time series showed that the two water tables are both sensitive to rainfall. This implies that water infiltration within the ground affects both reservoirs even if lags are observed. These new observations and analysis bring insight into the mechanisms of water infiltration at depth within clayey landslides and confirm the role of fissures connected to shear surfaces as preferential infiltration paths.

However, despite the original observations conducted in this work, several issues remain to understand the relationships between water infiltration and landslide activity. First, the balance between run-off and infiltration is required to establish more robust correlations. It would allow to discriminate the quantity of water actually infiltrating into the ground and reaching both water tables. Second, the impact of pore pressure increase on deformation, especially at the vicinity of shear surfaces, still needs to be assessed. A joint monitoring of these two parameters would allow to quantitatively evaluate the physical phenomenon leading to landslide acceleration (dilatancy and/or fluidization). The experimental monitoring of these parameters at depths of several metres within active landslides remains challenging.

Data Availability

Time series used for this work are available at https://www.researchgate.net/project/The-Charlaix-landslide-French-western-Alps. They are also available upon request.

Conflicts of Interest

The authors declare that they have no conflicts of interest.

Acknowledgments

This work benefited from a financial support from the VOR program (Université Grenoble Alpes) and from ISTerre. Grégory Bièvre is part of LabEx OSUG@2020 (ANR10 LABX56). This work is part of the project Anticipation and Mitigation of Slow Landslides in the Alps (MLA3) that benefited from the support of the European Union and Provence-Alpes-Côte-d’Azur region through the European Regional Development Fund (ERDF) interregional operational program “Massif des Alpes” (POIA). The authors thank two anonymous reviewers whose comments helped to improve the manuscript.