1. Introduction

Seismic site-response analysis is routinely used to quantify the influence of surficial soil layers on ground motion propagation. Typically, one-dimensional (1D) frequency-domain (i.e., equivalent linear) or time-domain (i.e., nonlinear) analyses are performed in a total stress framework. However, downhole array and surface acceleration records available at sites that experienced level-ground cyclic liquefaction-for example, 1964 Niigata (Kawagishi-cho apartment buildings; Ishihara and Koga, 1981); 1987 Superstition Hills (Wildlife site [^{1}])-illustrate that porewater pressure (PWP) generation and concurrent strain-softening can significantly alter the acceleration time histories and response spectra [^{2}]. When these sites are analyzed using a total-stress, equivalent-linear constitutive formulation, poor agreement generally is observed between predicted and field measured response [^{3}-^{6}].

In order to improve the predicted seismic site-response, researchers have proposed PWP generation models (e.g., Vucetic [^{4}], Matasovic [^{5}], Seed et al. [^{7}], Martin et al. [^{8}], Dobry et al. [^{9},^{10}], Green et al. [^{11}], Green [^{12}], Polito et al. [^{13}], Ivsic [^{14}]) and effective stress-based soil constitutive models (e.g., Wang et al. [^{15}], Elgamal et al. [^{16}], Park and Byrne [^{17}], Jefferies and Been [^{18}], Boulanger and Ziotopoulou [^{19}]) that can be used for seismic site-response analysis. For example, Booker et al. [20] first proposed an effective stress-based model (named GADFLEA) to study the generation and dissipation of PWP.

In this context, the objective of this study is to evaluate the performance of three coupled PWP and stress-strain constitutive models, used in granular soils, and further recommend them for subsequent application in site-response analysis in effective stresses. For this evaluation, the authors collected 12 bins cyclic simple shear tests (25 individual tests) [^{21}] available in the referred literature. To facilitate comparisons, the laboratory tests were binned based on similar soil response. The PWP coupled stress-strain constitutive models were then evaluated using the laboratory test data, and statistics (i.e., residuals) were used to identify the most effective model.

The three PWP coupled stress-strain constitutive models evaluated corresponded to the:

2. Methodology

In this paper, the authors collected 12 bins cyclic simple shear tests (25 individual tests) [^{21}] available in the referred literature. Those tests were reproduced using three PWP coupled stress-strain constitutive models and the framework of the models used in here is presented next. The results of the numerical simulation were compared with the data obtained in the Laboratory tests. The comparisons with the laboratory tests were evaluated using statistics (i.e., residuals) to identify the most effective model.

2.1. Coupled porewater pressure (PWP) and stress-strain constitutive models

In response to undrained cyclic loading, saturated soils generate excess PWP, which can be separated into transient and residual components. The transient excess PWP results from changes in the mean normal stress during cyclic loading [^{22},^{23}] and thus do not affect the soil effective stress. The residual excess PWP develops from the tendency for volume decrease during undrained cyclic loading, and therefore result in a decrease in effective stress and stiffness. This residual excess PWP can be defined as the PWP in excess of hydrostatic conditions when the cyclic stress is zero [^{9},^{11}]. Most PWP generation models were developed to represent residual excess PWP and they can be incorporated in 1D site-response analysis studies to account for changes in effective stress and stiffness related to PWP generation. A brief description of coupled PWP and stress-strain constitutive models is subsequently presented.

2.1.1. GMP model coupled with hyperbolic stress-strain model

Building on the work by Hardin and Drnevich [^{24}], Matasovic [^{9}] proposed two degradation indices that introduce excess PWP-induced material softening into a simplified hyperbolic constitutive model. The modulus degradation index (δ_{G}) and stress degradation index (δ_{τ}) reduce the shear stress mobilized during the loading-unloading process as a result of PWP increase [^{9}], and are defined as:

where r_{u} = excess PWP/σ'_{vo} or Δu/σ'_{vo}; and ϑ = dimensionless exponent generally equal to 3.5 [^{9}] based on matching the stress-strain hysteresis loops over a wide range of r_{u}-values for Santa Monica Beach sand, Wildlife Site sands A and B, Heber Road point bar (PB) and channel fill (CF) sands. The advantage of the degradation indices is that they can use r_{u} values defined by any PWP generation model.

The modified hyperbolic model presented by Phillips and Hashash [^{25}] was further modified to incorporate the degradation indices previously defined. Moreno-Torres et al. [^{26}] proposed the following equations to compute shear stresses (τ) during loading and unloading-reloading, respectively, corresponding to a given strain.

Loading:

Unloading - Reloading:

where γ_{c} = given shear strain; γ_{rev} = reversal shear strain; τ_{rev} = reversal shear stress; γ_{m} = maximum shear strain; γ_{r} = reference shear strain; t = dimensionless exponent; β^{'} = dimensionless factor; δ_{G} = modulus degradation index; δ_{τ} = stress degradation index; F(γ_{m}) = reduction factor; and G_{0} = initial shear modulus.

Based on the energy-based models proposed by Green et al. [^{11}], Green [^{12}], Polito et al. [^{13}], Davis and Berrill [^{27}], [^{28}] and Berrill and Davis [^{29}], proposed the following empirical expression, termed the GMP model:

where W_{s} = energy dissipated per unit volume of soil divided by the initial mean consolidation stress (i.e., normalized unit energy), and PEC = a calibration parameter termed the pseudoenergy capacity. The GMP model was developed using results from stress-controlled, undrained cyTxC tests performed on nonplastic silt-sand mixtures that ranged from clean sands to pure silts. For each test, W_{s} can be computed as:

where n = number of load increments imposed during cyclic loading on the stress-strain curve until r_{u} = 1 [^{11}]; (_{i} and (_{i+1} = applied shear stress at load increment i and i+1, respectively; and (_{i} and (_{i+1} = shear strain at load increment i and i+1, respectively. Simply stated, Eq. (4) employs the trapezoidal rule to compute the area bounded by the hysteretic stress-strain loop.

For the GMP model, It is used the model parameter correlations reported by Polito et al. [^{13}], which are defined as:

where PEC = pseudoenergy capacity; D_{r} = Relative Density; and FC = fines content.

2.1.2. Elgamal model

The numerical framework of this model uses two-phase (fluid and solid) fully-coupled FE (Finite Element) formulation [^{30},^{31}] based on Biot’s theory [^{32}]. The saturated soil system is analyzed such as two-phase material. The *u-p* formulation (displacement of the soil skeleton, *u*, and pore pressure *p* are the primary unknowns) is a simplified numerical framework of this theory [33]. This model was implemented for the simulation of two-dimensional (2D) and 3D response scenarios [^{30},^{31},^{34},^{35}].

The *u-p* formulation is defined by an equation of motion for the solid-fluid mixture, and by an equation of mass conservation for the fluid phase that incorporates equation of motion for the fluid phase and Darcy's law [^{33}]. The two governing equations are presented in the following finite element matrix form [^{33}]:

where [M] is the total mass matrix, {Ü} is the second derivative of displacement vector, [B]^{T} the strain-displacement transposed matrix, {σ^{'}} the effective stress tensor, Ω is the body volume, dΩ is derivative of the body volume, Q is the discrete gradient operator coupling the solid and fluid phases, [Q]^{T} is the discrete transposed gradient operator coupling the solid and fluid phases, {p} is the pore pressure vector, [S] is the compressibility matrix, {0} is a cero vector, and [H] is the permeability matrix. The vectors {f ^{s}} and {f ^{p}} represent the effects of body forces and prescribed boundary conditions for the solid-fluid mixture and the fluid phase, respectively. Those Equations are integrated in the time domain using a single-step predictor multi-corrector scheme of the Newmark type [^{30},^{33}].

For cyclic and seismic loading scenarios, multi-surface plasticity pressure independent (Von-Mises) and pressure dependent models are able to reproduce the cyclic shear stress-strain response characteristics of the soil conditions. A soil constitutive model pressure-dependent for frictional soils [^{30},^{31},^{34}] was developed based on the Prevost [^{36}] original framework. The focus of the formulation was placed on simulating the liquefaction-induced shear strain accumulation mechanism in clean medium-dense sands [^{31},^{37}]. Special attention was given to the deviatoric-volumetric strain coupling (dilatancy) under cyclic loading, because this causes increased shear stiffness and strength at large cyclic shear strain excursions (i.e., cyclic mobility). The main modeling parameters [^{38}] include standard dynamic soil properties such as low-strain shear modulus and friction angle, as well as parameters to control the dilatancy effects (phase transformation angle, contraction, and dilation), and the level of liquefaction-induced yield strain (γ_{y}).

2.1.3. PM4SAND model

The sand plasticity model presented follows the basic framework of the stress-ratio controlled, critical state compatible, bounding-surface plasticity model for sand presented by Dafalias and Manzari [^{39}]. The Dafalias and Manzari [^{39}] model extended the previous work by Manzari and Dafalias [^{40}] by adding a fabric-dilatancy related tensor quantity to account for the effect of fabric changes during loading. The fabric-dilatancy related tensor was used to macroscopically model the effect that microscopically observed changes in sand fabric during plastic dilation have on the contractive response upon reversal of loading direction.

**
Basic stress and strain terms
**

The basic stress and strain terms for the model are described. The basic model is supported on effective stresses and the conventional prime symbol is skipped from the stress terms for convenience. All stresses are effective for the model considered. The stresses are represented by the tensor σ, the principal effective stresses σ_{1}, σ_{2}, and σ_{3}, the mean effective stress p, the deviatoric stress tensor (s), and the deviatoric stress ratio tensor (r). The implementation was simplified by forming various equations and relationships in terms of the in-plane stresses only. This procedure limits the implementation to plane-strain applications and is not correct for general cases. This has the advantage of simplifying the implementation and improving computational speed by reducing the number of operations. Expanding the implementation to include the general case should not affect the general features of the model. Consequently, the relationships between the various stress terms can be summarized next:

where σ_{xx} is the stress in the plane x, σ_{yy} is the stress in the plane y, σ_{xy} is the shear stress in the plane xy and I is the unitary matrix.

The strain model is represented by a tensor (ε) (strains) that can be separated into the volumetric strain (ε_{v}) and the deviatoric strain tensor (e). The volumetric strain is represented as:

where ε_{xx} is the strain in the plane x, ε_{yy} is the strains in the plane y and the deviatoric strain tensor (e) is represented by:

In incremental form, the deviatoric (de) and volumetric strain (dε_{v}) terms are decomposed into an elastic and a plastic part as they are presented next.

where de^{el} is the elastic deviatoric strain increment tensor, de^{pl} is the plastic deviatoric strain increment tensor, dε_{v}
^{el} is the elastic volumetric strain increment and dε_{v}
^{pl} is the plastic volumetric strain increment.

3. Results

**3.1. Prediction of PWP using coupled PWP and stress-strain constitutive models**

The authors used 12 cySS bins for the prediction exercise (Table 1). Model parameters for each of the models were determined using recommendation and formulation for that purpose. Figs. 1 - 3 compare the PWP time histories computed using the GMP model coupled with hyperbolic stress-strain, Elgamal and PM4SAND models to the measured responses of cySS tests in bins S3B, S6B, and S8B, respectively. In general, the GMP model coupled with hyperbolic stress-strain reasonably predict PWP generation in loose to medium-dense specimens (Fig. 1 (a) and (b), respectively). As expected, the model more poorly predict PWP generation in dense specimens (Fig. 1 (c) ). The Elgamal and PM4SAND models fairly predict PWP generation in loose to dense specimens (Figs. 2 through 3, respectively) for different CSR magnitudes (i.e., different FSliq).

Fig. 4 presents the residuals for measured and predicted r_{u} values (at each cycle) for the 12 cySS bins (25 individual cySS tests) in the prediction dataset. These data illustrate that the coupled GMP with stress-strain constitutive model yields small residuals (±0.25), but contain some bias. The Elgamal and PM4SAND models yield medium residuals (±0.5) with high bias concentrated in the underprediction side for Elgamal model and PM4SAND presents bias in the under and overprediction sides. However, all models qualitatively capture initial liquefaction.

**3.2. Evaluation of predicted stress-strain behavior using coupled PWP and stress-strain constitutive models**

The authors used the same 12 cySS bins considered for the PWP prediction exercise to evaluate the stress-strain behavior of the models. Model residuals were calculated using computed and measured shear stresses when the measured shear strain, γ_{c} = 0%.

Figs. 5 through 7 compare the stress-strain and excess r_{u}-strain behavior computed using the GMP model coupled with hyperbolic stress-strain, Elgamal and PM4SAND models to the measured responses of cySS tests in bins S3B, S6B, and S8B, respectively. In general, the three models reasonably predict stress-strain and r_{u}-strain behavior in loose to dense specimens for strains less that γ_{c} = 5%. As expected, the coupled GMP model with stress-strain constitutive model more poorly predict stress-strain and r_{u}-strain behavior when dilation becomes more pronounced (γ_{c} > 5% and r_{u} > 0.65) and significant modulus degradation occurs. Instead, Elgamal and PM4SAND can capture the initial part of dilation spikes and the degradation of the modulus in a good manner when γ_{c} > 5% and r_{u} > 0.65, but the both codes are not able to support more than a couple of cycles after r_{u} > 0.98 and the codes stop the calculation.

Figs. 5 through 7 compare the stress-strain and excess r_{u}-strain behavior computed using the GMP model coupled with hyperbolic stress-strain, Elgamal and PM4SAND models to the measured responses of cySS tests in bins S3B, S6B, and S8B, respectively. In general, the three models reasonably predict stress-strain and r_{u}-strain behavior in loose to dense specimens for strains less that γ_{c} = 5%. As expected, the coupled GMP model with stress-strain constitutive model more poorly predict stress-strain and r_{u}-strain behavior when dilation becomes more pronounced (γ_{c} > 5% and r_{u} > 0.65) and significant modulus degradation occurs. Instead, Elgamal and PM4SAND can capture the initial part of dilation spikes and the degradation of the modulus in a good manner when γ_{c} > 5% and r_{u} > 0.65, but the both codes are not able to support more than a couple of cycles after r_{u} > 0.98 and the codes stop the calculation.

Fig. 8 presents the residuals determined for measured and predicted τ/σ'_{vo} values (CSRs) at each cycle for the 12 cySS bins (25 individual cySS tests) in the prediction dataset. These data illustrate that the coupled GMP with stress-strain constitutive model and PM4SAND yield relatively low residuals (< ±0.25) for r_{u} < 0.65 although the models generally overestimate the values of r_{u}. For high r_{u}-values (r_{u} > 0.65), the residual increases (±0.35), although the bias decreases for both models. Instead, Elgamal model yields very low residuals (< ±0.05) through all r_{u}-values, and no bias is presented with means that this model can reproduce the strain-stress behavior in a good shape.

4. Conclusion

In this paper, three coupled PWP and stress-strain constitutive models (termed GMP, Elgamal and PM4SAND) are described. Those constitutive model uses different approach to incorporate porewater pressures (PWP) predicted by PWP generation models available in the literature. The evaluation process consisted of a prediction phase, the coupled GMP and stress-strain constitutive model exhibited relatively low residuals compared with Elgamal and PM4SAND model (more advanced constitutive models), but the bias in the coupled GMP and stress-strain constitutive model is more reduced than the other models considered here.

Given the ability of the coupled GMP and stress-strain constitutive model to reasonably predict PWP generation in the prediction phase over the widest range of r_{u} and over all values of D_{r}, it was expected a reasonable stress-strain prediction. The modified hyperbolic constitutive model coupled with the GMP model qualitatively captures the stress-strain behavior for r_{u} ≤ 0.65, but the inability of the hyperbolic constitutive model to reproduce dilation at r_{u} values exceeding approximately 0.65 reduces the accuracy of the coupled GMP and stress-strain model compared with the Elgamal and PM4SAND models can reproduce the dilation process in better way. The advanced models predict PWP fairly reasonably until initial liquefaction with bias more pronounced. Even, the found limitations that were exposed, all the models can be used in site response analysis with good results in reproduction of the soil response in an earthquake event when liquefaction is an issue.