Hydraulic hysteresis of natural pyroclastic soils in partially saturated conditions: experimental investigation and modelling

In many geotechnical applications, especially in the study of weather-induced landslides, a reliable soil hydraulic characterization in unsaturated conditions is required. Currently, the experimental techniques that neglect the hydraulic hysteresis represent the greatest limitation to landslide forecasting. In this paper, a procedure to obtain an unsaturated soil hydraulic characterization on natural pyroclastic samples is proposed and verified. The approach enables the evaluation of the soil hydraulic properties along the main drying path and wetting/drying cycles to fully quantify the effects of the hydraulic hysteresis. Pyroclastic soil samples collected at a test site at Mount Faito in the Campania region (southern Italy) were tested. The experimental investigation consisted of a sequence of testing phases: a constant-head hydraulic conductivity test, a forced evaporation test followed by several wetting–drying cycles, and a drying test in a pressure plate apparatus. The hysteretic model proposed by Parker and Lenhard (1987) was adopted to fit the data, while inverse modelling of the forced evaporation tests allowed to derive the model parameters. Therefore, the main drying and wetting branches and the soil response to drying and wetting cycles from any reversal point were reproduced with the model, which suitably described the hysteretic behaviour of the pyroclastic soil under all conditions and along all paths.


A
Cross-sectional area of the soil sample CV b l Coefficient of variation of parameter b l d Distance between the tensiometer tips G s Specific gravity i Hydraulic head gradient between top and bottom of the sample K s ð Þ Soil hydraulic conductivity K sat Saturated hydraulic conductivity K 0 Hydraulic conductivity at the null suction K d sat Saturated hydraulic conductivity along the drying branch K w sat Saturated hydraulic conductivity along the wetting branch k-P-S Hydraulic conductivity -pore pressuresaturation hysteretic model ' Fitting parameter of the hydraulic conductivity function m Fitting parameter of the k-P-S model n v Fitting parameter of van Genuchten equation n Porosity R 2 Coefficient of determination S Effective degree of saturation s Matric suction s bot Matric suction measured by the bottom tensiometer in the ku-pf apparatus S s;b l z; t ð Þ Dimensionless scaled sensitivity of the matric suction to generic parameter b l s top Matric suction measured by the top tensiometer in the ku-pf apparatus

Introduction
Several geotechnical problems in unsaturated soils require an in-depth understanding and modelling of the hydraulic relation between matric suction and water content, which is known as the water retention curve (WRC) [36]. These curves play an important role in the assessment of unsaturated soil property functions such as the hydraulic conductivity, volume change and shear strength [70,73]. For example, the proper determination of the WRC is essential for modelling shallow landslides triggered by rainfall infiltration into the soil slope [8]. The analysis of the phenomena governing landslide triggering requires both (i) monitoring of hydrological and meteorological processes [45,55] and (ii) hydraulic and mechanical characterization of natural soils [39,40,66]. Over the last two decades, many experimental studies have focused on the triggering mechanisms of flowslides and debris flows in the pyroclastic soils in Campania, southern Italy [7,59]. In particular, representative testing sites of different regional geological contexts were instrumented to monitor the weather conditions, soil matric suction, and soil water content to investigate the field hydraulic properties of pyroclastic soils [12,13,15,21,33,45,[48][49][50][51][52][53][54].
The characterization of the hydraulic soil properties using a single branch (wetting or drying) of the WRC and the hydraulic conductivity function (HCF) is insufficient for practical purposes. Indeed, the hydraulic soil behaviour is hysteretic, and infinite different scanning paths can be followed in drying-wetting cycles depending on the reversal point, i.e., the suction at which the transition from drying to wetting occurs or vice versa [22]. When transient boundary conditions are involved, the hysteresis in the hydraulic properties can strongly affect the water flow regime [41]. Neglecting hysteresis can cause considerable errors in the mass flux calculations and soil shear strength determination. Nevertheless, hysteresis is generally ignored in most practical applications, and the characterization of hysteretic hydraulic properties and their implementation in computer codes is uncommon [18]. This negligence occurs most likely because the determination of hysteresis with classical investigation methods is laborious and time-consuming due to the large amount of experimental data required for model calibration, complexity of the numerical analysis, and uncertainty of the ability of the existing models to describe the behaviour of natural soils.
The first experimental apparatus to study the hysteresis was developed by Fireman [20]. Over the last 30 years, several efforts have been conducted in different areas such as hydrology, agronomy and geotechnical engineering to experimentally analyse the hydraulic hysteresis at the element scale (soil specimen) and medium scale (physical model) [17,36,41,56,69]. Subsequently, to model the hydraulic hysteresis, physically based models [38,61] and purely empirical models became available [28,46]. While the hysteresis aspect in the water content-suction, i.e., h(s), and hydraulic conductivity-suction, i.e., K(s), relationships is notable, the hysteresis effect in the hydraulic conductivity-water content, i.e., K(h), relationships is generally considered negligible [38].
Scanning paths in pyroclastic soils were recently documented through a systematic comparison of the coupled volumetric water content and matric suction measurements at identical depths in instrumented slopes in the Campania region [4,21,49,51]. Therefore, the existing research on the hysteresis of the hydraulic properties of natural pyroclastic soils is mainly based on field measurements collected under transient boundary conditions (meteorological conditions).
In this framework, the experimental activity described in the present study attempts to determine the parameters necessary for the hydraulic modelling of infiltration phenomena in silty pyroclastic slopes by adopting a hysteretic model useful to set up a reliable tool for early warning systems against the triggering of flowslides.
With respect to previous studies on the hydraulic hysteresis of pyroclastic soils performed on site, we present a complete laboratory procedure of experimental techniques devoted to calibrate a hysteretic model. The laboratory study consisted of well-controlled flow experiments with monotonically changing boundary conditions to reduce the effects of the soil heterogeneity and minimize the noise and bias, which commonly affect the water content and matric potential measurements [16]. The procedure was composed of a sequence of tests that provided experimental data on the saturated soil hydraulic conductivity, soil response to evaporation and imbibition over time, and matric suction level at the residual water content. Thereafter, the experimental data were adopted to calibrate Parker and Lenhard model [46], which described the soil hydraulic behaviour along the (i) main drying and wetting branches and (ii) scanning paths. In particular, the model calibration has been achieved by performing inverse analysis of the boundary value problem (i.e. the evaporation and the imbibition phases of the lab test) which consists in a series of transient phenomena during which both the suction and the permeability vary. Thus, a hydraulic conductivitypore pressure-saturation hysteretic model (K-P-S model) suitable for natural pyroclastic soils was established, and its predictive capacity was proven to be very satisfactory.
This procedure which results from consecutively applying different experimental techniques to the same soil specimen enables a large number of tests and is not excessively time-consuming, so it is suitable to correctly manage the inherent variability of the considered soils. Indeed, the available hysteretic models in the literature and those implemented in commercial codes were conceived for sedimentary soils and generally validated on artificial samples [5,11,60,64]. Here, the model is calibrated and validated on natural undisturbed soils samples that are characterized by macro-pores, i.e., a meta-stable structure. In addition, it is proven that: (i) the main drying curve and only one drying-wetting cycle performed within the experimental procedure proposed, are enough to completely determine the parameters of K-P-S model; (ii) using discontinuous influx of water that detect a scanning wetting phase occurring step-by-step can be used to determine the fitting parameters of the main wetting curve that usually represent a challenging task as a continuous and slow influx of water into a sample is hard to perform experimentally.

Samples and soil physical properties
Undisturbed soil cores were collected from a shallow layer at the test site located at an altitude of approximately 850 m on a north-facing slope of Mount Faito (40°40 0 32.29 00 N, 14°28 0 23.35 00 E) in southern Italy [17]. On average, the soil, originating from the deposition of pyroclastic material resulting from 79 AD eruptions of Mount Vesuvius, is characterized by a very high porosity (approximately 0.70) and low dry soil unit weight (approximately 8.2 kN m -3 ). The porosity and dry unit weight of Mount Faito soil are similar to those reported for other pyroclastic soils in the Campania region [47]. The soil was partially saturated at the site throughout the year [17], and the volumetric water content measured for the specimens sampled at the site during the dry period (July 2016) varied between 0.10 and 0.15. The grain size distribution and soil physical properties are shown in Fig. 1 and listed in Table 1, respectively. The soil was classified according to the Unified Soil Classification System [1] as sandy silt and sometimes as silty gravel with sand. The complete geotechnical characterization of this soil has already been reported in the literature [17,21]. The soil cores were collected horizontally from the walls of a road cut at the site. In the present study, a total of 11 specimens were cut from the original soil cores using a metallic sleeve with a cutting edge. The metallic sleeve had a diameter of 70 mm and a height of 60 mm. The wall sleeves contained two holes at 15 mm from each end to install the necessary sensors, which are described later. The soil specimens always remained inside the sleeves throughout the testing, thus allowing the soil to remain undisturbed. The preservation of the soil structure resulting from the deposition process is necessary to allow a better compatibility with field data.
On each specimen, a constant-head hydraulic conductivity test under saturated conditions was performed in a permeameter, followed by evaporation-imbibition cycles in a ku-pf apparatus. Finally, the water content in the highmatric suction range (up to 1 MPa) was determined by drying soil specimens in a pressure plate. The procedure proposed in this study does not include measurement of the soil volume variation. In this regard, it is well known that some soil types, such as expansive and collapsible soils, and generally most of the fine-grained soils undergo significant volume change due to suction change at constant net stress [58]. In the last twenty years, some, and valuable efforts within geotechnical field, have been done to determine experimentally and to model the soil water retention curve in deformable soils, mainly in compacted and reconstituted soils [23,24,27,35,58,67,68]. This volume change can markedly affect the water-retention response of the soil and ignoring such dependencies can lead to misinterpretation of Water Retention Curve (WRC) data and to extraction of imprecise soil parameters [74]. However, previous experimental studies have demonstrated that the volume changes due to water content variations in the examined pyroclastic soils are non-significant [40,43]. In this regard, there are evidences from oedometric tests performed on undisturbed specimens of pyroclastic soils collected in the Campania region in [43]. In particular, the soil volume strain measured in oedometer tests by submerging the soil specimen and keeping the vertical net stress constant and equal to that operative at site, resulted on average equal to 0.010%. Therefore, settlements due to decrease in suction (resetting to zero) were negligible. The negligible amount of volumetric collapse and the incompressible behavior under hydraulic path at constant vertical stress exhibited by the pyroclastic soil specimens collected in Campania is well documented also in [6,40,44,65]. This incompressible behavior exhibited by volcanic ashes was also experimentally observed by [19,29]. In these works, limited volumetric deformation on volcanic ash from the sides of the Irazú Volcano in Costa Rica were measured, during complete drying paths for suction increments at low vertical stresses.
Due to the insignificant changes in volume upon drying and wetting observed in the above-mentioned studies on volcanic soil, the experimental procedure discussed in this work is considered suitable for ashy soils and for all soil types for which the volume variations expected are negligible. Therefore, the experimental data will be presented in terms of the volumetric water content instead of the gravimetric water content; in fact, because of the null volumetric deformations, it is possible to indifferently use both variables.

Hydraulic conductivity under constant pressure head
The experimental procedure started with the soil specimen saturation followed by a hydraulic conductivity test. The apparatus consisted of a rigid-wall permeameter with two pressure sensors to measure the water pressure at the bottom and at the top of the cylindrical soil specimen, in addition to three water reservoirs with pressure regulators and two burettes to measure the water volumes flowing in and out of the soil sample (Fig. 2a). The rigid wall permeameter, wherein the specimens were vertically installed, contained two O-rings to prevent water leakage and two porous stone discs protected by filter paper to prevent small soil particles from being dragged into the porous stones (Fig. 2b).
The soil saturation process was performed by applying a pressure of 5 kPa at the bottom of the sample and allowing water to flow through the soil specimen until no air bubbles exited the upper part of the specimen. Then, a pressure of 10 kPa was applied in the upstream reservoir connected to the permeameter bottom, while a pressure of 5 kPa was applied in the downstream reservoir connected to the permeameter top. Hence, distilled water flowed upwards through the specimen. The volumes of water flowing in and out of the specimen were monitored over time using the two graded burettes. When a steady condition was reached (i.e., the incoming and outcoming flow rates became equal), the specimen was considered ready for testing. At least four-flushing cycles were repeated per specimen, where each flushing cycle corresponded to the emptying of the burette to measure the inward flow of water.
The comparison between the porosity of the specimens and the estimated volumetric water content after the hydraulic conductivity test revealed that the specimens were not fully saturated because some air bubbles remained entrapped during the flushing. However, the matric suction that was measured immediately afterwards was close to null. Therefore, the measured values consisted of the hydraulic conductivity at the null suction (K 0 ), which is used here as an approximation of the saturated hydraulic conductivity. The hydraulic conductivity at the null suction (K 0 ) was estimated with Darcy's law for one-dimensional flow (Eq. 1), where A is the specimen cross-sectional area, i is the hydraulic head gradient between the bottom and the top of the specimen, and DV is the average water volume flowing through the sample over the Dt time range. The hydraulic conductivity was calculated for each flushing cycle, and the average of the estimated values was adopted as the final value of hydraulic conductivity at the null suction.

Wetting and drying cycles in the ku-pf apparatus
After the hydraulic conductivity test was performed, the specimens contained still in the steel sleeve were moved in a ku-pf MP10 apparatus (Umwelt-Geräte-Technik GmbH) (Fig. 3a). In particular, the equipment consisted of a starshaped sampler changer, whose top view is shown in Fig. 4a. Each changer arm held a basket wherein the specimen was vertically installed (Fig. 3b); it was equipped with two mini-tensiometers connected to a conditioning unit (Fig. 3c). The mini-tensiometers measured the matric suction in the range of 0-80 kPa at a resolution of 0.01 kPa. The bottom part of the specimen was sealed using a plastic paraffin film (Parafilm M) to prevent water evaporation and drainage. The cling film and metallic cap were placed on the upper part of the specimen to prevent water evaporation. Figure 4b shows a side view of the basket holding the sample, which indicates the location of the mini-tensiometers, height and width of the soil specimen and position of the mini-tensiometers. The mini-tensiometers were saturated by flushing water through the porous stone to remove any trapped air bubbles with a needle. The calibration was made using the setup in Fig. 3b, where the mini-tensiometers were inserted in a sealed chamber to which negative pressure can be applied using a pump with an analogical dial. A two-point calibration procedure was adopted by imposing a pressure of 0 kPa (i.e., atmospheric pressure) and a negative pressure of 50 kPa (measured in the pump dial).
Two small holes were drilled on the side of the specimen (at 15 and 45 mm from the top) with a guide (Fig. 3d) to ensure the horizontal installation of the two mini-tensiometers. The specimen was placed on the structure in Fig. 3d by aligning the holes of the mould with the guides. An empty tube was pressed through the guides into the soil, and its content was emptied. All removed soil was weighed on a precision scale to correctly estimate the water content over the new sample volume. The mini-tensiometers were carefully installed in the holes in the specimen to ensure good contact between the porous tip and the surrounding soil.
The basket was installed in the rotating changer arm and weighed every 10 min on a precision scale with a resolution of 0.01 g to register all variations in water mass (Fig. 3e). The arm placed the basket above the balance, which was installed on a lifting mechanism. Simultaneously with the measurements of the soil weight, the matric suction was measured by the two mini-tensiometers, and the vertical hydraulic head gradient was estimated according to Eq. 2. Here, s top and s bot are the matric suction values measured by the top and bottom mini-tensiometers, respectively; d is the distance between the mini-tensiometers (d ¼ 30 mm); c w is the water specific weight (c w ¼ 10 kN m À3 ). The hydraulic gradient obtained by Eq. 2 was defined so that a positive value corresponded to an upward water flow.
After the installation of the mini-tensiometers, a waiting period was necessary to stabilize the matric suction measurements, during which the sample remained sealed. The stable conditions were reached when the matric suction measurements detected hydrostatic pressure conditions in the specimen, i.e., when the difference between the top and bottom measurements was 0.3 kPa.
The first drying phase was begun by removing the top sealing cap and cling film from the specimen and enabling water to evaporate (upward water flux). The decrease in weight of the specimen due to water evaporation and the suction were recorded until the top mini-tensiometer approached 80 kPa. The drying phase was stopped by covering the top part of the specimen with cling film and a metallic cap. While the specimen was sealed, the water redistributed inside it and reached a new hydrostatic condition. After the first drying, the specimens were subjected to a series of wetting steps that consisted of repeatedly pouring a fixed quantity of distilled water onto the top surface of the specimen and allowing it to infiltrate and redistribute inside the soil. The top cap and cling film were only removed to add water at the beginning of each wetting step. If the matric suction values measured by both tensiometers exceeded 15 kPa, 5 g of water was added to the top of the specimen. However, to obtain more refined data near the air-entry value (AEV) of the WRC, 3 g of water was added if the measured matric suction was below 15 kPa. The threshold value of the matric suction, below which the amount of added water was reduced, was selected as the upper boundary for the expected values of the AEV, which was estimated according to the grain size distribution of the tested soil. Each new wetting step was initiated after a stable hydrostatic pore water pressure distribution was attained as indicated by the measured matric suction values. The duration of a wetting step was at least 2 h, so that sufficient data points were recorded for each step to easily model the process. In particular, 2 h were generally sufficient when the suction was lower than 20 kPa and the soil was close to saturation; thus, the water infiltration and redistribution processes were very fast.
Three or two cycles of drying and wetting phases were repeated several times for each specimen. The cycles always ended with a drying phase.

Water content measurement in the high-matric suction range
A significant change in slope of the WRC occurred in the matric suction range corresponding to the residual saturation; hence, a large increase in matric suction corresponds to a relatively small decrease in water content [69]. For the tested soil, this region of the WRC corresponds to matric suction values well above the measurement range of the ku-pf mini-tensiometers (i.e., the matric suction is notably higher than 80 kPa). Therefore, a pressure plate apparatus (Fig. 5) was used to obtain the water content by applying matric suction of 600, 850 or 1000 kPa through the axis translation technique. In the present work, only one point of the WRC per tested specimen was determined with the pressure plate apparatus, as suggested by Nicotera et al., Vanapalli et al. [39,69]. Thus, although the test in the pressure plate takes a long time, the entire experimental technique continues to save time.
At the end of the last drying phase in the ku-pf apparatus, specimens were removed from the baskets and placed in the pressure plate. The bottom surface of the specimen was put in contact with the porous stone of the pressure plate to enable water exchanges. The pressure plate was sealed, and an air pressure of 600, 850 or 1000 kPa was applied to the chamber. The porous stone was crossed by water but not by air, so the specimens in the chamber could drain until the suction in the soil was in equilibrium with the applied pressure in the chamber. The weight of the specimen was regularly measured. The test was stopped when the changes in soil weight were negligible. The final water content was determined using the gravimetric method: the specimens were weighed before and after being oven-dried for 24 h at 105°C.

Hysteretic K-S-P model
The hysteretic model adopted herein is implemented in HYDRUS-1D software [62] that numerically solves the Richards equation with the implementation of standard Galerkin-type linear finite element schemes. The code can simulate the water movement in one-dimensional variably partially saturated media by solving the Eq. [73]: where t is time, z the spatial coordinate, k(h(s)) the hydraulic conductivity function, c w the unit weight of water, h is the volumetric water content, and s is the matric suction. This model assumes that the soil porous medium is rigid (no volume deformations). However, within the framework of behavior of slope constituted of ashy soil, the hydraulic modelling of infiltration phenomena using the Richards equation correctly takes into account the amount of water stored in the soil sample as proved by [9,10,14,26,42,56,57]. In particular, in [42] a good match between the volumetric water content measured at physical prototype filled by pyroclastic soil collected from a site in Campania, and the value estimated by the model was found. Greco et al. [26] reproduced satisfactorily the measurements collected at the monitored pyroclastic slope in Cervinara (Campania) by using a numerical model solving the Eq. 3. The hysteretic model adopted in this study, indicated as the PL model, was proposed by Lenhard et al. [32] as a simplification of the Lenhard and Parker [31] and Parker and Lenhard [46] models for two fluids. The main drying and wetting branches of the WRC are described by the van Genuchten [25] equation, but each branch is characterized by its own fitting parameters, as expressed in Eqs. 4 and 5. Superscripts d and w refer to the main drying and wetting curves, respectively, h is the volumetric water content and s is the matric suction. According to the most practical observations, in this model, the fitting parameters n v and h r are assumed to be equal for both branches. Differences in the contact angle between the solid and water and irregular pore geometry effects on the hysteretic loop are considered via the fitting parameters a d and a w , under the assumption of a d a w . The difference between h d s and h w s with h w s h d s is due to the air entrapment upon re-wetting. This model assumes that: (i) the air phase is immobile when it becomes entrapped by the water phase; (ii) the air pressure is the atmospheric pressure; (iii) the WRC and HCF are unimodal; (iv) the soil porous medium is rigid (no volume deformations); (v) the ''dynamic'' effects in WRC, which are defined as an apparent dependence of soil hydraulic properties on the flow dynamics, are negligible. Therefore, the soil hydraulic parameters are steady and do not depend on the flowrate.
The scanning paths are scaled from the main branches following the method of Parker and Lenhard [46]. In particular, drying scanning curves are scaled from the main drying curve, and wetting scanning curves from the main wetting curve. The scaling procedure is well documented in [56,71]; for more details about the entire model, the readers can refer to [30,72]. The HCF is described by Eqs. 6-7 [25,37], where K s ð Þ is the soil hydraulic conductivity, l and m are the fitting parameters, m ¼ 1 À 1=n v , S is the effective degree of saturation (Eq. 7), and K sat is the saturated hydraulic conductivity.
Consistently, this formulation of the HCF was modified by Lenhard and Parker (1987) to account for air entrapment, so the reader can refer to the original paper [31]. This model does not consider the hysteresis aspect in K(h) but only in K(s). An analogous hysteretic procedure can be applied to the unsaturated hydraulic conductivity function K(s); the main branches of the hysteresis loop are determined by using Eqs. 4 and 5 depending on the investigated phase, i.e., drying or wetting, to calculate h in Eq. 7.

Inverse analysis
HYDRUS-1D allows to perform an inverse analysis of the boundary value problem. In this case, the evaporation and the wetting phases of the lab test have been simulated. Once choosing the K-P-S model as the hydraulic constitutive model, the code is able to compute all the parameters of the constitutive model that allow the best agreement between the matric suction calculated and the matric suction measured during the test, by means of an optimization algorithm. Therefore, the modelling of a series of transient phenomena (evaporation and wetting phases) during which both the volumetric water content and the permeability varies, provide the best-fitting parameters of both WRC and HCF curves.
In particular, the fitting of experimental data with the hysteretic model was accomplished by means of inverse analysis conducted in two phases due to the large number of parameters. The first phase involved the parameter vector h d s ; h r ; a d ; n v ; l À Á associated with the main drying branch, and the second phase considered the parameter vector h w s ; a w À Á associated with the main wetting branch. The data sets were preliminarily filtered to remove outliers and data due to tensiometer malfunction. Subsequently, inverse analysis was performed using the HYDRUS-1D software where the objective function used to fit the data was minimized via the Levenberg-Marquardt nonlinear minimization method [34].
The fitting of the main drying curve followed the procedure defined by Nicotera et al. [39] after modifications. The initial condition was fixed as the hydrostatic pressure distribution estimated from the initial suction measurements. The water flow occurring within the specimen was vertical, with a null flux at the bottom and an upward flux at the upper boundary equal to the specimen weight change over a given time range. The value of K sat in Eq. 5 was assumed equal to K 0 in Eq. 1. The data sets and respective weights in the objective function were composed of (i) the matric suction values at the tensiometer positions during the monitoring time, with a weight of 1; (ii) pair s; h ð Þ obtained from the pressure plate, with a weight of 5; and (iii) pair s; h ð Þ corresponding to the AEV, with a weight of 5. The AEV was identified as the point of maximum curvature on the WRC, obtained by coupling the mean measured matric suction to the average water content estimated by the variation in soil weight recorded during the evaporation test. The adopted weights were chosen to compensate for the large number of suction measurements that, for this reason, greatly influence the results. The initial estimation and range of the fitting parameters were set as suggested by Nicotera et al. [39]. Parameters l, n v and h r are the same for both branches [32].
Only h w s and a w were fitted for the wetting branch. In particular, h w s was allowed to vary between 0 and h d s , and a w was allowed to vary between 0 and 100 kPa -1 . The initial estimation of h w s was obtained from the water content observed at the end of the wetting phase, while the initial value of a w was assumed to be equal to 2a d , as suggested in previous studies [28]. The boundary conditions adopted for the wetting phase reproduced the variations in water content in the specimen by imposing constant water fluxes at the upper boundary. Each wetting step consisted of an initial imbibition step within a very short time period (10 min), which increased the sample weight due to the added water. During the rest of each wetting step, an equalization and water redistribution process occurred. However, the upper boundary of the system was not perfectly sealed, and any water loss was registered by the ku-pf apparatus. Therefore, a constant evaporation flux at a very low rate was applied at the top boundary, while at the lower boundary, a null water flux was imposed. The data set used to fit the main wetting parameters consisted of the matric suction values measured by both minitensiometers with a weight of 1.
One of the purposes of this research was the identification of an accurate but expeditious experimental procedure that could be applied to study a large number of soil specimens within a reasonable time. Hence, it was decided to determine a compromise between the determination accuracy of the hysteretic model parameters and the minimization of the number of experimental data required to estimate them. This issue was investigated by determining the values of the coefficient of determination R 2 obtained by comparing all experimental data to the numerical simulation results of the whole test with the hysteretic model calibrated on the basis of only one and two drying-wetting cycles.

Sensitivity analysis
Before the model calibration phase, a sensitivity analysis was conducted to verify whether small changes in the fitting parameters resulted in large deviations of the model estimations [2]. The analysis in the present work followed indications by Nicotera et al. [39], who studied the model sensitivity to the fitting parameters along the main drying branch. Each fitting parameter was individually disturbed by adding or subtracting one standard deviation. The differences in matric suction estimated with the optimal values of these parameters and the disturbed ones were quantified through the function S s;b l z; t ð Þ, which is the dimensionless scaled sensitivity of the matric suction to the generic parameter b l of the parameter vector b. S s;b l z; t ð Þ was calculated according to: where CV b l is the coefficient of variation of parameter b l , r es is the estimated standard deviation of the matric suction measurement error (r es ¼ 0:112 kPa), z is the depth at which the matric suction was measured, t is the time, and Ds m z; t; b ð Þ is the matric suction difference between the best fitted and disturbed models. The disturbance adopted in the present work was equal to the standard deviation of the best fitted parameters calculated based on all tested specimens.
The dimensionless scaled sensitivities can be adopted to compare the importance of different observations to the estimation of a single parameter b i . Observations with large S s;b i values are likely to provide more information on parameter b i than those associated with small S s;b i values. Hence, the experimental setup should be designed such that the measurements yield the most information on the unknown parameters to be optimized, i.e., the measurements that are more sensitive to any changes in the unknown parameters [63].

Experimental data
The measured hydraulic conductivity at the null suction (K 0 ) varied between 1.93 9 10 -7 and 7.04 9 10 -6 m s -1 and spanned a range typical of silty sand ( Table 2). In some specimens, the hydraulic conductivity very slightly increased with increasing number of flushing cycles, but it generally fluctuated around an average value. Therefore, the average of the measured values over all flushing cycles was adopted as the final value of the hydraulic conductivity for each specimen.
However, these hydraulic conductivity values, which were determined at the end of the saturation phase, should be only considered a lower limit of the saturated permeability (K sat ) because the full specimen saturation was not always attained.
The estimated volumetric water content estimated in the pressure plate at matric suction values varying between 650 and 1000 kPa ranged from 0.129 to 0.164.
The test results for specimen no. 2 are shown in Fig. 6 as representing the observations obtained in all tests. Specimen no. 2 was subjected to three drying and two wetting phases. The matric suction increased during the drying phase with decreasing water content. The wetting phase consisted of a sequence of steps initiated as an abrupt increase in water content, which resulted in a sudden drop in matric suction measured by the top tensiometer. The pressure distribution inside the soil specimen tended to be hydrostatic (Fig. 6b), where the matric suction measured by the top mini-tensiometer should be approximately 0.3 kPa higher than that measured by the bottom minitensiometer. The specimen weight during each wetting step remained constant, but the water content did not remain uniform inside the soil specimen.
Experimental data were visualized via the WRC (i.e., volumetric water content versus matric suction), as shown in Fig. 7a (specimen 2) and Fig. 8a-d (specimens 1, 2, 3, and 10). The drying branches of the WRC were obtained by associating the average water content value to the mean matric suction value measured by the two mini-tensiometers. Due to the discontinuous nature of the wetting phase, only one point of the branch was derived per wetting step. The selected point per wetting step corresponded to the instant at which the gradient first reached a null value after water was added. As an example, to obtain point (A, B), as shown in Fig. 7a, which represents the sixth wetting step of the first wetting phase (wet 1), the null gradient value was identified, as shown in Fig. 7d. The mean matric suction value measured by the two mini-tensiometers (point B, Fig. 7c) and the volumetric water content (point A, Fig. 7b) corresponding to the same instant were then determined. This procedure was repeated for each wetting step.
The hysteretic features of the hydraulic behaviour of the tested soil are quite evident in Figs. 8a-d. These results indicate the following: (i) The hysteresis amplitude depends on the reversal point. For example, in Fig. 8a, the third cycle (E-   Fig. 8b. (ii) The suction corresponding to the knee along the wetting path is lower than that on the main drying curve (AEV) (see Fig. 8c). (iii) The amount of entrapped air can be estimated by observing the difference in the maximum volumetric water content measured between the main drying and wetting paths, h d s and h w s , respectively: this seems to be notable in specimen 1, low in specimen 2 and negligible in specimens 3 and 10 ( Fig. 8a-b). However, the amount of entrapped air is larger if the soil porosity of each specimen is considered instead of h d s . The reason is that full saturation is very difficult to attain. Fig. 6 Experimental measurements from the ku-pf apparatus of: a suction (top and bottom tensiometers) and volumetric water content (h) over time; b enlargement of a wetting step for specimen 2 Fig. 7 Experimental data on WRC plane: a comparison of the main loops obtained from the inverse analysis and experimental data for specimen 2. The method to obtain the points of the scanning wetting paths is schematized on the right-hand side: b volumetric water content; c suction; d vertical hydraulic head gradient of one wetting step 3.2 Determination of the parameters of the hysteretic K-P-S model Figure 9 shows the comparison of the fitted WRC and HCF to the experimental data set used for the model calibration of the main drying branch of specimen no. 2. The fitted model is consistent with experimental data. Additionally, the value of the coefficient of determination (R 2 ) was approximately equal to 1 for the tests on all specimens considered ( Table 3). The fitted values of the model parameters for the main drying branch of the WRC and the HCF reveal small variations among the specimens ( Table 3). The soil behaves like a coarse-grained material as a d varies between 0.06 and 0.17 kPa -1 . The residual volumetric water content (h r ) ranges from approximately 0.08-0.11. The volumetric water content at saturation (h d s ) varies between 0.55 and 0.59, but a lower value equal to 0.45 was obtained for specimen 2. The value of parameter n v , which affects the slope of the WRC, varies from 1.51 to 1.76. Parameter ' ranges from 0.90 to 3.00, exhibiting a very high variability.
Similarly, the fitting of the main wetting branch of the WRC provided an R 2 value approximately equal to 0.99 ( Table 4). The h w s =h d s ratio, which is an indicator of the air entrapment during wetting, ranges from 0.92 to 1.00, thus showing that as high as an 8% decrease in the maximum water content can be expected upon a wetting process produced by downward infiltration (Table 4). However, if h w s is compared to the soil porosity, n, instead of h d s , an average value of 0.75 is obtained for the h w s =n ratio with  [4]. Although the wetting process in the laboratory and in situ is directed downwards, a higher fraction of air could remain entrapped in the voids in situ, which cannot escape through the soil surface due to the higher velocity of the water infiltrating from the ground surface. The a w =a d ratio ranged from 1.73 to 2.18, which agrees well with the values reported by Kool and Parker [28], and these values are 2.08 ± 0.46 on average, with values of 1.88 ± 0.40 for undisturbed soil and 2.29 ± 0.47 for compacted soil. Finally, given that the relationship between the volumetric water content and hydraulic conductivity was considered non-hysteretic [28], the K w s values in Table 4 were calculated as K w s ¼ K d h w s À Á after the fitting of h w s . The datasets from one or two drying-wetting cycles were employed to calibrate the hysteretic part of the model. In the case of specimen 2, the model parameters were calculated considering one or two drying-wetting cycles. Thereafter, two R 2 values were determined by comparing all experimental data in terms of the suction (i.e., based on three drying and two wetting branches) to the results obtained from the model calibrated on only one dryingwetting cycle (R 2 = 0.93) and to those obtained from the model calibrated on two drying-wetting cycles (R 2 = 0.92). By comparing the two R 2 values, it was found that (i) considering a wider data set does not enhance the ability of the numerical model to reproduce the experimental data and (ii) the experimental determination of only one dryingwetting cycle greatly reduces the time required for each test.
The sensitivity analysis on the simulation of one cycle required the disturbance of each fitting parameter associated with the main drying and wetting branches. The disturbance was equal to the standard deviation of each parameter presented in Tables 3 and 4. However, h w s only decreased by r h w s , and it would not exceed h d s , while h w s ¼ h d s was adopted as the maximum. The sensitivity of the average matric suction head estimated at the top and bottom of the soil sample is reported as a function of the elapsed time of one cycle, as shown in Fig. 10a-h. The sensitivity associated with the main drying parameters increases over time up to approximately 100 h of testing (the drying phase). The sensitivity to n v initially decreased and then rapidly increased as evaporation proceeded.
A comparison of the dimensionless scaled sensitivities corresponding to the drying phase suggests that the matric suction head measurement has quite a different importance in the estimation of the model parameters. In particular, parameters n v , h r , h s and a d exhibit a more notable relation with the matric suction head measurements. These results are consistent with those obtained in similar soils by Nicotera et al. [39]. Indeed, it is interesting to note that the sensitivity of the modelling parameters of the main drying curve is nonnegligible during the wetting phase but approaches a value of one. This demonstrates their influence on the determination of the parameters corresponding to the main wetting curve. The model sensitivity along the wetting phase to the main drying parameters tends to decrease or remain constant. The disturbed parameter vector leads to deviations in the model-estimated suction values that propagate from the drying phase into the wetting phase. This observation also provides insights into why the adoption of a hysteretic model is important for the correct estimation of the groundwater regime. A disturbed estimation of h d s results in the propagation of any deviations in the estimated values across the entire suction value range. The remaining parameters (n v , h r , and a d ) exhibited a decreasing sensitivity with decreasing suction. The analysis results in Fig. 10g, h reveal that S s;b , which is related to h w s and a w , is null or very low up to approximately 200 h of testing, as expected because this period corresponded to the drying phase and matric suction equalization before the first addition of water. After 200 h, the matric suction was sensitive to variations in h w s and a w because S s;b exceeded 1, as shown in Fig. 10f, g, respectively.
Finally, a comparison of the main loops obtained for each specimen via inverse analysis and the respective experimental wetting and drying cycles is shown in Fig. 8: in all cases, the main loop contains all experimental data.

Considerations on the representativeness of the model parameters
The reproducibility and representativeness of the model parameters were tested by comparing the experimental data of the drying-wetting cycles for specimen no. 2 to the numerical simulations performed by considering the parameters obtained from another test (specimen no. 10). The variability within the tested soil affects the simulation results and reveals how poorly estimated parameters can compromise the simulation accuracy. The complete sequence of the drying and wetting phases was simulated for specimen no. 2 by applying a water flux to the upper boundary equal to that registered during this test. The initial suction profile within the soil sample was equal to the measured initial matric suction distribution in specimen no. The values of all model parameters of specimen no. 2 were considered in the simulation of the test, i.e., this case corresponds to the best-fit solution.
The simulation of the test on specimen 2 is plotted in the WRC plane in Fig. 11. The adoption of the correct parameters associated with the main drying curve enabled a good fitting of the experimental data (case 2, Fig. 11b), which was much better than that with the reference parameters (case 1, Fig. 11a). This result demonstrates that the correct estimation of the main drying curve is necessary to correctly characterize the soil hydraulic behaviour.
In case 2, adopting the fitting parameters of the wetting phase of specimen no. 10 caused an increase in the differences between the simulation and experimental data, especially in the range of 1-10 kPa during the wetting phase, which propagated into the following drying phase.
The parameter combination of case 3 did not produce a notable improvement in the overall model response, but adopting the correct value of h w s resulted in better estimations in the low-suction range: the simulated scanning wetting path (wet 1) is closer to the experimental data. In specimen no. 2, the h w s =h d s ratio is equal to 0.97, while in specimen no. 10, the ratio is equal to 1 because no air entrapment occurred.
The simulation resulting from the best-fit solution (i.e., case 4) provided excellent estimates of the matric suction along the scanning wetting paths (Fig. 11d) as a consequence of the adoption of a correct value for parameter a w . Indeed, specimen 10 presented the highest value of a w , which led to a wider hysteresis loop. Parameter a w , controlling the hysteresis amplitude, was required to capture the soil behaviour along the scanning wetting paths. Although the a w =a d ratios of soil specimens 2 and 10 are very similar, parameter a w for specimen 10 is much higher than that for sample 2. Considering all tests performed, the a w =a d ratio varied between 1.73 and 2.39 (Table 4) with an average of 2.18. Therefore, as a rule of thumb, once the main drying curve is experimentally obtained, the main wetting curve can be inferred by assuming a h w s =n value of 0.75 and a a w =a d value equal to the mean ratio from the literature, e.g., that for Campanian pyroclastic soils is 2.18. In general, to expeditiously evaluate a w , typical a w =a d values can be established for each soil.

Conclusions
Undisturbed soil samples of Campanian pyroclastic soils were subjected to a sequence of laboratory tests to investigate their hydraulic properties within a certain matric suction range (i.e., from 0.1 to 1000 kPa), which extends from saturation up to the residual water content in sandy and silty soils. An experimental procedure which consists of consecutively applying different techniques to the same soil specimen is proposed and used in this research as an extension of the method recommended by Nicotera et al. [39]. This approach allows to study the soil behaviour upon wetting and investigate the hysteretic effects. The hysteretic K-P-S model of Parker and Lenhard proved suitable to fit the experimental data via inverse analysis of the boundary value problem consisting of the wetting and drying phases measured in ku-Pf apparatus. In particular, the parameters of the van Genuchten equation for the main drying and main wetting curves were separately fitted. The inverse analysis results revealed that the experimental data collected along main drying curve and one drying-wetting cycle are sufficient to obtain excellent estimations with the adopted model, which can suitably reproduce the soil response along the drying and wetting cycles from any reversal point. The experimental procedure captured the effects of the entrapped air and the AEV reduction along the wetting branches. Different levels of accuracy are possible, depending on the scope of the application. For example, a faster determination of the parameters is possible because only the direct identification of the parameters associated with the main drying curve is essential for a satisfactory hydraulic characterization. Parameter a w , required to determine the hysteresis amplitude of the scanning paths, is closely related to a d and can be obtained through its correlation with this parameter. In this case, the hydraulic response of the model under a low matric suction can be improved through the direct determination of h w s . Otherwise, this parameter may also be obtained through an empirical correlation with h d s or with the soil porosity, n. In fact, the a w =a d and h w s =n ratios both remain relatively homogenous within the same lithotype.
The fully calibrated hydraulic hysteretic model provided in this work could be adopted to investigate rainfall-induced landslides in unsaturated pyroclastic slope and to set Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons. org/licenses/by/4.0/.