LEAP-UCD-2017 Centrifuge Test Simulation at UNINA

Within the framework of the LEAP-UCD-2017 exercise, Type B simulations of centrifuge tests were conducted assuming a hypoplastic constitutive model for sand. Differently from the most common elastoplastic approach, the hypoplasticity does not decompose the strain rate into elastic and plastic parts and does not use explicitly the notions of the yield surface and plastic potential surface. The process followed to calibrate the constitutive model is presented in detail. The initial state of stresses in the analyzed mesh, the key parameters used in the dynamic simulation phase, and a comparison of the simulation with some experimental results are reported. All the simulations were performed using the model parameters calibrated by using the laboratory test data. Finally, a sensitivity analysis of computed displacement to soil density and ground motion intensity show the in ﬂ uence of such factors on the seismic soil response of lique ﬁ able soils.


Introduction
The aim of this paper is to describe the process followed to calibrate the constitutive model, to obtain the initial state of stresses in the analyzed soil-deposit mesh, to properly document the key parameters used in the seismic simulation phase, and to briefly present a comparison of the simulations with some experimental results of the LEAP-UCD-2017 exercise.
The content of the paper focuses on the following features of the performed numerical simulations: (1) the constitutive model used in the simulations, which is quite different from the most popular approaches that are currently adopted for studying liquefaction; (2) the calibration process, which is described in details; and (3) a sensitivity study. The set of model parameters were calibrated on both triaxial and simple shear tests. The same parameters are used in the simulation of centrifuge tests, which are carried out at the prototype scale. Type B simulations section includes a description of the discretization, the numerical integration scheme, and numerical damping used in the analysis.
A sensitivity study (simulations NS1-NS7, exploring the sensitivity of displacement to soil density and ground motion intensity) has been carried out with the same parameters and using the same numerical model adopted in the Type B simulations.
The paper is divided into four sections. After this Introduction, Sect. 22.2 covers the essential features of the constitutive model, the final model parameters, the calibration philosophy, and the assumptions used in the calibration process. It also presents a comparison between the predicted and experimental cyclic strength curves. Section 22.3 discusses the main features of the numerical analysis platform used in the simulation, the model geometry and the discretization details, the boundary conditions of the numerical model, the solution algorithm employed, and some assumptions used in the reported analyses. Section 22.4 presents a comparison between the results of the Type B predictions and those of one of the centrifuge experiments. Finally, Sect. 22.5 shows the main results of the sensitivity study to explore the role of the soil density and ground motion intensity on slope behavior.

Hypoplastic Constitutive Model
The constitutive model used in the simulation exercise is the hypoplastic model with intergranular strain concept (Von Wolffersdorff 1996;Niemunis and Herle 1997), implemented in the finite element code PLAXIS 2016 (Mašín 2010).
Hypoplasticity is a particular class of incrementally nonlinear constitutive models, developed specifically to predict the behavior of soils. The basic structure of the hypoplastic models has been developed during the 1990s at the University of Karlsruhe. In hypoplasticity, unlike in elastoplasticity, the strain rate is not decomposed into elastic and plastic parts, and the models do not use explicitly the notions of the yield surface and plastic potential surface. Still, the models are capable of predicting the important features of the soil behavior, such as the critical state, dependency of the peak strength on soil density, nonlinear behavior in the small and large strain range, and dependency of the soil stiffness on the loading direction.
Hypoplastic models can well represent deformations due to rearrangements of the grain skeleton. However, applications of the hypoplasticity to cyclic loading with small amplitudes revealed some shortcomings. In the basic hypoplastic model (Von Wolffersdorff 1996), the stiffness during loading and reloading is almost identical, which leads to an unrealistic large accumulation of deformation in cyclic loading. Therefore, for simulations of cyclic or seismic loading, it is necessary to enhance the model with an extension, considering an increased stiffness in a small strain range. Thus, the applied hypoplastic constitutive equation corresponds to the model by von Wolffersdorff (1996) with the intergranular strain extension according to Niemunis and Herle (1997). In Plaxis, the constitutive model is integrated using explicit adaptive integration scheme with local sub-stepping (Mašín 2010).
The parameters of the hypoplastic model (Von Wolffersdorff 1996) are critical friction angle ϕ c , critical void ratio at zero pressure e c0 , void ratio at the maximum density at zero pressure e d0 , maximum void ratio at zero pressure e i0 , granular hardness h s , stiffness exponent n, and the exponents α and β. The parameter α affects the peak friction angle in dependence of relative density, while the soil stiffness is affected by the parameter β. Moreover, the intergranular strain extension (Niemunis and Herle 1997) requires further five parameters: R, m R , m T , β r , and χ. The stiffness increase after a change in the loading direction of 180 and 90 is controlled by the parameters m R and m T , respectively, with the loading direction referring to the direction of the strain rate. If m R ¼ 0, the intergranular strain concept is switched off, and the problem is simulated using the basic hypoplastic model (Von Wolffersdorff 1996).
The reduction in the shear modulus with increasing shear strain is controlled by the parameters R, β r , and χ.
Finally, in addition to the initial void ratio, e 0 , further six state parameters may be specified, δ 11 , δ 22 , δ 33 , δ 12 , δ 13 , δ 23 , which are the initial values of the intergranular strain tensor δ ¼ (δ 11 , δ 22 , δ 33 , 2δ 12 , 2δ 13 , 2δ 23 ). They describe the accumulation of permanent deformation. Their initial values have been set to zero. A more detailed description of the model can be found in the referenced publications. Table 22.1 reports the input parameters used in the hypoplastic model, as listed in Plaxis 2016 (Brinkgreve et al. 2016). The model requires overall 21 input parameters. Overall slope of compression curve n Curvature of compression curve e d0 Densest particle packing at zero mean stress e c0 Critical particle packing at zero mean stress e i0 Loosest particle packing at zero mean stress α Peak friction angle as function of relative density β Shear stiffness as function of relative density m R Initial shear stiffness for initial and reverse condition m T Stiffness upon neutral condition R Size elastic range β r Rate stiffness degradation χ Rate stiffness degradation SV: e 0 Initial void ratio δ ij i, j ¼ 1,2,3 Initial value of the intergranular strain tensor

Model Calibration
The approach used in the calibration of the model parameters is hereafter explained. The calibration made use of the results of the provided cyclic triaxial stresscontrolled tests and additional element test data (i.e., monotonic triaxial tests, simple shear tests). Shear strength parameters are computed from the monotonic triaxial test data, available on the NEES Hub (https://nees.org/dataviewer/view/1064:ds/1189). Drained triaxial compression tests, carried out by Vasko (2015) on loose and dense specimens, were used to define the critical state line in the plane q:p 0 and the constant volume friction angle, ϕ 0 c . As well known, the evaluation of critical state conditions in triaxial tests is a very complex issue, being such a test intrinsically affected by a number of experimental limitations (localization, bulging, shear stresses on the rough porous stones, difference between local and external displacements, etc.). One of the best ways to evaluate the final state is therefore by looking to dilatancy trend at the end of the tests. Based on all the elaborations of the available experimental data and considering that the hypoplastic model assumes the critical state as an asymptotic state at infinite strains, in this case the best fit of this parameter is the following: The critical void ratio at zero mean stress, e c0 , has been assumed, by extrapolation of experimental values as e c0 ¼ 0:746 ð22:2Þ Hence, the other two parameters representing the initial values at p 0 ¼ 0 of maximum and minimum void ratios have been set equal to e d0 ¼ e min ð22:3Þ e i0¼1:2 e c0 ð22:4Þ The parameter α was computed as a function of the void ratios, e, e c , e d , related to the peak mean stress, according to the equation proposed by Herle and Gudehus (1999): where the peak ratios are function of the peak friction angle, ϕ p , while the parameter a is function of the critical friction angle, ϕ c , Finally, the parameter β was calibrated on the results of the available drained triaxial tests on dense specimens (Vasko 2015).
The parameters h s and n were calibrated by interpolating the G 0 :p 0 relationship proposed by Robertson et al. (1995) for dense Ottawa C-109 sand ( Fig. 22.1) with the equation where the constant of proportionality m R has been set to 5 and ð22:11Þ .1 Shear wave velocity measurements, V S , vs. mean effective stress, p 0 during isotropic consolidation for Ottawa C-109 sand (modified after Robertson et al. 1995) In Eq. (22.10), a shear wave velocity of 250 m/s was assumed for a reference value of mean stress, p 0 ¼ 100 kPa (Fig. 22.1). Then, the stiffness at small strains, G 0 was computed as After that, the parameters h s and n were calibrated to fit the G 0 :p 0 relationship. According to Mašín (2015), the parameter m T was considered equal to a fraction of m R : The implemented model requires also the definition of the parameter p t , which is a shift of the mean stress due to cohesion. For the basic hypoplastic model, p t is set equal to zero, but a nonzero value of p t is needed to overcome problems with stressfree state (Mašín 2010). In order to identify the parameters R, β r , and χ, simulations of strain-controlled cyclic shear tests were carried out, and the secant shear modulus G was evaluated from them. The numerical results have been compared with the experimental data for (G/G 0 , γ) reported in literature for Ottawa sand 20-30 (Alarcon-Guzman et al. 1989). The values of R and χ were assumed as R ¼ 0.0001 and χ ¼ 1.0. Hence β r ¼ 0.4 was determined (Fig. 22.4).
The cyclic triaxial tests were simulated by using the "soil test" tool of Plaxis 2016 (Brinkgreve et al. 2016). The prescribed cyclic stress ratio, CSR CTX , was imposed for a number of cycles, N, large enough to induce liquefaction. The model parameters that were finally obtained are summarized in Table 22.2. They were adopted for the simulation of the centrifuge tests.

Liquefaction Strength Curves
For each group of tests with the same void ratio, the relative density was defined as the mean value of the relative densities calculated by the different experimenters (see Table 22.3). For the sake of simplicity, the obtained mean values were rounded down. Table 22.3 shows that the available laboratory tests were performed on specimens of dense sand. Figure 22.5 reports an example of the stress-strain and stress path curves simulated by the model for a cyclic stress-controlled test (e 0 ¼ 0.515, p 0 0 ¼ 100 kPa, and  (2017) and compares them with the experimental data. The simulated stress-strain cycles are quite different from the experimental results, while the prediction of the stress path is satisfied especially during the extension phase of the test. Liquefaction has been assumed when a computed axial strain of 2.5% was achieved. It is worth noting that such a value is affected by ratcheting that is a common issue for many constitutive models among those implemented in commercial codes. The liquefaction strength curves, obtained from the simulated cyclic stresscontrolled triaxial tests, are plotted and compared with the experimental results ( Fig. 22.6).

Analysis Platform
The simulations of the centrifuge tests have been carried out by using Plaxis 2D (Brinkgreve et al. 2016). This is a 2D commercial finite element method (FEM) code that includes several constitutive models. Many of them are available as user-defined model, such as the hypoplastic model with intergranular strain concept, implemented and made available in the platform by Mašín (2010). The main reason to use Plaxis rather than other platforms is that, although not specifically oriented to solve boundary value problems in earthquake geotechnical engineering, this numerical code is quite well widespread in the community of geotechnical practitioners. Hence, it was interesting to check the possible benefit of a rigorous validation of numerical simulation procedures implemented in Plaxis through experimental data, in order to apply those procedures to a boundary value problem involving soil liquefaction.
Unfortunately, modelling the pore pressure dissipation during a dynamic analysis was not possible in Plaxis at the time when simulations were carried out. Therefore, pore pressure buildup needed to be calculated by imposing undrained conditions during shaking. This is obviously a limitation of the selected calculation method, since it affects the prediction of pore pressure buildup.

Solution Algorithm and Assumptions
The Newmark time integration scheme is used in the simulations where the time step is constant and equal to the critical time step during the whole analysis. The proper critical time step for dynamic analyses is estimated in order to model accurately wave propagation and reduce error due to integration of time history functions. First, the material properties and the element size are taken into account to estimate the time step, and then the time step is adjusted based on the time history functions used in the calculation. During each calculation step, the Plaxis calculation kernel performs a series of iterations to reduce the out-of-balance errors in the solution. To terminate this iterative procedure when the errors are acceptable, it is necessary to establish the out-of-equilibrium errors at any stage during the iterative process automatically. Two separate error indicators are used for this purpose, based on the measure of either the global equilibrium error or the local error. The "global error" is related to the sum of the magnitudes of the out-of-balance nodal forces. The term "out-of-balance nodal forces" refers to the difference between the external loads and the forces that are in equilibrium with the current stresses. Such a difference is made nondimensional dividing it by the sum of the magnitudes of loads over all nodes of all elements. The "local error" is related to a norm of the difference between the equilibrium stress tensor and the constitutive stress tensor. It is made nondimensional dividing by the maximum value of the shear stress as defined by the failure criterion. The values of both indicators must be below a tolerated error set to 0.01 for the iterative procedure to terminate. In general, the solution procedure restricts the number of iterations that take place to 60, in order to ensure that computer time does not become too high.
A full Rayleigh damping formulation has been considered in the simulation, and the coefficients α RAY and β RAY are equal to 0.025 and 0.63 Â 10 À3 , respectively.

State of Stresses and Internal Variables of the Constitutive Model in the Pre-shaking Stage
All the centrifuge models are subjected to centrifugal accelerations that are increased from 1 g to a designated value during the spin-up process. The soil specimen is subjected to an increasing centrifugal acceleration that affects the state of stresses and the subsequent seismic response. The initialization of stresses and internal variables of the constitutive model after the centrifuge spin-up and just before the start of shaking was addressed in the simulation through a static analysis aimed at computing the initial effective stress state in the soil profile.

Results of the Dynamic Analysis
Relevant results of the simulations are hereafter discussed, while all the data of the numerical analysis can be found in the platform DesignSafe (www.designsafe-ci. org). As described in Kutter et al. (2019), the centrifuge model is monitored with accelerometers and pore pressure transducers as reported in Fig. 22.8. Figure 22.9 shows the comparisons among simulated and predicted time histories of acceleration, pore pressure, and horizontal displacement at the surface for the UCD1 centrifuge test. A reasonable prediction of acceleration and pore pressure is provided by the numerical simulation, while simulated displacements strongly underestimate the experimental time history.
Since the horizontal displacements result to be the most critical aspect to be reproduced, the simulated displacements at the surface are compared to the experimental ones in all the other centrifuge tests (Fig. 22.10). The comparison shows that the amplitude of the simulated time histories is quite similar to the experimental one in centrifuge test KyU3 and Ehime 2 ( Fig. 22.10c, d).

Sensitivity Analysis
The considered key input parameters of the sensitivity study are relative density and ground motion intensity, as described below. Table 22.4 shows the characteristics of each analysis.   The sensitivity study incorporates Type C simulation of the RPI-1 test (simulation NS-1), Type B simulation of RPI-2 (NS-6) test, and an almost Type A simulation of RPI-3 (NS-2). NS-1 to NS-3 has been used to deduce the sensitivity of simulations to relative density. NS-1, NS-4, and NS-5 allowed deduction of the sensitivity of simulations to motion intensity. NS-6 and NS-7 allowed assessment of the influence of superimposed high frequencies on simulation results. Figures 22.11,22.12,and 22.13 report the comparison between computed and measured acceleration and pore pressure time histories related to the RPI01 centrifuge test. Experimental acceleration time histories are well-reproduced until liquefaction condition is attained; after that, the simulations lead to an underestimation of the observed trend. Underestimation is also related to the prediction of the pore pressure in the center of the model, while pore pressure closer to the boundaries of the laminar box are better simulated. Table 22.5 summarizes the results of the analysis in terms of computed horizontal displacement at the middle point of the model surface. It must be noted that shaking was done imposing undrained conditions and a consolidation phase was run only after the end of the shaking; hence the analyses did not provide information about the duration of liquefaction after the end of shaking. In the following sections, the effects of relative density and motion intensity are shown.

Influence of Relative Density
With reference to the influence of relative density, three different analyses were carried out with a relative density of 65% (NS-1), 50% (NS-2), and 75% (NS-3), as reported in Fig. 22.14.

Influence of Motion Intensity
With reference to the influence of input motion, three different analyses were carried out with a maximum acceleration of 0.15 g (NS-1), 0.25 g (NS-4), and 0.11 g (NS-5), as reported in Fig. 22.16. As expected, significant differences can be observed in the simulated acceleration time histories along the vertical profile ( Fig. 22.16), while higher pore pressure are generated for increasing values of the maximum acceleration ( Fig. 22.17).

Conclusion
This paper described the model calibration and numerical simulations of centrifuge test for LEAP-UCD-2017, carried out at the University of Napoli Federico II (Italy). The adopted hypoplastic model allowed a relatively easy calibration of the model parameters compared to other advanced models. Cyclic resistance curves are well- predicted although with some drawbacks on the shear strain histories, where "ratcheting" effects can be identified. Modelling pore pressure dissipation during dynamic analysis was not possible in the version of Plaxis available at the time of simulations. Therefore, pore pressure buildup needed to be calculated by imposing undrained conditions during shaking. Nevertheless, in most cases, the buildup of excess pore pressure was reasonably captured. However, the magnitude of calculated horizontal displacements was in most cases unsatisfactory.
The sensitivity study highlighted the crucial role of relative density and the intensity of the applied input motion on the liquefaction behavior of the slope.