Analysis of landslide triggering conditions in the Sarno area using a physically based model

Rainfall is recognized as a major precursor of many types of slope movements. The technical literature reports both study cases and models of landslides induced by rainfall. Subsurface hydrology has a dominant role since changes in the soil water content significantly affect the soil shear strength. The analytical approaches used are very different, ranging from statistical models to distributed and complete models. The latter take several components into account, including specific site conditions, mechanical, hydraulic and physical soil properties, local seepage conditions, and the contribution of these to soil strength. This paper reports a study using a complete model, named SUSHI (Saturated Unsaturated Simulation for Hillslope Instability), to simulate the role of subsurface hydrology in rain-induced landslides, on a case of great interest both in terms of its complexity and its severity. The landslide-prone area in question is located in Campania (southern Italy), where disastrous mudflows occurred in May 1998. The region has long been affected by rainfallinduced slope instabilities, which often involve large areas and affect many people. The application allows a better understanding of the role of rainfall infiltration and suction changes in the triggering mechanism of the phenomena. These changes must be carefully considered when assessing hazard levels and planning mitigation interventions regarding slope stability.


Introduction
The problems and the damage caused by landslides have become increasingly complex, accounting each year for huge property damage in terms of both direct and indirect costs.Social and economic losses due to landslides can be reduced by effective planning and management.This includes limiting of development in landslide-prone areas, following appropriate construction regulations, using physical measures to prevent or control landslides, and setting up early warning systems.To tackle landslide problems, it is necessary to develop a better understanding of trigger mechanisms, propagation and impact.
As a result of rainfall events and subsequent infiltration, the soil moisture can be significantly changed with a decrease in suction in unsaturated soil layers and/or an increase in pore-water pressure in saturated layers.As a consequence, in these cases, the shear strength can be reduced enough to trigger a failure.
The severity of events is also influenced by the heterogeneity of hydraulic and geotechnical properties.The complex hydrological responses of natural slopes are strongly influenced by infiltration into unsaturated soil, surface runoff, slope-parallel flow through perched aquifers, subsurface flows from upstream areas, the effect of vegetation and flows through fractured bedrock.All these issues affect the predictive ability of the models and sometimes impact on the difficulty of interpreting the results.
In addition, the shear strength contribution from soil suction above the groundwater table is usually ignored if the major portion of the slip surface is below the groundwater table.However, negative pore-water pressures can no longer G. Capparelli and P. Versace: Analysis of landslide triggering conditions in Sarno area be ignored in situations characterized by a deep groundwater table and shallow failure surface (Lu and Godt, 2013).
Analytical approaches differ in terms of the spatial scale range adopted, the available data and the description of the processes occurring in the slope.
Empirical models that directly analyze the rainfall by identifying threshold values are very popular.These values are assumed on the basis of historical data and are shown in an intensity-duration plot as proposed by Caine (1980); they provide a lower limit of rainfall associated with the occurrence of landslides, shallow landslides and debris flows (Guzzetti, 2008).
Other types of rainfall thresholds (Glade, 2000) consider the effect of previous rainfall events as more important than the rainfall recorded on the day the landslide occurred.Usually this type of approach is related to the study of deepseated landslides.
Combined with rainfall forecasting models, both stochastic and meteorological, empirical models are an essential tool to support the prediction of landslides in early warning systems (Sirangelo et al., 2003;Versace et al., 2003;Capparelli and Tiranti, 2010).
However, these models do not provide any information on the hydrological processes involved in a landslide area and do not improve our understanding of landslide dynamics.
Complete models can help in understanding triggering mechanisms since they attempt to reproduce the physical behavior of the processes involved at hillslope scale, using detailed hydrological, hydraulic and geotechnical information (Montgomery and Dietrich, 1994;Wu and Sidle, 1995;Iverson, 2000;Baum et al., 2002;Rigon et al., 2006;Arnone et al., 2011).They develop an analysis over wide areas and usually produce a susceptibility map characterizing the landslide-prone zones according to a stability index.They are generally composed of a hydrological and a geotechnical module.
Hydrological modules differ considerably from each other.The geotechnical module in most cases makes use of the limit equilibrium method under the assumption of an infinite slope.Only in a few cases is more detailed analysis carried out, as in Lu and Godt (2013), where the authors provide a quantitative treatment of rainfall infiltration, effective stress and their coupling and roles in hillslope stability by introducing a unified effective stress modeling framework, linking soil suction to effective stress.
The hydrological approach proposed in SHALSTAB (SHALlow Landslides STABility) (Montgomery and Dietrich, 1994) assumes a constant infiltration rate, neglects soil moisture above the water table, does not take into account the transient response to rainfall and considers the groundwater flow parallel to the slope.These assumptions are too restrictive, for example, when pore-water pressure responds very quickly to transient rainfall and its redistribution has a large component that is normal to the slope.Wu and Sidle (1995) also combined the infinite slope limit equilibrium equation with a subsurface flow model based on the kinematic wave approximation, also taking into account the vegetation root strength.An enhanced version of this model is proposed by Dhakal and Sidle (2004), who investigate the influence of different rainfall characteristics on slope stability.Iverson (2000) developed a flexible framework by modeling a one-dimensional linear diffusion process in saturated soil, using an analytical solution of Richards equation.The model is valid for hydrological modeling in nearly saturated soil.According to this hypothesis, which is used to find an analytical solution to pressure heads, the infiltration capacity is assumed to be equivalent to the saturated hydraulic conductivity, instead of considering it as variable with time during the rainfall event.In addition, Iverson considers the ground surface of the hillslope subject to uniform rainfall.
In order to take into account the variability of rainfall intensity and duration, dynamic or quasi-dynamic models have been introduced.These include TRIGRS (Transient Rainfall Infiltration and Grid-based Regional Slope-stability model) (Baum et al., 2002), which gives a more precise description of slope hydrology but requires a large number of parameters.The model performs transient seepage analysis, using the linearized solution to the Richards equation proposed by Iverson (2000) and extended to the case of impermeable bedrock located at a finite depth.
A more recent version -TRIGRS-unsaturated (Baum et al., 2008) -predicts the pore-water pressure regime in unsaturated/saturated conditions, coupling the simple analytic solution for transient unsaturated infiltration proposed by Srivastava and Yeh (1991) with the original TRIGRS equation.
Most of the aforementioned approaches rely on the restrictive assumption of a steady-state subsurface flow, which can affect the predictive capability of the models both in terms of accuracy and timing of the prediction.All models must always be validated by checking, for example, the accuracy of the simulation with the experimental data available from real cases.The results, in fact, can sometimes be very different when different models are applied to the same event.This is described in Sorbino et al. (2010), who illustrate how, by applying three different physically based models (SHALSTAB, TRIGRS and TRIGRS-unsaturated) to the same set of geoenvironmental cases, different results were obtained.The results reveal the advantages and limitations of each model in landslide forecasting.
These types of spatial, distributed modeling are clearly well-suited for shallow landslides but in deeper landslides their effectiveness is hindered by the higher complexity of the phenomena (van Westen et al., 2003).
In this work a complete model, named SUSHI (Simulation for Saturated Unsaturated Hillslope Instability; Capparelli and Versace, 2011), is applied to a very complex case in order Hydrol.Earth Syst.Sci., 18, 3225-3237, 2014 www.hydrol-earth-syst-sci.net/18/3225/2014/ to improve the understanding of the slope failure mechanism during rainfall infiltration.SUSHI takes into account several components, such as specific site conditions, mechanical, hydraulic and physical soil properties, local seepage conditions and their contribution to soil strength.
It is composed of a hydraulic module to analyze the subsoil water circulation due to the rainfall infiltration under transient conditions and a geotechnical module, which provides the slope stability using the limit equilibrium equation.
The hydraulic process is highlighted by the implementation of a finite difference scheme, which solves the Richards equation used to represent saturated/unsaturated flows within a hillslope.The temporal and spatial distributions of moisture content in the subsurface are calculated, in order to evaluate the different contributions such as the downslope and vertical components of flow in the hillslope caused by unsteady rainfall.
The model was also developed for cases with strongly heterogeneous soils, irregular domains and variable boundary conditions in space and time.After a brief introduction of the model, the paper describes the analysis and results obtained for the volcaniclastic covers of Sarno, where catastrophic mudslides occurred in May 1998.
Mudslides have been analyzed in several papers focusing on the most significant, hydrological and geotechnical features and on the models for triggering mechanisms and landslide propagation.
There is still no commonly accepted interpretation of the triggering mechanism of landslides such as ones occurring in similar geomorphological contexts in Campania.
To summarize, the objectives of the presented application are: -to formulate a likely interpretation of the mudslide triggering mechanism, highlighting the main factors affecting slope stability; -study the rainfall infiltration dynamics, checking whether pore-water pressure, just before the landslide events, reached levels never attained before; -verify the model's ability to analyze a real case, characterized by very marked spatial heterogeneity.
Using a physically based model designed and created ad hoc is mandatory for such purposes, as unlike many commercial models, it allows interaction with the source code in such a way as to adapt it to the specific case by modifying the spatial and temporal resolution as well as the parametrization of the various processes and variables analyzed.It is thus possible to make replacements, to verify the accuracy of the results and to improve its performance by adapting it to specific cases.2 Study area

Description of the study site
The case study proposed in the paper is located in Campania (Italy), where catastrophic flowslides and debris flows in pyroclastic soils are frequent.A brief list of some recent events is reported in Table 1, which also includes information on the size of the landslide.Pizzo d'Alvano is a NW-SE oriented morphological structure, consisting of a sequence of limestone, dolomitic limestone and, subordinately, marly limestone dating from the Lower to Upper Cretaceous age.The slopes are mantled by very loose pyroclastic soils which are the result of the explosive activity of the Somma-Vesuvius volcanic complex, both as primary airfall deposits and volcaniclastic deposits, according to the mode of transport and deposition (Rolandi, 1997).Airfall deposits were dispersed along directions ranging from N-NE to S-SE, according to the prevailing wind direction, and covered a wide area reaching distances of up to 50 km.Pumiceous and ashy deposits belonging to at least five different eruptions have been identified.From the oldest to the youngest, they include: Ottaviano (8000 years BC; E-NE dispersion direction), Avellino (3800 years BC; E-NE dispersion direction), AD 79 (E-SE dispersion direction), AD 472 (N-NE dispersion direction) and AD 1631 (N-NE dispersion direction).The deposits are affected by pedogenetic processes determining paleosoil horizons during resting phases of the volcanic activity.The total thickness of the pyroclastic covers in these areas ranges between a few decimeters, along the steepest slopes, to 10 m, near the uppermost flat areas.The general structure of the soil progressively adapts itself to the morphology of the calcareous substratum thus showing complex and variable geometries (Rolandi, 1997).

Description of the landslide events
On 5 May 1998, a huge number of mudflows were triggered on the slopes of the Pizzo d'Alvano Massif (Fig. 1), involving an extension area of around 60 km 2 , a volume of 2 000 000 m 3 (40 % derived from the eroded materials along the channels) and leading to 159 victims and huge damage to the towns of Sarno, Quindici, Siano and Bracigliano.
These landslides were classified as very rapid to extremely rapid soil slip/debris flows (Ellen and Fleming, 1987), which traveled downslope and then propagated in highly urbanized areas.
A characteristic element is represented by the run-out distances, which ranged from a few hundred meters up to distances of over 2 km (Revellino et al., 2004) and speeds which, at the toe of channels, were estimated to be in the range of about 5-20 m s −1 .
Many similar phenomena have afflicted various other parts of the world, (Japan in 1985, California in 1973and 2005, Brazil in 1967, Venezuela in 1999), in some cases involving similar pyroclastic soils.
Although the triggering mechanisms are different and sometimes the soils involved are not always similar, the common feature seems to be the presence of particles with a high porosity and a very low degree of cementation; these change suddenly due to the action of an external agent (such as an earthquake or more often a rainfall event), which produces a rapid increase in pore-water pressure.What made the May 1998 landslides unique and what made the events even more tragic, was that they all occurred at the same time over a large area, and an extremely high amount of material was involved.
The flowslides of Campania have been analyzed in several papers which indicate the most significant geomorphological, hydrological and geotechnical features of the slopes involved, and propose models for the triggering mechanisms and propagation of landslides.Cascini et al. (2008) argue that the instabilities in Sarno were caused by a combined effect of water infiltrated into the surface layers plus water from a temporary spring in the bedrock.Calcaterra et al. (2000) discuss the role played by groundwater circulation inside both the pyroclastic deposits and the karst cavities of the underlying limestone bedrock.
Soil water circulation is important due to the typical stratification of the pyroclastic covers involved, where one or more layers of pumice, with high permeability at saturation, and layers of paleosoils, with lower saturated permeability, are present.
When persistent rainfall events occur, it encourages subsurface runoff, which may predispose the slope to instability in limited areas.
The role played by the possible interaction of the unsaturated cover with the underlying groundwater, through an impervious soil-bedrock interface, has also been analyzed, by Greco et al. (2013), for a slope similar to those in Sarno.
Vegetation and plant roots and their possible relations to slope instability have also been analyzed.Some authors emphasize the highly dynamic nature of the soil-vegetation as a whole, where the hydrological processes can be greatly affected by the dynamics of the vegetation.Changes in vegetation cover can produce equally rapid effects on the soil and water regime (Mazzoleni et al., 1998).

Description of field surveys and rainfall events
In the years following the landslide events, many field surveys have been carried out in order to assess hydraulic and geotechnical soil characteristics of the site.
In order to assess the influence of soil suction on the triggering mechanism, suction measurements were performed along the Tuostolo Basin (Sarno area), very close to areas that collapsed in May 1998 (Cascini and Sorbino, 2004), using quick-draw portable tensiometers and jet fill tensiometers.These measurements were taken at three sites (Fig. 2), at different depths from the ground surface.Site 1 was located in an area not affected by the landslides in 1998; sites 2 and 3 were located in landslide source areas.
A significant spatial and temporal variability of soil suction can be noted (Fig. 3), essentially related to the differences between the sites, the depths at which the measurements were carried out and also to local factors which induce changes at the end of the dry season when the acquired data show very high suction levels (up to 65 kPa).
The data in Fig. 3 also show that slope equilibrium is guaranteed by cohesion rather than by friction.In fact, in situ and laboratory investigations show that the pyroclastic materials involved have friction angles of between 32 and 38 • and effective cohesion ranging from 0 (excluding reworked pumices) to 4-5 kPa (ashes).The pyroclastic deposits covering the affected areas often rest on slopes with an inclination greater than 40 • , much larger than their friction angle.In such geomorphological conditions, soil suction, which increases soil shear strength, is a major contributor to slope stability.
These considerations support the interest in the suction assessment in these pyroclastic soils as possibly the main predisposing cause.The availability of models able to simulate the circulation of water in these complex terrains can provide useful tools to better understand these phenomena.
The rainfall data were collected by the rain gauge of Santa Maria La Foce, around 2 km from the triggering landslide area, which in 1998 recorded a total rainfall value of 900 mm (Fig. 4).
This station is located at 192 m a.s.l., lower than the landslide source areas (700 m a.s.l.).It is therefore likely that the rainfall along the landslide source areas could be greater than that measured by the rain gauge located at a lower altitude.Available data reveal that, the rainfall occurring in May 1998 was not extreme, in fact it was characterized by a return period lower than five years, but the period when it occurred was unusual.In fact, the monthly and maximum daily values of rainfall height observed in April and May 1998 are significantly higher than the mean rainfall and than the maximum daily rainfall values over the period 1967-1997 (Table 2).

SUSHI model framework
The model combines two modules: HydroSUSHI, which studies subsoil water circulation, and GeoSUSHI, which assesses slope stability.

G. Capparelli and P. Versace: Analysis of landslide triggering conditions in Sarno area
HydroSUSHI analyzes subsoil water circulation in a spatial two-dimensional domain which can be characterized by irregular soil stratigraphy.
Concerning the GeoSUSHI module, a stability analysis is performed to better understand how negative pore-water pressures (or matric suction) increase the shear strength of the soil.
The infiltration analysis is carried out using the Richards equation (1931), expressed in terms of water potential to facilitate applications on layered soils and transient flow regimes for both saturated and unsaturated conditions.
By adopting a Cartesian orthogonal reference system Oxz, with z axis positive downwards, the governing differential equation is where K (ψ) [L/T] is the hydraulic conductivity which depends on pressure head ψ [L] for unsaturated soils (neglecting soil anisotropy).The right-hand side of Eq. ( 1) is written to simulate water flows in both unsaturated and saturated zones, thus avoiding the use of different algorithms for the resolution of parabolic and elliptic equations.C (ψ) = dϑ/dψ [L −1 ] is the specific soil water capacity in the unsaturated zone, which represents the rate at which soil absorbs or releases water when there is a change in pressure head; S S [L −1 ] is the specific volumetric storage, accounting for soil deformation, which is a feature that not many models have; effective saturation S e [ψ] = (ϑ − ϑ r ) / (ϑ s − ϑ r ), where θ is the water content, θ s is the porosity and θ r is the residual water content, which can be computed using the soil water retention curve.The saturated flow equation is simply a special case of the Richards equation, in which the conductivity and storage terms are not functions of pressure head.This module was upgraded by incorporating a method to describe the evapotranspiration process although this component usually produces only secondary effects when slope mobilizations occur during very rainy periods (Capparelli and Versace, 2011).Since the case study here proposed took place during a very rainy spring, the effects related to evapotranspiration were neglected.

Model structure
In the HydroSUSHI module the finite differences method (FDM) with a fully implicit method is adopted.It is well known that the Richards equation only allows analytical solutions where simplifying hypotheses and/or particular boundary conditions are introduced (Iverson, 2000;Srivastava and Yeh, 1991).
Finite difference algorithms which deal with either variably saturated or fully unsaturated conditions have been proposed by Freeze (1978) and Vauclin et al. (1979).In this approach, the continuous problem domain is discretized so that the dependent variables are considered to exist only at discrete points.Figure 5 shows an example of the spatial discretization created by a regular mesh ( x; z).The size to be assigned to x and z should be selected on the basis on the complexity of the stratigraphy.The layers need to be faithfully reproduced, so as to guarantee a realistic representation of water flow exchanges.
With reference to the generic node with coordinates x = x 0 + i x, z = z 0 + j z according to the finite difference method, Eq. ( 1) can be written as where C SU = [C (ψ) + S e (ψ) S s ], the subscripts i ± 1/2, j and i, j ± 1/2 indicate quantities evaluated at the spatial coordinates (x 0 + (i ± 1/2) x, z 0 + j z) and (x 0 + i x, z 0 + (j ± 1/2) z), t is the time step, and the superscripts (k) and (k + 1) indicate quantities referring to time instants t = t 0 + k t and t = t 0 + (k + 1) t.
To solve Eq. ( 1), boundary conditions along the edges of the integration domain must be specified.A general form of the boundary conditions for the considered PDE (partial differential equation) can be written as (McCord, 1991) In particular (9) The value K (x, t) will depend on the values of ψ (x, t) at the point x at time t and on the nature of the K (ψ) curve for the surface soil at x.As previously mentioned, in the GeoSUSHI module the stability analysis is carried out by evaluating the influence of negative pore-water pressures as well.
It may be a reasonable assumption to ignore negative porewater pressures for many situations where most of the slip surface is below the groundwater table.However, for situations where the groundwater table is deep or where there is concern over the possibility of shallow failure surface, negative pore-water pressures cannot be ignored.
The procedure adopted here is an extension of conventional limit equilibrium methods adapted to unsaturated soils as suggested by Fredlund and Rahardjo (1993).The shear strength of an unsaturated soil can be formulated in terms of independent stress state variables (σ − u a ) and (u a − u w ) as follows: where the subscripts f indicate quantities evaluated on the failure plane at failure, and τ ff is shear stress, c is effective cohesion, (σ n − u a ) f net normal stress state, u a pore-air pressure, ϕ effective friction angle, (u a − u w ) matric suction and ϕ b angle indicating rate of increase in shear strength relative to the matric suction.This last term is evaluated using the expression proposed by Vanapalli et al. (1996).Equation ( 10) is an extension of the shear strength equation for a saturated soil.As the soil approaches saturation, the pore-water pressure, u w , approaches the pore-air pressure u a and matric suction (u a − u w ) goes to zero.The general limit equilibrium method (GLE) provides a general theory whereas other methods can be viewed as special cases.It is well known that the elements used in the GLE method to derive the safety factor (FS) are the summation of forces in two directions and of the moments about a generic point Figure 6.Picture of the case study area with some topographic details. (Fredlund and Rahardjo, 1993).FS is defined as the factor by which the shear strength parameters must be reduced in order to bring the soil mass into a state of limiting equilibrium along the assumed slip surface.
Here the calculations for evaluating the stability of a slope are performed by dividing the soil mass above the slip surface into vertical slices.The mobilized shear force at the base of a slice can be written using the shear strength for an unsaturated soil: where S m is the shear force mobilized at the base of the slice; and β is the sloping distance across the base of a slice.

Input data and slope scheme
The SUSHI model was applied to the mudflow occurring in the Tuostolo Basin, highlighted by the red square in Fig. 1, which destroyed the village of Sarno.
The actual geometry of the mudflow is the result of the coalescence of further landslides that took place over a period of 6-8 h.
The landslide mobilized a volume of about 92 000 m 3 of volcaniclastic materials resting over carbonate bedrock, including the eroded material within the channel.It developed from an altitude of about 725 m to a subvertical limestone wall (a morphological frame) situated at an altitude of about 500 m (Fig. 6).
Most of the landslide occurred in May 1998.It started from the point where morphological discontinuities are located, also represented by topographic variations or anthropogenic discontinuities such as tracks.
To define the dynamics of the water circulation in the subsoil, the solution requires a description of the investigated domain, the soil water characteristic curves, the permeability functions, the mechanical properties of the involved soils, as well as the boundary and initial conditions.Surveys and studies carried out by also using information available in the literature indicated the presence of alternating layers of pumice with a composition and thickness related to the characteristics of the eruptions and to the distance from the eruptive centers.
This sequence comprises both primary airfall and volcaniclastic deposits.The primary deposits are composed of alternating layers of pumice, with interbedded paleosoils.At the base of this sequence, above the bedrock, there is a layer of dark-red clayey ashy soil (regolith) with rare limestone fragments.
By using the available topographic maps showing the top surface of the ground before the events, the soil cover thickness and their distributions were obtained.
At the main scarp, the average thickness of the pyroclastic cover was about 4 m.From top to bottom under a topsoil formed by humified ashes including roots and organic matter (about 90 cm thick), the following layers were identified: (A) an upper layer (60 cm) of coarse pumice; (B) a layer (70 cm) of paleosoil; (C) a horizon (60 cm) of finer pumices; (D) a layer (80 cm) of paleosoil; (E) a bottom layer (40 cm) of weathered dark-red clayey ash in contact with the fractured limestone bedrock (Fig. 7).
To determine the mechanical and hydraulic properties of the cover involved, undisturbed specimens were collected, both from the investigated area and from other triggering areas on the Pizzo d'Alvano slopes.Table 3 reports the mean values of the physical properties of the various materials.
The hydraulic properties of the ashy soils in saturated conditions were investigated by conventional permeameter tests.In the unsaturated conditions, suction-controlled oedometer was used.
The experimental data were fitted by the expression proposed by van Genuchten and Nielsen (1985).The values of the parameters estimated for the various layers are given in Table 4.
As shown in Fig. 8, the SWRC (soil water retention curve) are extremely variable, butt in all cases they are typical of coarse soils with low air-entry pressure, low residual water content and steep slope within the transition zone.
Due to a lack of direct measurements, representative values have been considered as the specific storage S S , ranging from 3 × 10 −3 (m −1 ) for plastic clay to 1.3 × 10 −4 (m −1 ) for dense sand (Freeze and Cherry, 1979).
The variable boundary conditions were provided using both Dirichlet and Neumann conditions.At the top (i.e., on the ground surface) a flux boundary condition equal to the rainfall infiltration capacity was applied; the runs enabled the infiltration rate to be defined step by step for each node of the domain.At the bottom (i.e., at the contact between the pyroclastic cover and the bedrock) no flux was imposed, since the bedrock was assumed to be impervious.
Similarly, for the upslope, left side, a Neumann condition of no water flow was fixed, since the morphology of the analyzed area reasonably leads to the hypothesis of coincidence between the surface and underground watershed.These contributions of fluxes from upstream may be assumed to be equal to zero.For the right-hand side of the downslope, along the morphological frames, two different boundary conditions were imposed using a Neumann or a Dirichlet condition depending on whether or not saturation occurred.
For the slope section shown in Fig. 7 the mesh was made up of 130 000 nodes, using regular quadrilaterals with lengths and heights, respectively, equal to x = 0.20 m and z = 0.05 m.
The initial conditions were defined in a non-arbitrary way, thanks to the data provided by the tensiometers, which were located, as previously mentioned, very close to the selected study area.Such information was very useful for setting the initial conditions.Constant distribution suction throughout the domain was initially hypothesized by selecting the following values: ψ (x, z; t = 0) = 3; 4; 5; 8; 10; 14 [kPa]. (12) By starting a simulation with no rain, a warm-up was performed for each of these values, in order to redistribute the water content all over the domain.The warm-up was stopped when the standard deviation values for each spatial node were less than 10 −5 [m].The distribution thus obtained was considered as being representative of an equilibrium condition.
The profiles obtained were compared with the available in situ evidence recorded by tensiometers at the end of summer, which were in considerable agreement with the warmup results.By comparing these profiles, a strong similarity was evident with the distribution performed with (x, z; t = 0) = 6 kPA.This pore-water pressure distribution was set as the initial condition for simulating the evolution between 1 October 1997 and 5 May 1998.

Results and discussion
The period analyzed (from 1 October 1997 to 5 May 1998) was characterized by a total rainfall of 891 mm, with greater values of rainfall intensity having occurred between the end of October and December 1997.
Given the large extension of the domain investigated, the results presented refer to the conditions reached in two zones considered as representative of the domain: one in the upslope part, at Z = 720 m a.s.l.(hereafter referred to as "section A"), the other at the toe of the slope, at Z = 520 m a.s.l., ("section B").
The temporal variations in suction profiles are given in Figs. 9 and 10 for sections A and B, respectively.
For the sake of clarity, in the two figures, the various computed suction profiles are not all plotted on the same graph.The upper panel shows the first three computed soil profiles, for around the end of October, November and December 1997.The middle panel shows another two profiles, evaluated around the end of January and February 1998.To let the reader compare the new profiles with the old ones, the latter are also plotted in the same panel in light gray.Finally, the bottom panel shows the suction profile computed on 3, 4 and 5 May 1998, during the rainfall event which led to a catastrophic flowslide at the investigated slope.Also in this case, the previous suction profiles are plotted in light gray to help the reader follow the evolution of soil suction with time.
Figure 9 shows that no water table is ever computed along the upslope area (section A), since the values of the pressure head are always negative.This situation is consistent with the morphological characteristics of the investigated site, where the steepness of the slope prevents water from accumulating and lets it move through the layers towards the foot of the slope.
The situation is very different at the toe of the slope (section B, Fig. 10), where the lower layers soon reached saturated conditions, already in November 1997, much before May 1998.From the former date onwards, the lower layers remained saturated for the entire simulation period.
The results of the analysis do not lead to a clear interpretation of the triggering phenomenon, but suggest that the  However, a big difference lies in the fact that, on 5 May 1998, the calculated vertical suction profiles of water content show high values, not far from saturation, in the upper layers as well.
The wetting front generated in early May by the rainfall infiltration did not reach the lower layers and did not contribute to a further increase in pore-water pressure at these depths but led to the shallowest layers becoming close to saturation.
This result is even clearer when analyzing the evolution of the pore pressure distribution along the slope throughout the considered period.Figure 11 shows the pore-water pressures evaluated at 3 and 0.7 m below the ground surface.The first depth is representative of the behavior of a relatively deep layer, which reaches saturation during the first months of the rainy season, while the second represents the conditions of the upper layers well.
In the lower layers (Fig. 11a), the pressure levels remain approximately the same with the rainfall in late April and early May, while it increases sharply in the upper layers (Fig. 11b) from 4 to 5 May.
Regarding the computed pore pressure, a slope stability analysis was carried out to simulate failure conditions and their relation to increasing soil water content.For several depths below the ground surface (0.3, 0.7, 1.8, 2.1, 2.9, 3.1 and 3.8 m), the average pore pressure was calculated, and the corresponding FS was estimated under the hypothesis of an infinite slope.This approach is the simplest limit equilibrium method and gives reliable results for slides where the longitudinal dimension prevails over the depth of the landslide, as with the landslide analyzed here.
The plots in Fig. 12 provide the time history of the simulated FS values, which help in understanding the evolution of the slope stability conditions.
The values in the lower layers are always indicative of stability; lower values, but nevertheless above 1, are due to the greater thickness of the soil cover and to higher values of  pore-water pressure.By contrast, in the more superficial layers the trends are more variable and show, at a depth of 0.7 m, a decrease in FS to 0.98, corresponding to instability conditions on 5 May, 1998.

Conclusions
Landslide prediction is a problem that has concerned the general population for a long time and more recently also the scientific community.The ability to predict landslide onset is closely related to the ability to understand the dynamics of the underlying processes affecting the event.
The application reported here confirms how it is possible to approach the interpretation of a natural but complex and difficult event such as the triggering of a landslide, via simulation with physically based mathematical models.
The analysis outlined a possible interpretation of the Sarno landslides that could explain the triggering mechanisms occurring in May 1998.The model's output suggests that the saturation of deep layers was not the only reason for slope instability, and that the reduction in suction across the shallowest layers may have been the actual cause of the triggering mechanism.
Monthly rainfalls occurring in April and May 1998, though not exceptional, were quite unusual for that period of the year.The rain occurring at the end of April caused the superficial layers to reach a high water content and, consequently, an increase in hydraulic conductivity.Thus, the heavy rainfall in the early days of May produced a steep wetting front which did not reach the deeper layers and considerably reduced the shear strength of the upper soil layers.
This interpretation is partly confirmed by the slope stability analysis, which reveals high values of the safety factor in the lower layers and values close to one in the upper layers with regard to the landslide date.
New applications to other cases occurring after May 1998 would certainly better define the critical conditions and could provide useful information for a possible early warning system.
Further analyses should be carried out in order to better evaluate the influence of the bedrock, of the road cuts located in the upper zone of the triggering areas and other factors that could have helped trigger the landslide.
By pursuing the overall objectives, the SUSHI model confirmed its capacity to (i) simulate, with sufficient detail, the phenomena induced by rainfall in soils characterized by complex stratigraphy and hydraulic properties; (ii) achieve reliable numerical solutions even when running very long simulations, thus enabling antecedent rainfall to be taken into account, and (iii) provide a convincing interpretation of the phenomenon, without the need to introduce other hydrological forcing.
Edited by: R. Greco

Figure 1 .
Figure 1.Overview of the Pizzo d'Alvano Massif and the area affected by the May 1998 mudflows.The red square delimits the area where the SUSHI model has been applied.

Figure 2 .
Figure 2. Topographic map and sites where suction measurements were performed.

Figure 3 .
Figure 3. Suction trends recorded and daily rainfall measured by Santa Maria la Foce rain gauge at 0.20, 1.00 and 1.60 m from topsoil.

Figure 4 .
Figure 4. Comparison of daily and cumulative rainfall data recorded at Sarno rain gauge in 1998.
where α (ζ ), β (ζ ) and B (ζ, t) are given functions evaluated on the boundary ∂G, the expression ∂/∂n is the normal derivative operator and ζ is the spatial local vector.By applying this general formulation for water flow modeling, the boundary conditions are as in Fig.5: along the basal impermeable boundary BC and vertical boundary AB a Neumann condition is considered, with flux equal to zero, then α (ζ ) = 0, β (ζ ) = K (ζ ) and B (ζ, t) = q (ζ, t) = 0.In terms of total hydraulic head h = ψ − z: downslope side DC we take into account the influence of the increasing subsurface flow, by considering both unsaturated and saturated layers.This is computed by adopting boundary conditions moving from the Neumann to Dirichlet conditions, with specified flux or pressure head, respectively, α (ζ ) = 1, β (ζ ) = 0 and B (ζ, t) = h (ζ, t).boundary AD a time-dependent rainfall intensity r [L/T] is applied.The boundary condition can be stated by considering the infiltration rate I (x, t) as K (x, t) ∂h (x, t) ∂z AD = B (ζ, t) = I (x, t) .

Figure 5 .
Figure 5. Nodal network implemented for development of FDM equation.

Figure 7 .
Figure 7. Geometric and stratigraphic characterization of the investigated slope.
[m s −1 ] 3.2 × 10 −5 1.0 × 10 −Parameters of van Genuchten and Nielsen model.The values of the bubbling pressure, or air-entry tension, ψ b , were determined through the graphic method proposed byFredlund and Xing (1994).Topsoiln = 1.6; m = 0.38 ψ b = 1.65 (kPa Layer A n = 1.71; m = 0.42 ψ b = 0.2 (kPa) Layer B n = 1.66; m = 0.40 ψ b = 2.5 (kPa) Layer C n = 1.8; m = 0.44 ψ b = 0.3 (kPa) Layer D n = 1.9; m = 0.47 ψ b = 2.5 (kPa) Layer E n = 2; m = 0.50 ψ b = 2.7 (kPa)saturation of the underlying layers at the toe of the slope was clearly not the cause of the instability of the slope, although it contributed to it.The suction values achieved on 5 May 1998 in the lower layers of the soil profile had been already predicted by the model for previous periods.

Figure 9 .
Figure 9. Computed suction profile for the upslope section (Section A).

Figure 11 .
Figure 11.Pore-water pressures performed at (a) 3 m and (b) 0.7 m below the ground surface.

Figure 12 .
Figure 12.Slope safety factor depending on pore-water pressures and soil mechanism properties at different depths from the ground surface.

Table 1 .
Features of some recent flowslides in the Campania Region(Versace et al., 2009).

Table 2 .
Comparison between mean monthly rainfall computed for the period1964-1997 and the year 1998 (rows a-b); and between maximum daily rainfall for the period1964-1997 and the year 1998 (rows c-d).Some values are in bold because they reflect the situation of the specific period analyzed.This allows to differentiate the values from other reported data.Jan Feb Mar Apr May Jun Jul Aug Sep Oct Nov Dec

Table 3 .
Average values of pyroclastic soil properties.