Modelling suction instabilities in soils at varying degrees of saturation

Wetting paths imparted by the natural environment and/or human activities affect the state of soils in the near-surface, promoting transitions across different regimes of saturation. This paper discusses a set of techniques aimed at quantifying the role of hydrologic processes on the hydro-mechanical stability of soil specimens subjected to saturation events. Emphasis is given to the mechanical conditions leading to coupled flow/deformation instabilities. For this purpose, energy balance arguments for three-phase systems are used to derive second-order work expressions applicable to various regimes of saturation. Controllability analyses are then performed to relate such work input with constitutive singularities that reflect the loss of strength under coupled and/or uncoupled hydro-mechanical forcing. A suction-dependent plastic model is finally used to track the evolution of stability conditions in samples subjected to wetting, thus quantifying the growth of the potential for coupled failure modes upon increasing degree of saturation. These findings are eventually linked with the properties of the field equations that govern pore pressure transients, thus disclosing a conceptual link between the onset of coupled hydro-mechanical failures and the evolution of suction with time. Such results point out that mathematical instabilities caused by a non-linear suction dependent behaviour play an important role in the advanced constitutive and/or numerical tools that are commonly used for the analysis of geomechanical problems in the unsaturated zone, and further stress that the relation between suction transients and soil deformations is a key factor for the interpretation of runaway failures caused by intense saturation events.


Introduction
Humanity faces unprecedented challenges related with population growth, increasing energy demand, and degradation of ecosystems [1].In this context, significant threats come from the deterioration of the climate to a point that vast communities are exposed to growing risks of extreme events [2].To address such crises, a key role is played by geotechnologies [3], i.e. interventions on near-surface and/or underground deposits for the purpose of environmental remediation (e.g., disposal of CO 2 and hazardous wastes) and risk management (e.g., protection of cities and/or coastal areas from sea level rise and other geohazards).Research on the multi-physics of geological systems, as well as on the couplings between mechanical, hydrologic, chemical and biological agents in the natural environment, is therefore a new frontier to be explored to devise innovative technological solutions able to face the abovementioned challenges [4].
From a scientific viewpoint, intriguing problems derive from the inherent complexity of soils and rocks.Geomaterials are indeed porous and heterogeneous media hosting multiple fluids which establish physico-chemical interactions with the solids.Any alteration of their in-situ state affects the strength and deformation properties, thus playing a role in applications such as the underground storage of hazardous products, infrastructure management and geohazard forecasting.In the latter case, the literature offers numerous examples of near-surface instabilities caused by environmental agents.A recent example is the catastrophic event that took place at Oso (USA) in 2014 [5].This case involved a sudden failure of unconsolidated glacial deposits with average inclination lower than 20°.The mobilized masses crossed at high velocity a 1-km floodplain, claiming the lives of 43 people and suggesting that liquefaction mediated by transient flow played a role in the unstable propagation of the mass [6], thus stressing the need to quantify the processes that may lead to such kinds of events.Dramatic evidences of destructive nearsurface failures also derive from other geological settings, such as unsaturated residual soils, volcanic sediments and airborne deposits [7][8][9], i.e. materials in which watersensitivity generates strong hydro-mechanical couplings leading to solid-fluid transitions upon saturation [10].
It can thus be argued that, unless the multi-physical processes controlling soil and rock failure are understood, limited progress can be achieved in assessing the spatial and/or temporal distribution of geologic hazards.In this paper particular emphasis is given to a specific form of multi-physical interactions, which involve changes of the saturation conditions within the pores, thus impacting the loss of unsaturated soil strength, as well as the alteration of suction transients in deformable ground.Specifically, the goal is to stress the role of solid-fluid interactions in failure events in shallow unsaturated slopes, as well as to discuss the role that a loss of mathematical stability may have in the constitutive/numerical tools that are widely used for geomechanical analyses in the unsaturated zone.

Second-order energy input in partially saturated soils
Previous works by Nova and coworkers [11,12] have elucidated the connection between second-order energy input to a solid and the singularity of the principal minors of the constitutive tensor.Such relation has been pivotal to correlate non-positive values of second-order work to the loss of control of loading variables, and hence to the failure of soil specimens.In hydro-mechanical contexts characterized by multiple fluids, the incremental response of a porous medium is controlled by both mechanical and hydraulic agents (e.g., pressure changes, fluid injection or extraction).As a result, additional control variables must be added to the constitutive law, requiring an adaptation of both second-order work input and control theory.Such an extension for unsaturated soils with two continuous fluid phases (which, in usual environmental contexts, coincide with water and air) has been discussed by Buscarnera and Nova [13] and Buscarnera and di Prisco [14], and it is rooted in the following expression of second-order work: where * 1 1 is a stress measure work-conjugate with the strains, ij ij , s s the suction rate, n the porosity, r S the degree of saturation, and a the air density.In most applications, only the first two terms of the second-order work are considered, thus implicitly assuming passive air phase (i.e., constant atmospheric pressure).Nevertheless, cases of relevance for failure events in the near-surface can be envisaged when, upon a substantial increase in degree of saturation, the air loses continuity because of air-entrapment, possibly causing undrained conditions for the air phase and simultaneous growth of both air and water pressures (Fig. 1).In the following, Equation (1) will be used to study coupled failure mechanisms under unsaturated conditions that may imply a simultaneous loss of control of variables associated with the two pore fluids.For this purpose, to express the second-order energy input as a function of variables commonly used in laboratory experiments, the expression of 2  d W is rearranged as follows:

Loss of Air Continuity
where ij is the total stress, w r e eS is the water ratio, and a m is the mass of air in the porous medium.This choice enables expression of the second-order energy input in terms of variables that are directly imposed during soil testing (i.e., total stress, ij , suction, s , and mass of air and/or water, i.e. w e and/or a m ), thus offering a convenient approach to link the energetics of hydromechanical perturbations to the constitutive description of the soil.

Controllability analyses
The second-order work input can be linked to constitutive laws to compute failure indices related with given hydromechanical controls (a strategy that is usually referred to as controllability analysis [12]).Here, various drainage conditions will be considered, with the goal to assess the role of the pore fluids on the onset of failures involving either uncontrolled water pressures or sharp variations of the pressures of both water and air.

Controllability conditions
Let us consider a constitutive law for a three-phasic soil written in incremental form.To define a simple tangent constitutive operator compatible with the second-order energy input, let us restrict the analysis to stress variables reflecting axisymmetric stress conditions, as follows: is a mean skeleton stress increment working for the volumetric strain rate, v v , while q q is the deviatoric stress increment, working for the shear strain rate, s s .The constitutive matrix in While the latter form of coupling can be included with no alterations of the general strategy of analysis [15], for the sake of analytical simplicity here the study will be restricted to uncoupled water retention responses (i.e., 0 ), thus focusing solely on the suction dependence of the solid skeleton (i.e., 0 ).The detection of constitutive singularities involves the analysis of specific loading scenarios.In particular, for a given choice of control variables, (i.e., variables that are directly imposed) and response variables, (i.e., variables work-conjugate to and obtained as an outcome of the soil response) it is possible to state: where X is a control matrix compatible with the second- order work input.It thus follows that a loss of control leading to undefined values of is associated with the following singularity condition [12,16]: To facilitate the selection of the control variables, it is convenient to perform a linear transformation replacing the variables in (3) with the two following sets of workconjugate variables: w a p q u u p p p q q q q q q w u a u u which are consistent with the expression of 2 d W given in (2).This choice enables rewriting of the constitutive laws in (3) as follows: where the following transformation is used: Hereafter, three control scenarios will be used to assess the impact of hydro-mechanical couplings on the constitutive singularities: (i) stress-control combined with drained conditions for both fluids (i.e., , , , w a p q u u p a p u u , , , , , p q u , , ); (ii) stress-control plus undrained conditions for water and drained conditions for air (i.e.p, , ,  p, , , ).Scenario (i) leads to the following singularity condition: which resembles the failure index that was obtained by di Prisco et al. [17] for shear failure in saturated ground and by Buscarnera and di Prisco [18] for unsaturated slopes at constant suction.By contrast, the scenario (ii) produces the following singularity condition:  (11) are corrective terms related with the coupling effects.The expression in (10) resembles that found by Buscarnera and di Prisco [10] for failure under water-undrained paths and/or water injection under the hypothesis of passive air phase.Finally, scenario (iii) provides the following loss of controllability condition: in which the potential for air pressurization contributes to the singularity of X via a term affected by the air-phase compressibility ( GG D ) and by the shearing response of the skeleton ( QQ D ).In the following, such indices will be computed with a constitutive model for unsaturated soils and will be used to assess the potential for coupled failure in volcanic soils involved in flowslide events [9].

Constitutive Models
An elastoplastic constitutive law based on that proposed Buscarnera and Nova [19] will be used for the following analyses.This model includes versatile yield surface and plastic potential expressions proposed by Lagioia et al. [20] to reproduce non-associated flow and a saturationdependent hardening rule based on earlier propositions by Jommi and di Prisco [21], where s p is the internal variable associated with the size of the yield surface, p v and p s are volumetric and shear components of the plastic strain, respectively, and p B , s and sw r are constitutive parameters.In particular, the parameter r sw controls the contraction of the yield surface that would occur during saturation (i.e., when the suction increment is negative), which plays a crucial role on soil instabilities induced by saturation [15].To define the soil E 2016 response for stress states inside the yield surface, the model is completed by an uncoupled hypoelastic law based on the following tangent bulk and shear moduli: which describe power laws linking the soil stiffness to the mean skeleton stress, 1 r w r a p p S u S u , and where 0 K and 0 G are reference values of bulk and shear moduli, respectively, while 0 p is a reference pressure.
Finally, to link the skeleton response to hydrologic processes, the model is integrated by a water retention curve based on the well-known expression proposed by Van Genuchten: as well as by the ideal gas law, which under isothermal conditions can be arranged as follows to relate changes in air density to the air pressure rate:

Calibration Strategy
To illustrate the use of controllability indices, this section makes reference to a class of volcanic soils involved in flowslide events that took place in the Campania region of Italy [9,22].While the site-specific properties of such soils vary across the region, some common features can be identified for the soils of the area, which include the suction-sensitivity of both compressibility and shear strength and their high potential for liquefaction.
To study their behaviour, here a baseline calibration is obtained based on data reported by [23] and [24] for soils retrieved from a specific site of the region (Table I).This calibration will be later applied to reinterpret flow failure processes detected at the laboratory scale on a similar pyroclastic soil from another site [22].
Figure 2 illustrates the calibration of the retention curve parameters with reference to data reported by [23] and [25].The range of retention curve variability reported by [23] for different classes of Campanian soil and stress conditions is also illustrated in the figure .The parameters controlling the shearing response of the Campanian ashes have been assessed on the basis of suction-controlled triaxial compression tests conducted at varying levels of suction and net confinement (Fig. 3).The good agreement between model predictions and data corroborates the validity of the average skeleton stress for the simulation of shear failure in Campanian soils.
To assess the non-associativity of the plastic flow, the liquefaction resistance of saturated samples subjected to undrained shear has been considered [26], identifying a set of parameters capturing accurately the lower bound of the instability line emerging from tests (Fig. 4).Finally, to have a first-order assessment of the ability of the model to reproduce plastic strains due to wetting, the suction-hardening characteristics have been calibrated to mimic the compression behaviour of Campanian soils under saturated and unsaturated conditions, as well as the occurrence of wetting-induced compaction in oedometric experiments conducted at varying levels of vertical net stress (Fig. 5).A satisfactory agreement between model predictions and data was found, corroborating the choice made for the suction-hardening law in Equation ( 13) and its applicability for the study of non-linear phenomena promoted by saturation paths.

Numerical Assessment of Unstable States
The set of parameters calibrated above provides a starting point to inspect the evolution of the state of stability of ashy soils during wetting.An example is the unstable behaviour observed upon suction reduction in pyroclastic soils under constant stresses (Fig. 6).Similar experiments have been reported by Olivares and Damiano [22] on soil samples from a neighbouring area.Such tests were based on an initial suction of 82kPa gradually reduced at a rate of 1.4kPa/hr to ensure drained conditions.The applied loads remained unchanged, thereby imposing a constant stress deviator under decreasing mean skeleton stress.Axial and volumetric strains were tracked throughout the test, along with the degree of saturation.As wetting progressed, a sharp increase in rate of soil deformation was observed in correspondence with a critical suction value close to saturated conditions, which was indicative of unstable behaviour and, ultimately, failure.
The material tested by [22] exhibits different hydraulic retention behaviour from that of the baseline calibration, but similar wetting paths can be simulated using the parameters in Table 1, with stress conditions as shown in Fig. 7b.The resultant volumetric compaction observed during these simulations is presented in Fig. 7a, along with points at which the controllability criteria given in ( 9), ( 10) and ( 12) are violated.It can be noted that the controllability indices computed with these parameters suggest that different hypotheses about the dominating control conditions at the onset of failure yield different times of potential instability.In particular, the first index to vanish upon saturation is found to be that associated with unstable pressurization of both water and air phases (detX WA ), shortly followed by that associated with a pressurization of the water phase (detX W ) under constant air pressure, and finally by the index associated with constant suction failure (detX S ).Most importantly, it must be noted that the saturation levels where these indices vanish are in proximity of states at which spontaneous undrained response and sharp acceleration of the soil deformation characteristics were found [22], thus corroborating the connection between vanishing controllability indices and unstable couplings between fluid flow and deformations.Such links will be explored in greater detail in the next section in light of the field equations that control suction transients in deformable unsaturated ground.

Implications for the stability of suction transients
The analyses discussed above point out the important role of hydro-mechanical couplings on an activation of unstable processes.Most notably, they stress the fact that such failures imply deformation-induced transient alterations of the pressure regime of both fluids.Hence, it is of interest to inspect the implications that constitutive instabilities of unsaturated soils have on the equations that are commonly used to track the spatial and temporal variations of the pressure regimes of deformable deposits.
To achieve this goal, an idealized one-dimensional flow condition is considered (Fig. 8) to derive a mass balance equation.For the sake of simplicity, only pore water flow is considered, and, in agreement with previous analyses, axisymmetric stress states are used to model triaxial compression tests characterized by transient saturation.
These hypotheses rely on simplifying assumptions similar to those used by di Prisco [27] to study pressure transients in saturated soils.In particular, it is assumed that water flow takes place only along the vertical axis of the sample, thus ignoring flow along horizontal planes.In addition, the state of stress, strain and suction is assumed to be uniform until the onset of transient effects, total stress control is hypothesized to be a good representation of the local constitutive response, and the presence of shear stresses at the bases of the sample is neglected.
For vertical flow, the mass balance equation can be written as follows: where z q is the Darcian velocity of water along the vertical axis and the storage of water volume over time is expressed in terms of water ratio, w e .By introducing the Darcy's law in the following form where w k is the suction-dependent permeability for water flow, w is the mass density of water and g is the gravitational acceleration, Equation ( 17) leads to: To write the mass balance equation (19) in terms of water pressures it is possible to use the constitutive equations of a deformable unsaturated soil.For this purpose, relations expressed in terms of tangent plastic compliance are adopted: where hypotheses similar to those adopted for the matrix in (3) are used by neglecting the compliance terms linked with the strain-dependence of the retention curve (i.e., 0 ).It is worth noting that, since only water pressure transients are considered, the effect of the air pressure is neglected in the constitutive law.This implies a hypothesis of passive air (constant air pressure) which implies negligible values for the last term in Equations ( 1) and ( 2).In other words, ( 20) is consistent with a second-order work expression accounting only for soil deformations and suction changes [14].
(a) (b) By using a transformation of variables similar to that used to pass from the variables in (3) to those in (6), it is possible to rewrite (20) in terms of increments of water ratio, as follows:  which was shown by [14] to be equivalent to (1) under constant air pressure and axisymmetric stress conditions.By using (21)  The differential form in Equation ( 24) is convenient to study the ill-posedness of the temporal evolution of a pressure field controlled by diffusive effects [28], and it has been used in early works by Rice [29] and di Prisco [27] to study pore pressure instabilities in fully saturated geomaterials.In particular, while the function b reflects forcing associated with external loading and advective effects linked with inhomogeneous permeability, the sign of the function A in (25a) is related with the spontaneous temporal evolution of the pore pressure field under lack of external forcing (i.e., homogeneous equation).In other words, the sign of A is associated with the parabolicity of the differential form (24), implying that a self-stabilizing diffusive evolution of the excess pore water pressures (and, hence, of suction transients) is possible only when A is positive, i.e. under the condition Violations of the condition above suggests the possibility of "backwards" diffusion processes, according to which small perturbations of the pore pressure field can quickly grow into singularities [27,29].Since ( 26) is a principal minor of the constitutive matrix in (21), it is interesting to explore the connection of this loss of well-posedness of the transient water mass balance equation in terms of controllability.Indeed, by writing the loss of control conditions as a function of compliance terms, it is readily apparent that the loss of control for a loading scenario characterized by , , / 1 w p q e e p p p / 1 p q , , / 1 e / 1 / 1 / 1 p, , p, , , p, , provides: implying that a loss of parabolicity of the water pressure field coincides with the loss of control for waterundrained loading.This finding suggests an additional interpretation of coupled failure modes at the material point level, as it shows that loss of controllability of pore pressures may imply alterations of the differential form of the field equations, thus impacting the robustness of the numerical techniques used for coupled flow/deformation analysis, as well as the interpretation of measurements of

Conclusions
This paper has discussed a set of analytical techniques to detect coupled failure events in unsaturated soil samples, i.e. states at which a simultaneous singularity of hydraulic (e.g., air and/or water pressure) and mechanical variables (e.g., axial and/or volumetric strains) is possible.For this purpose, energy statements defining a second-order work input in three-phase media have been linked to a set of constitutive relations for partially saturated soils, with the purpose to detect singularities of the tangent constitutive operator.A range of hydro-mechanical control conditions has been considered, encompassing constant suction and stress control (drained loading), constant water content under stress control (i.e., water-undrained conditions) and simultaneous lack of drainage for both the water and the air phase under imposed total stresses.It has been shown that the choice of control conditions impacts considerably the singularity of the control matrix, implying vanishing values of determinants that are met at different states during a saturation process.Such possibility was explored in detail by making reference to a particular volcanic soil involved in past rainfall-induced flowslides.In particular, after defining a baseline calibration of the parameters of a hydro-mechanical model reproducing first-order features of the selected soil, a constant shear wetting path was performed to interpret the unstable behaviour witnessed in similar tests reported in the literature.Conditions for loss of control were monitored during the path and points of controllability violation were identified.Most notably, such points of potential instability were detected to occur in proximity of the state at which failure accompanied by an increase of the rates of deformation and pore pressure growth were measured in laboratory experiments.While the unstable mode corresponding to simultaneous air and water-undrained conditions was shown to be the first to be activated upon wetting, it was also shown that waterundrained instability involving unstable pressurization of the water phase was activated shortly afterwards.Both failure modes were found to anticipate the onset of shearfailure at constant suction, suggesting that coupled controllability indices provide advantages to identify failure modes involving simultaneous changes of deformation and pore pressures.To further corroborate such findings, it was shown that loss of controllability conditions can also be associated with the loss of wellposedness of the field equations that govern the mass balance of water flow and, hence, the stability of suction transients in deformable unsaturated ground.As a result, the tools discussed in this paper find multiple applications not only in the detection failure mechanisms relevant for geohazard initiation, but also to assess the robustness of the analytical/numerical techniques that are routinely used to study coupled flow/deformation problems in unsaturated ground, as well as to conceive regularization strategies mitigating the loss of computational stability in presence of highly non-linear constitutive relations.

Figure 1 .
Figure 1.Schematics indicating transition from a continuous air phase characterized by drained conditions (i.e.0 a u 0 a u ) to states characterized by air entrapment and growth of air pressure due to lack of drainage (i.e.0 a m 0 a m (3) is built under the simplification that the pressure-density relation for the air is uncoupled from the other hydro-mechanical variables, thus enabling the use of basic state equations for perfect gases to define the term GG D .The other diagonal terms reflect either bulk or shear stiffness ( PP D and QQ D ) or the water retention behaviour ( HH D ).Cross couplings are also included in form of hydrologic effects on the stress-

Figure 3 .
Figure3.Simulated shear and volumetric responses at different levels of suction and net confinements (data after[23]).

Figure 5 .
Figure 5. Compression behaviour of saturated and unsaturated samples and simulations of wetting-collapse (data after[24]).

Figure 6 .
Figure 6.Schematic of constant shear stress path under suctioncontrolled wetting, with critical suction (s crit ) attained prior to full saturation.

Figure 7 .
Figure 7. Computed volumetric strains during wetting with points of loss of controllability for suction control (detX S ), constrained water mass (detX W ), and constrained air and water mass (detX WA ) (a) and consequences of switching to waterundrained shearing ( 0 w m 0 w m ) or shearing under simultaneous undrained conditions for air and water ( 0 a w m m 0 a w m m a m ) (b).

Figure 8 .
Figure 8. Schematic describing axisymmetric stress conditions in a triaxial compression test (a) and one-dimensional water flow (b).

C
is a modified compliance term defined as: role of isotropic compliance and hydromechanical coupling.It is worth noting that the variables in (21) are compatible with the following form of2

Table 1 .
Baseline model parameters for Campanian ash.
to express the rate of e w in terms of stress and suction increments, and by assuming that under constant air pressure