Liquefaction-induced lateral spreading may result in significant damage and disruption of functionality for structures and Slope Ground System. In this regard, finite-element simulations are increasingly providing a versatile environment in order to assess economical and effective countermeasures. Several systematic bidimensional FEM computations have been conducted to evaluate mitigation strategies under the action of an applied earthquake excitation. The presented study highlights the potential of computations in providing insights for analysis of liquefaction-induced lateral deformations. In the analysis, some specific assumptions are introduced and verified such as a nine-node quadrilateral elements, massive columns of soil with periodic boundary conditions, and a Lysmer-Kuhlemeyer dashpot used to model the finite rigidity of the underlying elastic medium. Moreover, the study aims to systematically explore the effectiveness of densification as a countermeasure and then evaluate the best extension comparing two scenarios.

1. Introduction

Lateral spreading refers to the development of large horizontal ground displacements due to earthquake-induced liquefaction. This phenomenon may result in significant damage and considerable replacement costs for existing buildings and civil engineering structures (quay walls, bridge piers, etc.) since it imposes notable lateral loads and may lead to widespread failures. Such adverse response is documented during several seismic events, such as the earthquakes of Niigata, Japan (1964, [14]), Dagupan City, Philippines (1990, [58]), Chi-Chi, Kocaeli, Turkey (1999, [9]), Taiwan (1999, [10]), and recent Tohoku, Japan (2011, [11, 12]). Particularly important is the dynamic slope stability under liquefaction lateral deformation into areas where potential structures can be of interest to landslides, such as dams and river or pier banks, especially if previously predisposing phenomena have acted in the static stress field [13, 14]. In this regards, several ground remediations have been developing, such as vibroreplacement [15], solidification (cementation) [16, 17], gravel drains, or stone columns [1821].

This paper aims to assess the seismic reliability and evaluate the ground improvement method of densification for an Italian real case study of a pier founded on a partially submerged layered ground slope strongly vulnerable to liquefaction. The earthquake response of the Ground System is simulated with a two-dimensional, advanced, nonlinear finite element model adopting the open-source computational platform OpenSees [22]. The performance of structural and geotechnical systems subjected to static and seismic loadings can be simulated thanks to the platform implementing a framework for saturated soil response as a two-phase materials, following the U-P formulation of Chan [23] and Zienkiewicz et al. [24].

The present study may be viewed as a further development of earlier and ongoing efforts [2528] to generate appropriate numerical models for simulation prediction of liquefaction-induced ground response using OpenSees. In particular, the main numerical challenge was modelling the boundary conditions to reproduce the elastic half-space under the soil system.

In the following sections, the employed computational formulation is described. The computational platform for slope-ground analysis is then presented focusing on the adopted materials, boundary conditions, and analysis assumptions in Section 3. Results of the conducted analysis for the free-field configuration are presented as well as the two scenarios countermeasures’ response, respectively, in Sections 4 and 5. Finally, insight and conclusions based on the effectiveness of the proposed remediation are drawn.

2. Computational Formulation

All simulations were conducted using the open-source computational platform OpenSees [22]. This platform allows for developing applications to simulate the performance of structural and geotechnical systems subjected to static and seismic loadings. Implemented in OpenSees [2528] is an analysis framework for saturated soil response as a two-phase material following the U-P formulation of Chan [23] and Zienkiewicz et al. [24], where U is displacement of the soil skeleton and P is pore pressure. This implementation is based on the following assumptions: (i)small deformation and rotation;(ii)solid and fluid density constant in both time and space;(iii)porosity locally homogeneous and constant with time;(iv)soil grains incompressible; (v)solid and fluid phases equally accelerated.

The soil constitutive model, Figures 1 and 2 implemented in OpenSees [2529], is based on the multisurface plasticity theory for frictional cohesionless soils proposed by Prevost [30], where is the effective confining pressure, the octahedral shear stress, and the octahedral shear strain. In particular, this constitutive model was developed for simulating the characteristics of cyclic mobility observed in saturated medium to dense cohesionless soil response [2528]. Within a multisurface plasticity framework, the model incorporates shear-induced contractive, perfectly plastic and dilative response phases implemented through an appropriate nonassociative flow rule motivated by experimental observations as to capture the involved phenomena. Emphasis is placed on accurately reproducing the development and accumulation of shear deformations. The hardening rule was also introduced to enhance numerical robustness and to increase efficiency. Finally, a model calibration procedure based on monotonic and cyclic laboratory sample test data was conducted.

3. Computational Platform for Slope: Ground Analysis

The 2D Ground System is a typical slope 1000 m long, 60–130 m high as shown in Figure 3. Such model is built with several layers of cohesive and cohesionless material, mainly of silty medium to fine sands (S), pliocene clays (C) as background layer, and a superficial layers of gravelly sands (G) located near the ground surface. OpenSees may implement a wise number of soil models including multiyield surface cohesionless (Drucker-Prager cone model) and cohesive (Mises or J2) ones. In order to model the different layers of the problem two models were adopted: (i)Pressure-Independent Multiyield (Table 1) for pliocene clays (C) and(ii)Pressure-Dependent Multiyield02 (Tables 2 and 3) for silty medium to fine sands (S).

Model implementation is based on several computational assumptions (regarding finite elements, boundary conditions, and analysis), shown in next sections.

4. Finite Elements

The 9_4_QuadUP elements used (Figure 4) in this study were developed in plain-strain deformation conditions following the Biot Theory of porous medium. Such elements allow to take into consideration the solid skeleton of the soil (all 9 nodes) and also the fluid phase (4 corner nodes). In particular, the 4 corner nodes have 3 DOFs (2 displacements and 1 pore pressure) while the interior nodes have only 2 DOFs (two displacements).

5. Boundary Conditions

Numerical simulation of dynamic Slope Ground System problems requires many efforts to reproduce the real wave propagation adopting realistic boundary conditions. At this aim, several assumptions such as one elastic half-space under the soil system, massive columns of soil with periodic conditions, one Lysmer-Kuhlemeyer dashpot [31] at the base of the model, and nodal mass simulating water conditions have been introduced. First of all, the entire site is underlain by an elastic half-space that was chosen to be consistent with the existence of bedrock below the slope site, as to allow the energy imparted by the seismic event to be removed from the site itself. The nodes on the base of the model are free to displace in horizontal directions and fixed against the vertical translation. Secondly, in the horizontal direction, the model represents a small section of a presumably infinite (or at least very large) soil domain. To ensure that free-field conditions exist at the horizontal boundaries of the model, the elements in these locations (indicated as red lines in Figure 1) are modelled significantly more massive than the interior mesh. For this purpose, their thickness was increased notably and the nodes on either side of these columns are tied together. Finally, to ensure that the critical portions of the model are not affected by the horizontal boundaries, the free-field columns were located sufficiently far away from the critical regions. Moreover, the dynamic excitation motion was assigned to a Lysmer-Kuhlemeyer [31] dashpot defined through a single zeroLength element. The first end of the dashpot element is fixed against all displacements, while the other end is connected to the soil node in the lower left-hand corner. The constitutive behaviour of the Lysmer-Kuhlemeyer dashpot in the horizontal direction is modelled by a viscous uniaxial material, that requires the dashpot coefficient, C. Following the method of Joyner and Chen [32], this coefficient is defined as the product of the mass density and shear wave velocity of the underlying medium (assumed to have the bedrock properties). In order to ensure that equivalent loading is applied to the model, the dashpot coefficient must be scaled by the area of the base of the model. Finally, the slope was considered completely submerged. Aimed to incorporate the dynamic effects of the water on the site without altering the effective stresses in the soil elements, a nodal mass is assigned manually to each node on the boundary of the mesh which is below water level. For the nodes on the level surface, the horizontal mass is set to zero and the vertical mass is set as the mass of the volume of water supported by the node.

6. Analysis

To control the various parts of the problems and to manage the wise quantities of results, the analysis was split into two consequent steps. The first one is gravity application and the second one is the dynamic analysis itself. Gravity application ensures that the distributions of pore pressure and effective stresses are appropriate for the site conditions prior to the application of a ground motion. Separate recorders were set up as to distinct the gravity analysis from any other postgravity results. Nodal displacements, accelerations, and pore pressures are recorded along with the elemental stresses and strains at each of the nine Gauss points. In order to achieve hydrostatic pore pressure conditions, gravity application analysis is divided into two parts: soil elements are considered to be firstly linear elastic and then elastoplastic. The elastic portion of the gravity analysis is run as a Transient Analysis with very large time steps, thus simulating a static analysis. Gravity is applied for 10 steps with a time step of 500, and 10 steps with a time step of 5000. Therefore, the plastic portion of the gravity analysis is run using smaller time steps to aid in convergence. Dynamic excitation analysis is developed using the method of Joyner and Chen [32]. The force time history was applied as a Plain Load Pattern at the Lysmer-Kuhlemeyer [31] dashpot. Two input motions (shown in Figure 5) are considered at the following hazard levels:(i)T-475: 5% of probability of exceedance in 475 years, representing service limit conditions;(ii)T-5000: 5% of probability of exceedance in 5000 years, both representing collapse limit conditions.

7. Free Field Response

Figures 6 and 7 show the results in terms of horizontal displacement at final time step. As expected, the site is strongly subjected to soil liquefaction and consequent approach fill settlement and lateral spreading. The zoomed views show that the main values are reached in correspondence with the superficial layers that slide on the deeper ones. Under service limit condition (T-475 motion) the displacements are around 0.70 m, while for collapse conditions (T-5000 motion) they grow up to more than 2.50 m. On the superficial layer, three main locations (Figure 8) and respective time histories were considered (Figure 9). Even if the time histories have similar shape, locations 1 and 2 present different values. In particular, for T-475, the displacements at the final time step are around 0.70 m for location 1 and around 0.40 m for location 2. For T-5000, these values grow to around 2.50 m for location 1 and more than 1 m for location 2. Location 3, the upper side of the slope, shows different time histories. This behaviour evidences that the upper sides need more time to slip if compared to the other parts of the slope. In particular, for T-475, the displacement at the final time step is around 0.65 m, while for T-5000 is around 2.80 m.

8. Countermeasures Response

This study aims to systematically explore the effectiveness of densification (Figure 10) consisting of the following:(i) construction of two series of steel 2.5 m diameter hammered piles (EI = 4.032 108 kNm2) in the downhill and uphill of the central zone. The piles are modelled with elastic beam column elements with equivalent flexural characteristics; (ii) densification of the areas that are resulted to be mainly subjected to liquefaction risk. In particular, such areas are evaluated with two different increasing extensions (named model 1 and model 2). Model 1: densification inside the two series of steel piles and for the first 5 meters out of section and 30 meters maximum depth (as painted in green colour in Figure 11). Model 2: the model 1 densification is extended to the superficial layer (at maximum 15 meters depth) out of section 2 for 150 meters (as shown in green colour in Figure 12). The densification technique consists of increasing the superficial layer (named as Fill in Figure 3) relative density to a value equal to 70% that considerably modifies the development of pore pressure. Such value was taken as a reference as the main goal to the technique itself in order to practically annulling liquefaction reproduction;(iii) refill with the same soil material as the original one (named Fill) in order to create an horizontal platform.

The comparison between the two scenarios and the free field response is discussed in correspondence with the two sections pointed out in Figure 13. Longitudinal displacement time histories for the two scenarios are drawn in Figures 14, 15, 16, and 17 for T-475 and T-5000. Scenario 1 is seen to be effective only in correspondence with section 1; while extending the densification in the superficial layer (scenario 2), the displacements in both sections decrease sensitively. These considerations are expressed numerically in Tables 4 and 5, where the two scenario displacements are compared with the free field response values for both sections. In Figures 18 and 19 the displacement time histories at the top of the two piles are compared for both considered motions. Finally, Figures 20 and 21 show the displacement and contour deformation at final time step for the two scenarios for the collapse limit conditions (T-5000).

9. Conclusions

The study conducted in this paper demonstrates OpenSees high potentialities in performing appropriate numerical simulations for predicting liquefaction-induced lateral deformation. The results confirm the assumptions concerning the reproduction of boundary conditions such as the elastic half-space under the soil system, massive columns of soil with periodic boundary conditions, and one Lysmer-Kuhlemeyer dashpot. Moreover, the assignment of nodal mass for each node on the boundary below water level allowed to reproduce the dynamic effects of the water on the site without altering the effective stresses themselves. Finally, splitting of the analysis into two consequent steps (gravity application and dynamic motions) with different computational choices in terms of time steps and force time history reveals its potentiality.

The analyses confirm the vulnerability of such a partially submerged slope in terms of approach fill settlement and lateral spreading due to liquefaction. In particular, the results offer a main reference in the evaluation of the countermeasures design. The extensions for ground improvement densification are compared in terms of displacements in correspondence with significant positions for both service and collapse conditions.

The analyses evidence the effectiveness of densification in reducing the original fill settlement and lateral spreading due to liquefaction. In this regards, comparing these results with economic evaluation can help to quantify the performance and risk of liquefaction using metrics that are of immediate use to both engineers and stakeholders.