Lianxu Zhou, Ph.D. ^{1} ; Michele Barbato, F.ASCE ^{2} ; and Aijun Ye ^{3}

Abstract: Scour can increase the earthquake-induced damage in pile group foundations. Quantifying the parameter sensitivity of the seismic performance for scoured pile group foundations is essential for the optimal design and retrofit of bridges located in seismic-prone regions. Such quantification requires numerical models that are computationally efficient and accurate in describing the mechanical behavior associated with the complex soil-foundation-structure interaction of these systems. This study proposes an efficient finite-element model (FEM) of pile groups based on a beam on the nonlinear Winkler foundation approach. This FEM uses asymmetric p-multipliers to describe the different soil resistance exerted on leading and trailing piles when applying cyclic lateral loads. The proposed FEM is validated by comparing the numerical response with the experimental measurements taken from a quasi-static test available in the literature and is used to perform an extensive parametric sensitivity analysis to quantify the response sensitivity to 11 structural and soil parameters. Tornado diagrams are employed to identify an importance ranking of these parameters on the seismic performance of scoured pile groups. The obtained results indicate that the proposed FEM is able to capture both the global and local structural responses of pile group foundations. The parametric sensitivity analysis shows that pile group foundations have considerable ductility capacity. Pile diameter and axial load ratio of piles are the most important parameters for the seismic performance of pile groups. Increasing the pile diameter is the most efficient approach to improve the seismic performance of a pile group when considering scour effects. The seismic performance of a scoured pile group deteriorates with increasing piles' axial load ratio. For a deep pile group foundation, seismic performance is very little sensitive to pile length and relative density of sand. Based on the results of the parametric analysis, recommendations are proposed for the seismic design of pile group foundations with scour effects. DOI: 10.1061/JBENF2.BEENG-5861. © 2023 American Society of Civil Engineers.

Author keywords: Bridge scour; Pile group effect; Finite-element model; Parameter sensitivity analysis; Tornado diagram; Ductility capacity; Soil-pile interaction.

Introduction

Reinforced concrete (RC) pile group foundations are widely used in bridge engineering practice (Fayyazi et al. 2012)²⁴. In addition to bearing vertical loads, pile group foundations can be affected by significant lateral loads produced, e.g., by earthquake ground motions. Existing seismic design specifications often require that the pile group foundations remain in their elastic behavior range under design-level earthquake excitations based on the capacity design philosophy (Mander et al. 1998)³⁷ or allow for the formation of plastic hinges at the pile-cap connection (AASHTO 2022)². However, pile damage is often unavoidable when strong earthquakes occur (Kawashima et al. 2009)³¹; (Wei et al. 2008)⁶⁰. In addition, scour has been reported as the main hazard causing bridge failures (Alipour and Shafei 2016)⁶; (Capers et al. 2013)¹⁹; (Wardhana and Hadipriono 2003)⁵⁹ by exposing pile group foundations and reducing their lateral load capacity. Recent studies found that scour makes pile group foundations more prone to earthquake-induced damage than columns in a pile-supported bridge system (Wang et al. 2014)⁵⁸, (Wang et al. 2015)⁵⁵, (Wang et al. 2019)⁵⁷. Therefore, the seismic design of scoured pile group foundations is a research topic of significant practical relevance. Extensive experiment- and/or simulation-based studies have focused on soil-pile interaction effects, including pile group effects and inertial/kinematic interaction of soil-pile-structure systems, in which the piles remained in the elastic range or exhibited a limited amount of plastic behavior (Boulanger et al. 1999)¹⁴; (Brown et al. 1988)¹⁷; (Hussien et al. 2016)²⁷; (Rollins et al. 2005)⁴⁸. To investigate the seismic failure mechanism and ductility capacity of RC pile group foundations, a few experimental studies adopted cyclic static loads imposed on pile group foundations to simulate earthquake loadings (Banerjee et al. 1987)¹¹; (Chai and Hutchinson 2002)²⁰; (Liu et al. 2020)³⁶; (Park and Falconer 1983)⁴⁵; (Wang et al. 2016)⁵⁶; (Zhou et al. 2021a)⁶².

Significant research efforts have been devoted to developing numerical approaches for simulating and predicting the behavior of pile group foundations subject to lateral loads. The lateral response of piles is commonly analyzed using the beam on a nonlinear Winkler foundation (BNWF) approach (Adeel et al. 2020)⁴; (Allotey and El Naggar 2008)⁷; (Heidari et al. 2014)²⁶; (Liu et al. 2020)³⁶; (Matlock and Ripperger 1956)³⁹; (Wang et al. 1998)⁵⁴; (Zhang and Hutchinson 2012)⁶¹, in which the soil-pile interaction is described using a p-y curve, where p denotes the soil resistance and y denotes the lateral displacement of the pile. Boulanger et al. (1999)¹⁴ developed a

nonlinear constitutive model based on a combination of three components in series: (1) an elastic spring in parallel with a dashpot to model radiation damping, (2) a plastic spring, and (3) a gap component consisting of a nonlinear closure spring in parallel with a nonlinear drag spring. This material constitutive model was implemented as the uniaxial material denoted PySimpleI in the Open System for Earthquake Engineering Simulation (OpenSees) platform, version 3.3.2 (McKenna 2011)⁴¹ and has been widely adopted by both the practicing and academic civil engineering communities (Brandenberg et al. 2007)¹⁶; (Hutchinson et al. 2004)²⁹; (Kramer et al. 2008)³⁴. This BNWF-based soil-pile interaction modeling approach has also been validated by a series of centrifuge tests (Boulanger et al. 2003)¹⁵, quasi-static tests (Hutchinson et al. 2005)²⁸; (Zhou et al. 2022b)⁶⁵, and shake table tests (Shang et al. 2018)⁵⁰. An aspect that has received significant attention for pile group foundations is the quantification and modeling of the so-called pile group effect (Brown et al. 1988)¹⁷, which corresponds to the reduction of the lateral capacity of a pile group with respect to the sum of the lateral capacities of the individual piles. This phenomenon is produced by the overlapping of the soil zones affected by the different piles, which is particularly evident for closely spaced pile groups. Due to the pile group effect, different pile rows provide different contributions to the overall lateral capacity of the pile group, with generally higher loads applied to the leading piles, which also tend to exhibit higher curvatures and higher ductility demands (Rollins et al. 2005)⁴⁸; (Wang et al. 2019)⁵⁷; (Zhou et al. 2021a)⁶². Therefore, developing an accurate modeling method for the pile group effect is a critical issue in predicting the performance and the ductile behavior of pile group foundations subject to seismic and lateral cyclic loads. Brown et al. (1988)¹⁷ proposed the p-multiplier method to simulate the pile group effect. In this method, the lateral soil resistance of each pile row at a given embedded depth is described using a p-y spring, in which the load p is reduced by a load reduction factor, f_{m} (i.e., the p-multiplier). Different p-multipliers can be used for different pile rows to reproduce the experimentally measured effect of different contributions between leading and trailing piles, providing a simple and widely adopted approach in engineering practice (AASHTO 2020)¹. However, for seismically excited pile groups, the loading direction continually changes during the seismic excitation, with piles frequently interchanging their condition between leading and trailing piles multiple times during any given seismic event. To address this issue, previous studies adopted an approximate approach using the average constant value of the p-multiplier for all piles in the group, called group efficiency factor or group reduction factor, to modify the soil resistance in front of the piles (Adeel et al. 2020)⁴; (Brown et al. 2001)¹⁸; (Lemnitzer et al. 2010)³⁵; (Liu et al. 2020)³⁶. This modeling approach has been shown to provide accurate estimates of the global force-displacement response of a pile group (Lemnitzer et al. 2010)³⁵; (Liu et al. 2020)³⁶. However, this uniform reduction factor cannot simulate the difference in lateral soil resistance among different rows in a pile group subjected to cyclic loads, which would be necessary to capture the curvature differences between piles in different rows.

In addition, experimental testing can only investigate a limited set of physical and modeling parameter combinations due to the high cost associated with each experimental sample. However, many different parameters can affect the behavior of pile group foundations, particularly when subjected to scour hazards (Blanco et al. 2019)¹³; (Song et al. 2022)⁵². To mitigate this issue, a numerical simulation based on finite-element modeling (FEM) can be used to investigate the effects of parameters for which direct experimental testing is unfeasible. Blanco et al. (2019)¹³ carried out a numerical parametric pushover analysis on the ductile behavior of RC pile group foundations. However, their study was based on a limited number of parameters and, in particular, did not investigate (1) the effects of cyclic loading, (2) different pile group configurations beyond a 2 \times 3 configuration, and (3) the relative importance of different parameters on the seismic performance of the pile group foundations. Therefore, a reliable finite-element modeling approach is a necessary complement to experimental investigations to understand the effects and relative importance of the different parameters that are expected to affect the performance of scoured pile group foundations subject to seismic actions. In addition, an extended parametric sensitivity analysis would significantly help engineers identify the critical parameters for improving the seismic performance of pile foundations with scour potential.

This paper proposes a practical and straightforward FEM approach, based on a BNWF model with asymmetric p-multipliers, to simulate the soil-pile interaction of pile groups with multiple rows of piles subjected to cyclic loading from seismic excitation. The proposed model is validated by comparing the experimentally measured and numerically simulated global and local response of a scoured RC pile group foundation in sandy soil. This study also performs a detailed parametric sensitivity analysis based on the newly proposed numerical model to identify critical structural and/or soil parameters affecting the seismic performance of scoured pile group foundations and to determine the sensitivity rankings of these parameters.

Novelty and Relevance

This paper proposed, for the first time, an asymmetric p-multiplier to better model the differences between leading and trailing piles in RC pile group foundations with multiple pile rows in the direction of the loading and subject to cyclic loading. The proposed approach could be easily extended to other types of piles. This study also performed, for the first time, a comprehensive parametric sensitivity analysis for scoured RC pile group foundations subject to cyclic loading. The sensitivity analysis results were reported in terms of the effects on the piles' performance, particularly in terms of damage levels exhibited by the piles after cyclic loading. The results presented in this paper could represent the basis for future improvements in the design and seismic retrofit of pile group foundations under the combined effects of earthquake loads and scour.

Numerical Modeling and Validation

Description of Experimental Test

Zhou et al. (2021a)⁶² carried out a series of quasi-static tests on RC specimens of scoured 2 \times 3 pile group foundations to investigate their ductile behavior during cyclic loading and their postearthquake vertical load-carrying capacity under different lateral damage states. Their experimental data for Specimen #3 (loaded to a maximum lateral displacement level of 100 \mathrm{~mm}) were used to validate the numerical model in this study. Fig. 1 shows the test layout. The RC pile group consisted of six circular piles with a diameter D = 0.12 \mathrm{~m} and a length H = 4.30 \mathrm{~m}. The pile head was connected together by a concrete cap with a dimension of 1.50 \times 1.00 \times 0.60 \mathrm{~m}. The center-to-center spacing of adjacent piles was 0.36 \mathrm{~m} (i.e., 3D). The specimen was positioned in the center area of a container with an inside dimension of 3.10 (\mathrm{length}) \times 1.50 (\mathrm{width}) \times 4.20 \mathrm{~m} (height) and embedded in homogeneous sand with an average relative density D_{r} = 55\%. The embedded depth was 3.70 \mathrm{~m} (30.83D), and the exposure length of each pile was equal to 0.60 \mathrm{~m}, representing a scour depth of 5D. Table 1 lists the property

Fig. 1. Quasi-static test overview: (a) schematic side view layout; (b) pile steel reinforcement; (c) full view of test; and (d) view of sand and aboveground piles (all units are in cm if not otherwise indicated). (Data from Zhou et al. 2021a.)

Table 1. Sand parameters in the experimental test from Zhou et al. (2021a)⁶²

of the sand used in the test. The total initial axial force applied on the piles was equal to 85.4~\mathrm{kN}, corresponding approximately to an axial load ratio, \eta, equal to 5\% for each pile. The axial load ratio is defined here as follows:

where P = axial (dead) load exerted on the individual pile; f_{c} = peak strength of the unconfined concrete (with f_{c} = 25.20\mathrm{MPa} for this experimental test); and A_{g} = pile gross cross-sectional area (with A_{g} = 0.0113\mathrm{m}^{2} for this experimental test). The lateral load was provided by an actuator, identified as Actuator #1 in Fig. 1(a). Fig. 2 presents the lateral loading protocol for Specimen #3.

Fig. 1(b) shows the steel reinforcement configuration for each pile. The longitudinal steel reinforcement ratio was 1.5\% and was provided by six 6-mm-diameter longitudinal steel rebars. The core concrete of the piles was spirally confined by 3.5-mm-diameter galvanized-iron-wires (GIWs) with a center-to-center spacing of 35~\mathrm{mm}, leading to a transverse reinforcement ratio of 1.215\%. Tables 2 and 3 summarize the mechanical parameters of the concrete and steel used in the pile group specimen.

Fig. 2. Loading protocol for Specimen #3 according to Zhou et al. (2021a)⁶².

Selection of p-Multipliers for Pile Group Effect

The p-multipliers have been typically obtained from full- or small-scale experimental quasi-static or centrifuge tests, or estimated using finite-element analysis, often based on three-dimensional (3D) models. Table 4 summarizes the values of the p-multipliers for pile groups in sandy soil reported in the literature. It is observed that the value of the p-multiplier mainly depends on pile layout, pile-row location, and the ratio of pile spacing to diameter, as reported in previous studies (AASHTO 2020)¹; (Adeel et al. 2020)⁴, whereas the relative density of sand seems to have a relatively small effect. The value of f_{m} increases with the increase of pile center-to-center spacing S, and the p-multiplier for leading piles is generally larger than that for trailing piles. Based on the collected data of f_{m} listed in Table 4, the mean values of f_{m} for first-, second-, and third-row piles of a three-row pile group with S = 3D are 0.75, 0.41, and 0.33 and their standard deviations are 0.058, 0.034, and 0.055, respectively. For a two-row pile group with S = 3D, the mean values of the p-multiplier for the first- and second-row piles are 0.84 and 0.56, and their standard deviations are 0.048 and 0.082, respectively.

Fig. 3 compares the p-multipliers for three-row pile groups obtained from the literature and reported in Table 4 with the values recommended by the AASHTO specifications as a function of the ratio S/D (AASHTO 2020)¹. The values suggested by AASHTO refer to pile groups with three or more rows in the load direction; they are equal to 0.8, 0.4, and 0.3 for the first, second, and third or higher row when S = 3D and equal to 1.0, 0.85, and 0.7 when S = 5D, respectively. A linear interpolation (shown in Fig. 3) was used to determine the p-multiplier for pile spacing contained between 3D and 5D. These values were found to be generally in good agreement with the p-multipliers obtained from the literature;

Table 2. Mechanical properties of concrete from Zhou et al. (2021a)⁶²

Table 3. Mechanical properties of steel reinforcement from Zhou et al. (2021a)⁶²

thus, they were used in the modeling performed in this study to estimate the p-multipliers. For the three-row pile groups with S = 2.5D, the p-multipliers were taken equal to the mean values obtained from Table 4, i.e., f_{m} = 0.66, 0.38, and 0.29 for the first, second, and third row of piles, respectively. For the two-row pile groups with a pile spacing S = 3D, the p-multipliers were taken equal to the mean values obtained from Table 4, i.e., f_{m} = 0.84 and 0.56 for the first and second row of piles, respectively.

Finite-Element Modeling

A finite-element (FE) model of the pile group foundation in sandy soil was built based on the BNWF approach to simulate the quasi-static test previously described. Fig. 4 illustrates the FE numerical model, which is developed and analyzed using the OpenSees software framework (McKenna 2011)⁴¹.

Modeling of Piles and Cap

The piles were modeled using displacement-based beam-column elements with distributed plasticity and fiber sections (Barbato et al. 2010)¹²; (Liu et al. 2020)³⁶. An FE mesh convergence analysis was performed in this study to determine an appropriate FE discretization for the piles. This analysis indicated that the FEMs with pile element lengths equal to 0.5D, 1.0D, and 1.25D produce a converged (i.e., almost identical) response for both global and local response quantities. Therefore, each pile was discretized into FEs of length equal to 0.12\mathrm{m} (i.e., 1D), except for the pile bottom element with a length of 0.22\mathrm{m}, as shown in Fig. 4(a). Each beam-column FE has five Gauss-Lobatto integration points. Different constitutive models were assigned to fibers corresponding to unconfined concrete, confined concrete, and longitudinal steel reinforcement. In particular, the axial stress-strain behavior of the concrete fibers was simulated using the uniaxial material Concrete01, which corresponds to the Kent-Scott-Park model with zero strength in tension (Scott et al. 1982)⁴⁹. This model has been shown to properly represent the stress-strain behavior of GIW-confined concrete (Zhou et al. 2021b)⁶⁴, (Zhou et al. 2022b)⁶⁵. Fig. 4(c) shows the backbone curves of this concrete model. The model parameters for confined and unconfined concrete are listed in Table 2. The axial stress-strain behavior of the longitudinal steel reinforcement fibers was modeled using the uniaxial material Steel02, which corresponds to the Menegotto-Pinto model with kinematic and isotropic strain hardening (Filippou et al. 1983)²⁵. The model parameters for the longitudinal steel reinforcement are given in Table 3. The pile cap was modeled using two elastic beam-column elements, and the cap bottom was connected with all six pile heads by elastic beam-column elements. The axial and flexure stiffnesses of these elastic elements were set equal to 1,000 times those of the pile elements to simulate an approximatively rigid link between the cap and the pile heads. The gravity load corresponding to the self-weight of each pile element was applied to the corresponding nodes, and the cap weight was applied to the cap center. A constant vertical load was imposed on the cap-top node to produce an axial load ratio of 5\% on each pile head section. All the degrees of freedom (DOFs) of the aboveground nodes were left unconstrained. For the belowground nodes, because the lateral loads were applied to the cap-center along the three-row pile direction only (i.e., along the global X-axis), the two translational DOFs in the XY plane were connected to zero-length elements representing the soil-pile interaction (which is described in the following section), the rotational DOF about the Z-axis was a free DOF, and the other three remaining DOFs (i.e., translation along the Z-axis and rotations about the X- and Y-axis) were fixed.

Modeling of Soil-Pile Interaction Considering the Pile Group Effect

This study introduced an innovative approach to model the pile group effect in soil-pile interaction systems subject to cyclic or dynamic loads. In fact, the approach commonly adopted in the literature describes the soil resistance to the lateral movement of a pile group through the use of a constant p-multiplier (i.e., the so-called group efficiency factor) equally applied to all piles of the group (Adeel et al. 2020)⁴; (Brown et al. 2001)¹⁸; (Lemnitzer et al. 2010)³⁵; (Liu et al. 2020)³⁶. This group efficiency factor is commonly calculated as the average value of the p-multipliers for different pile rows and is adopted because different piles alternate the roles of leading and trailing piles during cyclic loading or earthquake excitations. However, this constant reduction factor cannot correctly simulate the difference in lateral soil resistance among different rows in a pile group under cyclic loading, leading to inaccurate estimates of differential curvatures for piles in different rows.

This study proposed a new practical modeling method, which simulates the soil resistance in front of a pile at depth h by using two parallel springs consisting of (1) a common nonlinear p-y spring and (2) a nonlinear asymmetric spring, as illustrated in Fig. 4(e). In particular, the end nodes of each pile element below the ground surface were connected to the fixed nodes representing the soil site via zero-length elements. The load-displacement response of the zero-length element was described by a nonlinear p-y spring and a nonlinear asymmetric spring in parallel in the horizontal direction (X-axis) and by a nonlinear t-z spring in the vertical direction (Y-axis). The nonlinear p-y spring was modeled using the uniaxial material PySimple1 in OpenSees (Boulanger et al. 1999)¹⁴, which is commonly used to describe the force-displacement relation for soils acting on piles. The backbone of the p-y constitutive model for sand is given as follows (API 2007)⁸; (Chai and Song 2012)²¹:

Table 4. Values available in the literature of p-multipliers for pile groups with three or two rows in the load direction

Note: \mathrm{NS} = not specified; and - = not present.

where p = lateral resistance of soil at the embedded depth h; y = lateral deflection of the pile at depth h; p_u = ultimate resistance of the sand at depth h; A = loading factor, which is equal to 0.9 for cyclic loading; p_{us} and p_{ud} = ultimate resistances of soil in the shallow and deep regions, respectively; n_h = initial subgrade reaction modulus of sand, which can be obtained from API specification as a function of the sand friction angle; \gamma = soil weight density; and C_1, C_2, and C_3 = nondimensional coefficients that depend on the effective friction angle (API 2007)⁸; (Chai and Song 2012)²¹.

Fig. 3. p-Multipliers for three-row pile groups in the sand.

The nonlinear asymmetric spring is approximatively modeled using the uniaxial material QzSimple1 in OpenSees (Boulanger et al. 1999)¹⁴, which has a behavior similar to the constitutive model for the p-y spring on the compression side but has an asymmetric and smaller soil strength on the tension side, as shown in Fig. 4(g). The parameters of the backbone q-z curve are adjusted to approximately reproduce the same backbone curve used for the p-y spring, whereas the suction factor is set equal to zero. This nonlinear asymmetric spring is used to model the asymmetric value of the p-multiplier of any given pile when a pile group is subject to cyclic or seismic loads, i.e., when two different p-multiplier values need to be applied to the same pile in a given row that is switching from leading (corresponding to the larger p-multiplier value, f_{m,t}) to trailing pile (corresponding to the smaller p-multiplier value, f_{m,t}) as the load changes direction. The asymmetric spring is oriented so that the compression side coincides with the side in which the pile is nontrailing. It is noted here that a single asymmetric p-y spring (with two different p-multipliers for the leading and trailing directions) could be used to produce the same behavior obtained through the combination of the symmetric p-y spring and the

Fig. 4. Numerical modeling: (a) schematic illustration of the entire model; (b) fiber section discretization of the piles; (c) concrete model; (d) steel model; (e) soil-pile interaction modeling; (f) p-y spring model; (g) asymmetric spring model; and (h) ENT spring model.

asymmetric q-z spring proposed in this study. However, such a constitutive model is not currently available in OpenSees. Thus, the soil resistances of the two lateral parallel springs are given as follows:

where p_{sym}^{(n)} = \mathrm{soil} resistance for the nth row piles provided by the symmetric p-y spring; p_{asym}^{(n)} = \mathrm{soil} resistance for the nth row piles provided by the asymmetric spring; superscript n = 1,2,\ldots ,n_{\mathrm{max}} = pile row number; and n_{\mathrm{max}} = total number of rows. As a result, the required input parameters of the PySimple1 material in OpenSees, p_{\mathrm{ult}} and y_{50} (Blanco et al. 2019)¹³, are given by

where p_{\mathrm{ult}} = ultimate soil resistance provided by the symmetric p-y spring; y_{50} = soil displacement at 50\% of p_{\mathrm{ult}} and L_{t} = tributary length of the soil-pile contact associated with the given node. The required input parameters of the QzSimple1 material, q_{\mathrm{ult}} and z_{50} are given by

where q_{\mathrm{ult}} = ultimate soil resistance provided by the asymmetric q-z spring; and z_{50} = soil displacement at 50\% of q_{\mathrm{ult}}. It is noted here that, for cases in which the p-multiplier value for a pile row remains constant in the two opposite loading directions, the soil resistance corresponding to this pile row can be modeled more simply by using only the symmetric p-y springs with the appropriate value of the p-multiplier. For the numerical model of the quasi-static test considered in this study, which involves a pile group foundation with three rows of piles in the loading direction with a 3D pile spacing, the p-multiplier values are f_{m,l}^{(1)} = f_{m,t}^{(1)} = 0.8, f_{m,t}^{(1)} = f_{m,t}^{(3)} = 0.3, and f_{m}^{(2)} = f_{m,l}^{(2)} = f_{m,t}^{(2)} = 0.4.

The vertical soil-pile friction behavior was simulated using a t-z spring modeled with the TzSimple1 material in OpenSees (Boulanger et al. 1999)¹⁴. The corresponding input parameters t_{\mathrm{ult}} and z_{50} are given by (Mosher 1984)⁴⁴

where t_{\mathrm{ult}} = ultimate friction force at the soil-pile interface within the tributary length L_{t}, k_{0} = coefficient of lateral earth pressure at rest and is set equal to 0.4; \phi = friction angle of sand; z_{50} = displacement at which the friction force reaches 50\% of t_{\mathrm{ult}} and k = initial tangent stiffness and can be expressed as a function of the friction angle (Mosher 1984)⁴⁴. Finally, the FE model of the benchmark example used for validation describes the boundary conditions at the pile tips with vertical springs, the behavior of which is given by an elastic-no-tension (ENT) material with an initial stiffness of 1\times 10^{7}\mathrm{kN / m} as shown in Fig. 4(h). The modeling parameter values for the quasi-static test used in Eqs. (2)-(13) are given in Table 5.

Table 5. Modeling parameters for the sand soil in the benchmark example

Numerical Model Validation

The proposed FEM approach for pile group effect modeling under cyclic loading, referred to as the proposed cyclic model (PCM) hereinafter, was validated through a comparison with the experimental results available in Zhou et al. (2021a)⁶². To assess the performance of this approach with those commonly used in the literature, two additional FE models were built in OpenSees. The first additional FE model, referred to as ordinary monotonic model (OMM), uses ordinary symmetric p-\nu springs with p-multipliers equal to 0.8, 0.4, and 0.3 for leading, middle, and trailing piles, respectively, and is subjected to a monotonic pushover analysis with a maximum lateral displacement equal to 100\mathrm{mm}. These p-multiplier values are equal to those recommended in AASHTO (2020)¹ for a three-row pile group in sandy soil with a spacing of 3D. The second additional FE model, referred to as the ordinary cyclic model (OCM), adopts a constant p-multiplier (i.e., group efficiency factor) equal to 0.5 applied to all piles of the group and is subjected to a quasi-static cyclic loading. The value of the group efficiency factor was calculated as the average value of the p-multipliers for different pile rows, as recommended by Brown et al. (2001)¹⁸. All other modeling details were identical for the three considered FE models. The lateral load was applied to the cap center using displacement-controlled loading, and the nonlinear residual equations of equilibrium were solved using the Newton-Raphson algorithm with command Newton in OpenSees (Mazzoni et al. 2006)⁴⁰. It is observed here that all three modeling methods for the pile group effect are based on the BNWF model, which is commonly used in practical applications. In addition, the computational effort associated with both PCM and OCM is almost identical; in fact, the clock time for both cyclic analyses were approximately 1,280 and 1,180\mathrm{s}, respectively, when using a personal computer with an Intel Core-i7-10750H CPU @ 2.60 GHz and 32 GB RAM.

Fig. 5 presents the global force-displacement curve comparisons between the experimentally measured and numerically simulated results. In particular, Figs. 5(a and b) compare the global hysteretic force-displacement curves predicted by the PCM and OCM, respectively, with the experimental result, whereas Fig. 5(c) compares the experimental backbone curve corresponding to positive displacements with the lateral force-displacement curve predicted by the three different numerical models used in this study. It was observed that both the PCM and the OCM provide an overall very good agreement with the cyclic experimental results. The comparison of the global lateral force-displacement backbone curve results indicated that the three FE models used in this study provide almost identical results in terms of global response quantities, which are in good agreement with the corresponding experimental results. This result (i.e., negligible differences in global response quantities) was expected because the three pile group effect modeling approaches used in this study differ only in the way the lateral soil resistance is distributed among different pile rows, which has negligible effects on the overall soil resistance exerted on the entire pile group.

Fig. 5. Comparison of experimentally measured and numerically predicted lateral force-displacement results: (a) cyclic response predicted by PCM; (b) cyclic response predicted by OCM; and (c) backbone curve for positive displacements.

Fig. 6 compares the experimentally measured and numerically predicted pile curvatures. As reported in Zhou et al. (2021a)⁶², the experimental values of the curvatures along the piles were obtained from strain gauges (identified by circles in Fig. 6), except for the values at the pile heads for displacement values of 30 and 50~\mathrm{mm} which were obtained from the data measured by linear potentiometers (identified by crosses in Fig. 6). To compare consistent experimental and numerical values, the numerical values of the curvature were obtained as the cross-sectional curvatures at the location of the strain gauges for the experimental curvatures obtained by using the strain gauges, whereas they were obtained as the average curvature for the finite elements corresponding to the length of the linear potentiometers for the experimental curvature obtained by using the linear potentiometers. To quantify the accuracy of the different numerical results, Table 6 reports the normalized root mean squared error (Rizzo et al. 2018)⁴⁷ for different displacement levels, which is calculated as follows:

Fig. 6. Comparison of experimentally measured and numerically predicted pile curvatures. [Note: Exp. (s.g.) = experimental values obtained from strain gauges, Exp. (l.p.) = experimental values obtained from linear potentiometers.]

where \phi_{\mathrm{exp},i} and \phi_{FE,i} = experimentally measured and numerically predicted pile curvatures, respectively, at location i = 1 2, ..., N along the pile; N = total number of the experimental curvature data points collected along a given pile; and \mathrm{FE} = \mathrm{PCM} OCM,or \mathrm{OMM} = considered FE model. It is observed that the predicted pile curvatures by the PCM and OMM are almost identical (i.e., |\epsilon_{\mathrm{PCM}} - \epsilon_{\mathrm{OMM}}|\leq 2.6\% for all considered displacement levels and all pile rows) and present a good agreement with the experimental results (i.e., 2.1\% \leq \epsilon_{\mathrm{PCM}}\leq 10.2\% and 2.7\% \leq \epsilon_{\mathrm{OMM}}\leq 10.0\%). When the pile group was subjected to

Table 6. Normalized root-mean-squared error for pile curvatures at different lateral displacements

cyclic loads, the PCM provided better predictions of the pile curvatures than the OCM, with \epsilon_{\mathrm{PCM}} < \epsilon_{\mathrm{OCM}}, except for leading piles at 10 \mathrm{mm} displacement (for which \epsilon_{\mathrm{PCM}} = 8.9\% and \epsilon_{\mathrm{OCM}} = 5.6\%). In particular, the OCM underestimates the curvature in the belowground plastic hinge regions and overestimates the embedded depth of the plastic hinge for the leading piles, whereas it overestimates the curvature in the belowground plastic hinge regions and underestimates the embedded depth of the plastic hinge for the trailing piles.

Parametric Sensitivity Analysis of Pile Group Foundations in Sandy Soils

Parametric Study Matrix

This study presented results of an in-depth parametric analysis to quantify the seismic performance sensitivity of scoured pile group foundations. The parametric analysis was based on the validated FE modeling method for pile group response under lateral

loading and monotonic pushover analysis. Pushover analysis was preferred over cyclic analysis to reduce the computational effort and because the validation results presented in the previous section of this paper showed that both PCM and OMM produce practically identical global and local responses that are in good agreement with experimental results.

Eleven modeling parameters were investigated in this study: seven geometric parameters, i.e., pile configuration (P.C.), pile length (L_{p}), pile diameter (D), pile center-to-center spacing (S), scour depth (L_{a}), longitudinal steel reinforcement ratio (\rho_{l}), and transverse steel reinforcement ratio (\rho_{s}); one loading parameter, i.e., axial load ratio of the piles (\eta); two material parameters, i.e., concrete strength (f_{c}) and yield strength of the steel reinforcement (f_{y}); and one soil parameter, i.e., relative density of sand (D_r). A central composite design method was used to select the parameter combinations to be investigated, as shown in Table 7. For these combinations, C0 was selected as the central point, and every other combination was obtained by varying one parameter at a time, with three different values for each parameter. These parameter values were selected based on typical values encountered in practice and in the literature (Aviram et al. 2008)¹⁰; (Blanco et al. 2019)¹³; (Das 2002)²³; (Jones et al. 2002)³⁰. As a result, 23 cases are considered in this study.

Table 7. Parameter matrix and studied cases

Note: The boldface values for each parameter denote the middle, lower, and upper bounds, respectively.

The modeling details and the constitutive model of the cover concrete and the longitudinal steel rebars were identical to those used in the FE model of the quasi-static test. The elastic modulus and strain hardening ratio of the reinforcement were 201 GPa and 0.83\% (Zhou et al. 2022a)⁶³, respectively. To model steel reinforcement with spirals, the core concrete of the piles was modeled by using the uniaxial material Concrete04 available in OpenSees, which corresponds to the Mander model with zero tension strength (Mander et al. 1988)³⁸. Table 8 summarizes the values of the constitutive parameters used in the parametric study to model the core concrete. The modeling approach used here for the lateral and vertical soil-pile interaction was identical to that described for the OMM.

Ductility Development and Performance Limit States of Scoured Pile Groups

The piles in a scoured pile group form a frame-like structure due to the constraints imposed by the surrounding soil and the pile-cap connection. Multiple plastic hinges can develop on the pile shafts, both above and below ground levels, when subjected to lateral loads (Liu et al. 2020)³⁶; (Wang et al. 2016)⁵⁶; (Zhou et al. 2021a)⁶². The plastic hinge development can be used to define limit states for both load and resistance factor design and performance-based design and to assess the safety and need for retrofit of damaged bridges after an earthquake event. Based on the structural behavior identified in the existing literature (Blanco et al. 2019)¹³; (Liu et al. 2020)³⁶; (Wang et al. 2016)⁵⁶; (Zhou et al. 2021a)⁶², five performance limit states for a pile group are considered in this study: (1) first aboveground yielding (FAY), corresponding to the first yielding of any longitudinal steel reinforcement in the pile group, which for scoured piles is generally located at the interface between cap and leading piles; (2) first belowground yielding (FBY), corresponding to the easy-to-repair limit state introduced by Blanco

Table 8. Core concrete parameters used in the studied cases

et al. (2019)¹³, so called because the pile damage observed before reaching this level is limited to only the aboveground portion of the piles, which is easily accessible for postearthquake inspection and repair; (3) peak lateral strength (PLS), beyond which lateral strength degradation initiates; (4) severe structural damage (SSD), which is identified by the core concrete crushing or the rupture of the longitudinal steel reinforcement in any pile within the pile group, whichever occurs first; and (5) ultimate residual strength (URS), defined here as a 15\% reduction of the lateral strength from its peak value (Ataei and Padgett 2012)⁹; (Shen et al. 2021)⁵¹. A multiple-level displacement ductility index for a pile group is defined here as follows:

where \Delta_{i} = horizontal displacement of the cap center corresponding to any specific damage state i = 1,2,3,4, and 5 (i.e., FAY, FBY, PLS, SSD, and URS limit states, respectively).

Fig. 7 plots the backbone curves with markers to identify the displacements corresponding to the considered limit states for all 23 cases. For the 2 \times 3 pile groups, the FAY and FBY occur on the leading piles, consistent with the findings of previous experimental studies (Liu et al. 2020)³⁶; (Zhou et al. 2021a)⁶², whereas the FAY and FBY occur on the trailing pile for the 2 \times 2 pile group (i.e., case C1). The ductility development for the 3 \times 3 pile groups is very similar to that for the 2 \times 3 pile groups, with a lateral resistance equal to 1.5 times that of the corresponding 2 \times 3 pile groups, as the pile group effect is identical for these two cases. For all the considered cases, the SSD corresponds to the condition of concrete crushing at the head of the leading piles, which always occurs before the rupture of the longitudinal steel reinforcement. The SSD limit state always occurs before the URS limit state, except for cases C11 and C14. This phenomenon is observed because (1) for C11, the low longitudinal steel reinforcement ratio induces a rapid postpeak strength reduction, thus reaching the URS limit state before the core concrete can reach its ultimate strain; and (2) for C14, the high transverse reinforcement ratio significantly increases the ultimate strain of the core concrete, thus delaying its crushing limit state.

The pile length has a negligible impact on the lateral resistance of the pile group and the displacements corresponding to the different limit states. It is noted that this conclusion is valid here because the piles behave as deep foundations without uplift. As expected, the lateral strength of a pile group significantly increases with increasing pile diameters. More specifically, the peak lateral strength of the studied pile groups is equal to 1,511, 7,527, and 18,788\mathrm{kN} for D = 0.6, 1.2, and 1.8\mathrm{m}, respectively. The pile spacing has a slight effect on the peak lateral strength and the limit state displacements of a pile group. An increasing scour depth produces a reduction of the peak lateral strength of the pile groups and an increase of the displacement values corresponding to different limit states. In particular, the peak lateral strength of the pile groups is equal to 9,988, 7,527, and 5,860\mathrm{kN} for a scour depth equal to 3D, 5D, and 7D (with D = 1.2\mathrm{m}), respectively. Increasing the longitudinal steel reinforcement ratio can enhance the lateral resistance of a pile group and increase the displacement corresponding to the first three limit states (FAY, FBY, and PLS). By contrast, the transverse steel reinforcement ratio has a negligible effect on the global lateral load-displacement response of the pile group before PLS, but it becomes significant in controlling the strength degradation after the peak. An increasing axial load ratio slightly increases the peak lateral strength of a pile group; however, it also accelerates the crushing of the core concrete and triggers a rapid degradation of the lateral strength, which is an undesirable behavior for structures subjected to earthquake excitations. Increasing the concrete strength can improve the lateral strength of a pile group, but it results in a more rapid degradation of the lateral strength. By contrast, increasing the yield strength of steel reinforcement can delay the lateral strength degradation after peak, but it is not the most desirable approach for improving the lateral resistance of a pile group. The relative density of sand has very small effects on the peak lateral strength and the limit state displacements of a pile group.

Sensitivity Rankings of the Studied Parameters

The sensitivity rankings of the 11 parameters considered in this study in terms of seismic performance of a scoured pile group foundation were investigated by using the tornado diagram method (Barbato et al. 2010)¹². The seismic performance of a scoured pile group is assessed with respect to three aspects, i.e.: (1) global resistance to an earthquake excitation; (2) residual displacement corresponding to the FBY limit state, which is directly related to the postearthquake repair costs; and (3) ductility capacity. The parameter sensitivity for a given response R is quantified by the total relative swing, sw_{R}, which is defined as follows:

where R_{0} = value of the response quantity when the considered parameter is equal to the middle (reference) value; R^{(- )} and R^{(+)} = values of the response quantity when the considered parameter is assumed to equal to the lower and upper bounds, respectively, which are given in Table 7; and sw_{R}^{(- )} = \Delta R^{(- )} / R_{0} and sw_{R}^{(+)} = \Delta R^{(+)} / R_{0} = lower and upper relative swing, respectively, which are also referred to as one-side relative swings hereinafter. These calculated total relative swings are sorted from high to low values and plotted from top to bottom to form the corresponding tornado diagram.

The global resistance of a pile group was quantified by three indexes, i.e. (1) the lateral strength at the FAY limit state, referred to as yield strength of a pile group hereinafter, beyond which the yielding of a pile group initiates; (2) the strength enhancement coefficient of a pile group after yielding, denoted SE, and defined as follows:

where F_{\mathrm{FAY}} and F_{\mathrm{PLS}} = lateral strengths at FAY and PLS limit states, respectively, which represents the capacity of a pile group after yielding; and (3) the normalized strength degradation rate after peak, SDR, which is proposed in this study for the first time and is defined as the ratio between the postyield stiffness and the

secant stiffness at yielding of a pile group as follows:

where F_{\mathrm{URS}} = lateral strength at the URS limit state; and \Delta_{\mathrm{FAY}} and \Delta_{\mathrm{URS}} = displacement at the FAY and URS limit states, respectively.

Fig. 7. Sensitivity of lateral capacity to parameter variability: (a) pile configuration; (b) pile length; (c) pile diameter; (d) pile spacing; (e) scour depth; (f) longitudinal steel reinforcement ratio; (g) transverse reinforcement ratio; (h) axial load ratio; (i) concrete strength; (j) yield strength of reinforcements; and (k) soil relative density.

Fig. 8 shows the sensitivity rankings of the 11 parameters with respect to the global seismic resistance of the studied pile groups. The parameters below the horizontal dashed line have a small effect on the studied response with a maximum one-side relative swing of less than 0.05. It was observed that the yield strength of a pile group is most sensitive to the pile diameter (with a relative swing equal to 2.19), whereas the strength enhancement coefficient and normalized strength degradation rate were most sensitive to the axial load ratio (with a relative swing equal to 0.19 and 0.85, respectively). Increasing the pile diameter could significantly improve the yield strength, raise the strength enhancement coefficient, and delay the strength degradation of a pile group. Increasing the axial load ratio slightly enhanced the yield strength of a pile group, but it considerably accelerated the strength degradation and reduced the strength enhancement coefficient of a pile group, which is an undesirable effect for the scoured pile groups. The scour was also undesirable for a pile group because it significantly weakened the pile group (e.g., the yield strength decreased by 36% in this study when the scour depth increased from 3D to 7D), decreased the strength enhancement coefficient, and accelerated the strength degradation of a pile group. The transverse steel reinforcement ratio had a small effect on the yield strength and the strength enhancement coefficient of a pile group, but it could significantly delay the lateral strength degradation of a pile group after peak strength is reached (with a relative swing equal to 0.73). In addition, the global resistance of a pile group was mostly insensitive to the relative density of sandy soil and the pile length. Compared with the yield strength and the normalized strength degradation rate, the strength enhancement coefficient was relatively stable and showed less sensitivity to the studied parameters because its maximum one-side relative swing was only 0.12.

Fig. 8. Parameter sensitivity rankings for the global resistance: (a) lateral strength at FAY limit state; (b) strength enhancement coefficient of a pile group after yielding; and (c) normalized strength degradation rate after the peak.

Fig. 9 presents the parameter sensitivity ranking for the residual displacement of a pile group at FBY limit state (denoted as DR_{\mathrm{FBY}} hereinafter). The residual displacement was obtained by the following two analysis steps using the PCM: (1) load the model to determine the lateral displacement of the pile group at the FBY limit state and (2) unload the model until it reaches a zero-lateral force state. The displacement corresponding to the zero-lateral force state was regarded as the residual displacement of the pile group (Zhou et al. 2021a)⁶². As shown in Fig. 9, the top five parameters to which the residual displacement of a pile group at the FBY limit state is sensitive are (in decreasing order of importance) (1) axial load ratio, (2) yield strength of the steel reinforcement, (3) longitudinal steel reinforcement ratio, (4) pile diameter, and (5) scour depth in sequence, corresponding to the relative swings of 0.91, 0.54, 0.51, 0.49, and 0.26, respectively. Their increase leads to the increase of the residual displacement of a pile group, which is detrimental to the postearthquake recovery efforts. The pile length and transverse steel reinforcement ratio have less effect on the residual displacement of a pile group corresponding to the FBY state since the maximum one-side relative swing is less than 0.05.

Fig. 9. Parameter sensitivity ranking for the residual displacement at the FBY limit state.

Fig. 10 shows the sensitivity rankings of the different parameters with respect to the displacement ductility of a pile group at the FBY, PLS, and SSD limit states. The top four parameters to which the FBY displacement ductility of a pile group is sensitive are (in decreasing order of importance) P.C., scour depth, pile center-to-center spacing, and pile diameter. The other remaining seven parameters have a small effect on the FBY displacement ductility since their maximum one-side relative swing is less than 0.05. The 2 \times 3 and 3 \times 3 pile groups present an almost identical FBY displacement ductility since they have the same pile group effect. The axial load ratio affects the most the displacement ductility of a pile group at the PLS and SSD states (with a relative swing equal to 0.55 and 1.18, respectively), but it has a very small impact on the ductility of a pile group at the FBY state (with a relative swing smaller than 0.006). The higher the axial load ratio causes a lower pile group ductility corresponding to PLS and SSD states. The increase in the scour depth reduces the ductility of a pile group at the FBY and PLS states, but it has almost no impact on the ductility of a pile group at the SSD state. Increasing the transverse steel reinforcement ratio improves the ductility of a pile group at the SSD limit state; however, the transverse steel reinforcement ratio has a negligible impact on the ductility of the pile group at the FBY and PLS limit states. The displacement ductility of a pile group is insensitive to the pile length and the relative density of sand.

Considerations for Seismic Performance Improvements

Fig. 10. Parameter sensitivity rankings for the displacement ductility at (a) FBY; (b) PLS; and (c) SSD limit states.

The parametric sensitivity analysis results can be used to provide meaningful recommendations for the seismic design and retrofit of pile group foundations in sandy soils subject to scour effects, as summarized in Table 9. In this table, the plus sign (+) indicates that increasing the corresponding parameter value improves the seismic performance of the pile group (i.e., having a positive effect); the cross mark (\times) represents that the increase of the parameter value deteriorates the performance (i.e., having a negative effect); the horizontal line (- ) indicates that the response is not monotonically affected by the parameter; the slash mark (\backslash) indicates that the corresponding parameter has a generally small effect on the seismic performance (i.e., both lower and upper relative swings are less than 0.05). In addition, the marks corresponding to the three parameters affecting the most each considered response are highlighted in bold, and the most important and second most important parameters for each considered response are identified by three and two marks, respectively.

Table 9. Seismic performance variation of a pile group by increasing the parameter value

In general, increasing the pile diameter is the most efficient measure to improve the seismic performance of a pile group when considering scour effects because it can significantly increase the yield strength and the strength enhancement coefficient of a pile group and reduce the normalized strength degradation rate. Increasing the axial load ratio deteriorates the seismic performance of a scoured pile group because it (1) increases its residual displacement at the FBY state, (2) accelerates the lateral strength degradation after peak strength, and (3) degrades the strength enhancement coefficient of the pile group. Increasing the center-to-center spacing between piles up to 5D has a small positive effect on the seismic performance of a scoured pile group; however, it requires increasing the cap size, resulting in an increase in construction costs. For a deep pile group foundation, the pile length and the relative density of sand have a negligible effect on the seismic performance of the pile group. The transverse reinforcement ratio has a small effect on a pile group's seismic performance before the peak lateral strength is achieved; however, the transverse steel reinforcement must be sufficient to prevent buckling of the longitudinal rebars, which becomes important for large inelastic behavior. The increase in the scour depth generally deteriorates the seismic performance of a pile group.

Conclusions

This study investigates the seismic performance of scoured pile group foundations. A new practical finite-element modeling approach for the pile group effect is first proposed and validated by the experimental data available in the literature. The proposed approach can simulate the soil resistance difference among different pile rows in a pile group under lateral cyclic loads and provides more accurate local pile curvature results than the currently adopted methodology based on a constant pile group factor. It is also observed that both approaches provide almost identical results in terms of global response quantities, e.g., global lateral force-displacement response, which is in good agreement with existing experimental results. By using the validated model, an in-depth parametric analysis is performed to explore parameter sensitivity with respect to the seismic performance of a pile group foundation. The main findings of this parameter sensitivity analysis are as follows:

1. Scour significantly weakens the seismic capacity of a pile group. It reduces the lateral strength and displacement ductility of a pile group and increases the residual displacement of a pile group at the first belowground yielding (FBY) limit state (corresponding to an easy-to-repair condition).
2. Increasing the pile diameter is the most effective measure to improve the seismic performance of a scoured pile group.

Increasing the concrete strength, the yield strength of the steel reinforcement, and the longitudinal steel reinforcement ratio can increase the yield strength of a pile group. However, increasing the yield strength of the steel reinforcement and the longitudinal steel reinforcement ratio can also increase the residual displacement of a pile group at the FBY state. A higher axial load ratio slightly increases the lateral strength of the pile group foundation but also significantly increases its residual displacement, reduces its lateral strength enhancement coefficient, and accelerates its lateral strength degradation. The transverse reinforcement ratio has negligible effects on the yield strength and the residual displacement at the FBY limit state of a pile group; however, a higher transverse reinforcement ratio can decrease the degradation rate of the lateral strength. For a deep foundation in medium-dense sand, the pile length (for values higher than or equal to 30 times the pile diameter), pile spacing (for values smaller than or equal to five times the pile diameter), and relative density of sand (for values between 40\% and 60\%) have negligible effects on seismic performance of a scoured pile group foundation.

This study uses a set of quasi-static tests available in the literature to validate the proposed pile group effect modeling approach for pile group foundations subjected to lateral cyclic loads. It is recommended to validate this modeling approach also for pile group foundations subjected to dynamic conditions, e.g., seismic shaking. This study mainly focuses on scoured pile groups in homogeneous sand subjected to lateral loads. Further studies are needed to quantify the impacts of different soil conditions and profiles. In addition, future numerical studies should investigate the use of incremental dynamic analysis to develop fragility curves and achieve a better understanding of the ductile behavior and the failure mechanisms of bridges with scoured pile group foundations. It is also suggested to investigate whether the p-multipliers vary with the depth and lateral deformation level of the pile group foundation and how this potential variation affects the performance of a scoured pile group foundation.

Data Availability Statement

Some or all data, models, or codes that support the findings of this study (including the section analysis executable codes and the force-displacement data for the pile group specimens) are available from the corresponding author upon reasonable request.

Acknowledgments

The authors gratefully acknowledge the funding support of this work by the National Natural Science Foundation of China (Grant Nos. 51778469 and 52178155) and the State Key Laboratory of Disaster Reduction in Civil Engineering, Ministry of Science and Technology of China (Grant No. SLDRCE19-B-20). The first author is thankful for the financial support from the China Scholarship Council. Any opinions, findings, and conclusions expressed are those of the authors and do not necessarily reflect those of the sponsoring organizations. Special thanks to Professor Ross W. Boulanger at the University of California, Davis who provided significant suggestions on the modeling of the pile group effect under cyclic loading. All staff members in the Multi-Functional Shaking Table Laboratory of Tongji University at Jiading Campus are also gratefully acknowledged.

References