Engineering, Technology & Applied Science Research Vol. 8, No. 1, 2018, 2577-2584 2577 www.etasr.com Soomro et al.: 3D Numerical Modeling of Pile Group Responses to Excavation-Induced Stress Release ... 3D Numerical Modeling of Pile Group Responses to Excavation-Induced Stress Release in Silty Clay Mukhtiar Ali Soomro Department of Civil Engineering Quaid-e-Awam University of Engineering, Science & Technology Nawabshah, Pakistan eng.soomro@gmail.com Abdul Salam Brohi Department of Civil Engineering Quaid-e-Awam University of Engineering, Science & Technology Nawabshah, Pakistan eng.armaan@gmail.com Mohsin Ali Soomro Department of Civil Engineering Quaid-e-Awam University of Engineering, Science & Technology Nawabshah, Pakistan drmohsin@quest.edu.pk Daddan Khan Bangwar Department of Civil Engineering Quaid-e-Awam University of Engineering, Science & Technology Nawabshah, Pakistan skb_khan2000@yahoo.com Shahzad Ali Bhatti Department of Civil Engineering Quaid-e-Awam University of Engineering, Science & Technology Nawabshah, Pakistan shahzadalibhatti72@gmail.com Abstract—Development of underground transportation systems consists of tunnels, basement construction excavations and cut and cover tunnels which may encounter existing pile groups during their construction. Since many previous studies mainly focus on the effects of excavations on single piles, settlement and load transfer mechanism of a pile group subjected to excavation- induced stress release are not well investigated and understood. To address these two issues, three-dimensional coupled- consolidation numerical analysis is conducted by using a hypoplastic model which takes small-strain stiffness into account. A non-linear pile group settlement was induced. This may be attributed to reduction of shaft resistance due to excavation induced stress release, the pile had to settle substantially to further mobilise end-bearing. Compared to the Sp of the pile group, induced settlement of the single pile is larger with similar settlement characteristics. Due to the additional settlement of the pile group, factor of safety for the pile group can be regarded as decreasing from 3.0 to 1.4, based on a displacement-based failure load criterion. Owing to non-uniform stress release, pile group tilted towards the excavation with value of 0.14%. Due to excavation-induced stress release and dragload, head load of rear piles was reduced and transferred to rear piles. This load transfer can increase the axial force in front piles by 94%. Keywords-multipropped excavation; pile group; tilting; load- transfer I. INTRODUCTION In dense urban environments where buildings are closely spaced, deep excavation for basement construction and other underground facilities such as mass rapid transit station sand cut-and-cover tunnels is unavoidable. As these excavations are usually carried out close to existing buildings, a major concern is to prevent/minimize the damage to adjacent buildings and underground utilities [1, 2]. Thus, it is important for designers to estimate potential damage resulting from nearby deep excavations to existing piles. Most researches estimated the buildings settlement and tilting considering wall movements and ground surface settlement trough using empirical approaches. The performance of pile group depends on the stress state in the soil and surrounding sub-surface soil movements [3]. In addition, deep excavation in soft clay induces negative excess pore water pressure [4] which induces long term pile group settlement with the dissipation of excess pore water pressure. Therefore, it is essential to investigate the pile group response mechanism adjacent to deep excavation in soft soil. Authors in [1, 2] reported case studies in granular soil and alluvium residual soil respectively. They demonstrated that lateral soil movements due to excavation can be detrimental to nearby existing piles. In both the reported case studies pile toe levels were much deeper than the excavation level and they reported only the lateral pile behavior. Apart from field monitoring, a number of centrifuge tests were also conducted to investigate the response of single pile [5] and pile group in soft kaolin clay [6]. They concluded that the induced bending moment and lateral deflection of piles were highly influenced by distance from wall and pile head condition. In studies, lateral response of end bearing piles without initial applied load was reported. In reality, pile group in soft clay behaves as floating pile group and subjected to initial applied load from superstructure. In the presence of initial applied load, soil surrounding the pile group experiences higher stress level before the commencement of adjacent excavation. Authors in Engineering, Technology & Applied Science Research Vol. 8, No. 1, 2018, 2577-2584 2578 www.etasr.com Soomro et al.: 3D Numerical Modeling of Pile Group Responses to Excavation-Induced Stress Release ... [7-8] conducted centrifuge tests to investigate the effects of an unpropped excavation on the behavior of nearby single piles and pile groups in dry dense Toyoura sand. It was found that the distance from the pile to the retaining wall and pile head conditions had a large influence on the induced pile bending moment and lateral deflection. Authors in [9] reported the results of three centrifuge tests which were carried out to study the effects of a multipropped deep excavation in-flight on the behavior of single piles in dry Toyoura sand. Piles were laterally restrained in terms of rotation and deflection right at or above ground surface in the three different tests. It was concluded that lateral restraints imposed on the pile head have a significant influence on induced pile bending moment. Induced bending moment due to excavation can exceed the pile bending capacity. Authors in [10] developed design charts to compute the lateral behavior of a single pile adjacent to deep excavation in soft ground. They performed two staged analysis considering plane strain conditions and linear elastic soil model. Similar work was conducted by authors in [11] using finite element method. In both these studies the lateral response of single pile was investigated. The settlement behavior of pile and development of excess pore water pressure and consolidation settlement were not investigated. Authors in [12] proposed an analytical method to investigate the capacity reduction and settlement increase of a nearby pile during excavation. It was reported that pile settlement due to excavation depends on the percentages of end bearing and shaft friction of the pile, the soil movement pattern, and the distribution of the maximum shaft friction with depth. However, shaft resistance in these methods is calculated on the assumption that horizontal stress acting on the piles does not change during excavation. This assumption may not be valid and the pile settlement may be underestimated using the preceding methods, leading to a non-conservative prediction. Most of the previous studies focused on the lateral response of single pile foundation. The vertical response of floating pile group combined with working load adjacent to deep excavation in soft clay has not been studied. In this study, to investigate the settlement, tilting, load transfer and load redistribution of a elevated (2×2) pile group affected by nearby excavation, a three-dimensional coupled- consolidation finite element analysis (FEA) was carried out in saturated silty clay. II. THREE-DIMENSIONAL COUPLED CONSOLIDATION ANALYSIS It is well-known fact that the stress-strain relationship of soils is highly nonlinear even at very small strain. The stiffness of most soils decreases as strain increases and depends on the recent stress or strain history of the soil [13-14]. Owing to non- linear soil behavior, an excavation can cause reduction in the ground stiffness. Therefore, it is vital to investigate the pile responses to adjacent excavation in silty clay. To obtain a satisfactory numerical model of the pile group responses to excavation-induced stress relief, the analysis needs to take account of the small strain non-linearity of soil. To gain new insights into pile group responses to a nearby multipropped excavation in saturated silty clay, this study conducts a 3D coupled consolidation numerical analysis. Figure 1(a) shows the elevation view of the configuration of numerical simulation in which a multipropped excavation was carried out adjacent to a loaded elevated (2×2) pile group. The final depth of the excavation (He) was 10m. The embedded length (Lp) and diameter (dp) of the piles were 18m and 08 m respectively. The heads of all four piles were rigidly embedded into a 4.5m × 4.5m ×1.0m (length×width×thickness) pile cap. The pile spacing was 3.2dp. The pile cap was elevated by 1m from the ground surface. The modeled pile represents a cylindrical reinforced concrete (grade 40, reinforcement ratio=1) with a bending moment capacity of 800kNm. The clear distance between diaphragm wall and the front pile row was 3.0m (0.3He). The excavation was supported by 0.6m thick diaphragm wall. The ratio of wall penetration depth to excavation depth is typically 0.5-0.2 in engineering practice [15-16], and thus a value of 0.5 was adopted in this study. The retaining wall was supported by three levels of props, at 1.0, 4.0 and 7.0 m depth. The props were modeled as soft with axial rigidity of 81 × 103kNm [15]. Horizontal spacing of props was 10m. Figure 1(b) illustrates the plan view of the configuration of the numerical simulation. The length of the excavation is 12m. Due to symmetry, only half of the excavation was simulated. A monitoring section was selected at the transverse centerline of the excavation. In addition to this simulation, a pile load test (L) was conducted numerically in “greenfield” conditions (i.e., without excavation) to obtain the ultimate capacity of the pile in silty clay. Based on this, the working load was then calculated with a factor of safety of 3.0. The obtained working load was applied to the pile in the analysis simulating excavation. Ground water table Working load (3.44 MN) Silty clay Prop 1 (1.0 m) Prop 2 (4.0 m) Prop 3 (7.0 m) P la n e o f sy m m et ry Excavation level 1 (3.0 m) Excavation level 2 (7.0 m) Formation level (10.0 m) Embedded length = 18 m Diameter = 0.8 m Diaphragm wall Thickness = 0.6 m Depth = 14.0 m 3.0 m (2×2) Pile group (4.5m×4.5m×1m) Pile cap (a) Monitoring section 50m 20m P la ne o f sy m m et ry 0.6 m thick diaphragm wall Prop Prop 10m 6m 1.0 m thick pile cap (b) Fig. 1. Configuration of the numerical run (a) elevation, (b) plan view Engineering, Technology & Applied Science Research Vol. 8, No. 1, 2018, 2577-2584 2579 www.etasr.com Soomro et al.: 3D Numerical Modeling of Pile Group Responses to Excavation-Induced Stress Release ... III. FINITE ELEMENT MESH AND BOUNDARY CONDITIONS Figure 2 shows an elevation view of a finite element mesh used to analyze the excavation effects on the pile group. The size of the mesh for numerical runs is 50m×20m×40m. These dimensions were sufficiently large to minimize boundary effects in the numerical simulation as further increment of the finite element mesh dimensions did not lead to any change in the computed results. Regarding the element size in the mesh, it is found that further halving the adopted mesh size only leads to a change of computed results of no more than 0.2%, suggesting the mesh is sufficiently fine. Eight-noded hexahedral brick elements were used to model the soil, the pile and the diaphragm wall, while two-noded truss elements were adopted to model the props. Roller and pin supports were applied to the vertical sides and the base of the mesh, respectively. Therefore, movements normal to the vertical boundaries and in all directions of the base were restrained. The water table was assumed to be at ground surface. Initially, the pore water pressure distribution was assumed to be hydrostatic. Free drainage was allowed at the top boundary of the mesh. The excavation process was simulated by deactivating soil elements inside excavation zone. In the meantime, the truss elements representing the props were activated. Fig. 2. Three dimensional finite element mesh and boundary conditions IV. CONSTITUTIVE MODEL AND MODEL PARAMETERS USED IN FINITE ELEMENT ANALYSIS The basic hypoplastic model was developed to capture the non-linear behavior (upon monotonic loading at medium to large-strain levels) of granular materials [19-20]. The basic model consists of five parameters  r and ,,, ** cN  . The parameters N and * define the position and the slope of the isotropic normal compression line in the ln(1+e) versus lnp’ plane [21], where e is the void ratio and p’ is mean effective stress. The parameter * defines the slope of the isotropic unloading line in the same plane. φc is the critical state friction angle and the parameter r controls the large strain shear modulus. To account for the strain-dependency and path- dependency of the soil stiffness (at small strains), authors in [22] further improved the basic hypoplastic model by incorporating the concept of intergranular strain. The intergranular strain concept requires five additional parameters  RTr mmR and ,,,  : R controlling the size of the elastic range, βr and χ controlling the rate of stiffness degradation. The parameters mT and mR control the initial shear modulus upon 180° and 90° strain path reversal, respectively. In this study, all the model parameters for silty clay reported in [23] are adopted and summarized in Table I. The coefficient of lateral earth pressure at rest, Ko is estimated through (1) [24]:      sin0 sin1 OCRK (1) TABLE I. SILTY CLAY MODEL PARAMETERS ADOPTED [23] Description Parameter Effective angle of shearing resistance at critical state: ’ 33o Parameter controlling the slope of the isotropic normal compression line in the ln(1 + e) versus lnp plane, * 0.103 Parameter controlling the slope of the isotropic normal compression line in the ln(1 + e) versus lnp plane, * 0.015 Parameter controlling the position of the isotropic normal compression line in the ln(1 +e) versus lnp plane, N 1.31 Parameter controlling the shear stiffness at medium- to large- strain levels, r 0.3 Parameter controlling initial shear modulus upon 180 strain path reversal, mR 12 Parameter controlling initial shear modulus upon 90 strain path reversal, mT 12 Size of elastic range, R 2×10-5 Parameter controlling the rate of degradation of the stiffness with strainr 0.09 Parameter controlling degradation rate of stiffness with strain 0.7 Initial void ratio, e 0.7 Dry density (kN/m3) 1615 Coefficient of permeability, k (m/s) 1×10-9 The concrete pile, the diaphragm wall and the props were assumed to be linear elastic with Young's modulus of 35GPa and Poisson's ratio of 0.25. The wall thickness was taken as 0.60m. The concrete unit weight was taken as 24kN/m3. The parameters for the piles and the diaphragm wall are summarized in Table II. TABLE II. CONCRETE PARAMETERS ADOPTED IN FINITE ELEMENT ANALYSIS Description Parameter Young's Modulus, E 35 GPa Poisson's ratio,  0.3 Density,  2400 kg/m3 V. NUMERICAL MODELING PROCEDURE Each numerical analysis is modeled according to the following steps: Step 1: Set up the initial boundary and initial stress conditions (i.e., static stress conditions with varying K0 with depth). Step 2: Activate the brick elements representing pile group (modeled as “wished-in-place”). Engineering, Technology & Applied Science Research Vol. 8, No. 1, 2018, 2577-2584 2580 www.etasr.com Soomro et al.: 3D Numerical Modeling of Pile Group Responses to Excavation-Induced Stress Release ... Step 3: Apply the working load (determined from numerical pile load test) on the pile. Step 4: Allow excess pore pressure, which generated in result of application of working load on the pile, to dissipate. Step 5: Activate the brick elements representing the diaphragm wall. Step 6: Staged multipropped excavation was simulated. After excavating to 3m depth, the first level of props was installed at 1m below ground surface. Soil was then excavated to 7m below ground surface, followed by the installation of the second level of props at 4m depth. Finally, excavation was extended to the target level of 10m depth with installation of the third level of props at 7m. VI. INTERPRETATION OF COMPUTED RESULTS A. Determination of Working Load for the Pile Group The study’s objective was to investigate a pile group (subjected to a working) load responses to an adjacent excavation. Prior to the simulation of the excavation, it was necessary to determine the ultimate axial load carrying capacity of the pile group. The working load can then be obtained using a factor of safety (FOS) of 3.0. Therefore, a numerical pile load test was carried out on a different finite element mesh to obtain the load settlement relationship and the capacity of the pile group without excavation. The load applied on the pile group was gradually increased to 11MN over a period of 24h. The resulting pile load-displacement curve for the simulated pile group is shown in Figure 3. The ultimate axial load capacity was determined based on a displacement-based failure criterion proposed in [25]. This failure criterion is expressed as follows: ,max 1 0.045 2 h p ph p p p P L d A E    (2) where δph,max is the maximum pile head movement which defines the ultimate load, Ph is the pile head load, Lp is the pile length, Ep is the pile shaft elastic modulus, Ap is the cross- sectional area of the pile, and dp is the pile diameter. Based on the failure criterion, the ultimate bearing capacity of 10.34MN was calculated. With a factor of safety (FOS) of 3.0, the working load was determined to be 3.44MN. Owing to the applied working load, the pile group settled by 0.85%dp (Figure 3). B. Progressive Pile Settlement and Apparent Loss of the Pile Capacity During Excavation Figure 4 shows the incremental settlement (Sp) of the pile with different excavation stages. Construction stages of the excavation are indicated by the depth (i.e., h) from the ground surface. Sp and h are normalized by the pile diameter (dp) and the final excavation depth (He) respectively. The measured induced single pile settlement due to excavation in centrifuge modeling reported in [9] and computed excavation-induced settlement of a single pile reported in [27] are also shown in the figure for comparison. The centrifuge test was carried out to investigate the effects of a multipropped deep excavation (final depth of excavation=8m in prototype) in-flight on the behavior of single pile (diameter of the pile=1.25 m in prototype) in dry Toyoura sand (i.e. Dr=70%). Fig. 3. Computed load settlement curve from the pile group load test without excavation Fig. 4. Normalised pile cap settlement during excavation The embedded depth of the pile and clear distance between the pile and model wall is 20m and 3m (in prototype) respectively. The configuration of their numerical analysis (i.e., excavation depth, pile diameter and length and clear distance between pile and diaphragm wall) are the same as those of the numerical model in this study. It can be observed from the figure that induced settlement with increasing excavation depth was approximately bilinear. During the first two excavation stages (i.e., h/He=0.3 and h/He=0.7), the pile settled almost linearly. After that, the pile settlement occurred at an increased rate during last excavation stage (h/He=1.0). This may be attributed to reduction of shaft resistance. Due to excavation induced stress release, the pile had to settle substantially to further mobilise end-bearing. Similar settlement characteristics of a single pile due to excavation observed from the centrifuge test reported in [9] and the computed results in [27]. However, induced settlement of the single pile was larger than that of the pile group after excavation. The total settlement of the pile (including settlement due to working load and excavation) was Engineering, Technology & Applied Science Research Vol. 8, No. 1, 2018, 2577-2584 2581 www.etasr.com Soomro et al.: 3D Numerical Modeling of Pile Group Responses to Excavation-Induced Stress Release ... 19.5mm (i.e., 2.5% pile diameter). This value satisfies a reliability-based serviceability criterion (i.e. 56mm), which was developed based on 95 composite (i.e., reinforced concrete and steel) buildings subjected to settlement [26]. This conclusion may not be applicable to scenarios in which the ground conditions or excavation depth are different from those adopted in this study. C. Reduction of Factor of Safety After Excavation The pile capacity is often determined using displacement based criteria. The induced settlement can be regarded (conveniently but approximately) as an additional load on the pile cap. To investigate the reduction of an equivalent factor of safety (FOS) of the loaded pile group due to tunneling, an additional pile cap settlement of 12.75mm (1.6% pile diameter) after excavation is considered as the cause of the reduction in pile capacity when a displacement-controlled failure load acceptance criterion as given in (1) [25] is used. Due to the initial applied working load (pile cap displaced vertically by 6.8mm or 0.85 pile diameter) and the excavation effects, the pile group is settled by 19.5mm in total. It is possible to consider that the pile cap is subjected to an equivalent axial load of 7.3MN including the excavation effects, resulting in a 19.5mm pile head settlement (see the pile-load displacement curve in Figure 3). As the applied working load of the pile increases from the original of 3.44MN to an equivalent applied working load of 7.3MN, this means that an equivalent FOS of the pile drops from 3.0 to 1.4 due to the excavation effect. However, it should be pointed out that the calculated reduction of the equivalent FOS does not really mean that the ultimate load capacity is physically reduced due to the excavation. The reduction of the equivalent FOS just simply means that designers should be reminded that the serviceability limit state of the pile group can be violated as a result of additional excavation-induced pile settlement. D. Transverse Tilting of Pile Cap During Excavation Figure 5 illustrates computed transverse tilting of the pile cap during excavation. Each tilting is determined as a ratio of differential settlements measured at two edges of the pile cap to the distance between them. A positive value means the pile cap tilts towards the first tunnel and vice versa. It can be seen from the Figure that during the excavation, the transverse titling of the pile group increases non-linearly with tunneling stages. The rate of the tilting increases as during entire process of the excavation. This is because the front pile row experienced larger stress release that that of the rear pile row. The computed tilting of the pile cap during the first, second and final excavation stage (i.e. h/He equal to 0.3, 0.7 and 1.0) is 0.02%, 0.07% and 0.14% respectively. All the values do not exceed the allowable tilting limit (0.2%) suggested by Eurocode 7 for typical commercial and residential buildings. Owing to same stress release along longitudinal direction of the excavation, no longitudinal tilting of pile cap was observed. E. Computed Load Redistributions among Piles in the Group In order to study the load transfer mechanism during excavation, Figure 6 shows the changes in load at the pile head within the 2×2 pile group at various excavation stages. In the figure, the change in head load due to tunneling (p) is normalized by the computed head load of each pile before excavation (pi). Before excavation, the applied working load is equally shared by each pile. It can be observed from the figure that the head loads carried by the front row piles increased linearly with excavation stages. It was expected that head loads carried by the front row pile would decrease due to excavation- induced stress release as the pile cap tilted towards excavation. On the contrary, the head loads of piles P1 and P2 decreased due to substantial drag load resulting from excavation-induced soil movement (discussed in section F). To support the constant working load acting on the pile cap, the head load of pile P2 as well as that of piles P4 decreased as a result of load redistribution from piles P1 and P3. Owing to load redistribution among piles, the front row piles experienced the most significant increase of 32% in head load. 0.00 0.05 0.10 0.15 0.20 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0 T il ti n g : % Normalised excavation depth (h/He) Pile group Pile group Ground surface Forma tion level Wall P1 P2 P3 Fig. 5. Computed transverse tilting of the pile cap during excavation F. Changes in Axial Load Distribution Since load redistribution in the front row of piles (piles P1 and P3) and rear row of piles (piles P2 and Pile P4) is the same, only one pile is selected from each row to discuss axial load distribution along the pile length. Figure 7 shows the axial force distribution along the front pile (P1) and rear pile (P2) with normalized depth (Z/Lp) below the ground surface after completion of the excavation. The axial load distribution before excavation (after applying the working load) is also included in the figure as a reference. The applied working load on pile cap including weight of pile cap 4,000kN, was shared equally by each pile before excavation. The load taken by each pile is 1000kN. Before excavation, the pile carried approximately 61% of the working load (i.e., 607kN) with its shaft resistance and the remainder with its end-bearing resistance. The load distribution along the pile shaft started to alter when the excavation commenced. It can be seen from the figure that the computed load increased along the entire pile length of P1 after completion of the excavation. On the other hand, the load decreased along the upper portion (Z/Lp≥0.3) of P2 and increased along the remaining length of the pile at the end of the excavation. The maximum increment in the axial force 94% and 33% of that Engineering, Technology & Applied Science Research Vol. 8, No. 1, 2018, 2577-2584 2582 www.etasr.com Soomro et al.: 3D Numerical Modeling of Pile Group Responses to Excavation-Induced Stress Release ... working load was computed in P1 at Z/Lp=0.43 and in P2 at Z/Lp=0.51 respectively. By inspecting the axial load distribution it is observed that along the upper portion of both piles (Z/Lp≥0.47), negative skin friction (NSF) is mobilized resulting from soil movement due to excavation induced-stress release. This suggests that this portion of the pile is subjected to “dragload” by the surrounding soil. Consequently, the load transferred to the lower portion of the pile as well as load redistribution occurred within pile group. To maintain equilibrium, the pile had to settle (see Figure 4) to further mobilize the end-bearing and shaft resistance along the lower portion of the piles (Z/Lp>0.47). The percentage increment of the end-bearing of P1 and P2 were computed as 74% and 37% respectively. -40 -30 -20 -10 0 10 20 30 40 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0 N o rm al is ed c h an ge o f p il e h ea d lo ad ( Δ p /p i): % Normalised excavation depth (h/He) Pile P1 Pile P2 Pile P3 Pile P4 E x ca v at io n Pile P1 Pile P2 Pile P3 Pile P4 Fig. 6. Computed load-redistributions among piles (P1, P2, P3 and P4) in the group 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0 300 500 700 900 1100 1300 1500 1700 1900 N o rm al is ed p il e d ep th ( Z /L p ) Axial force: kN At working load Pile P1 Pile P2 E x ca v at io n Pile P1 Pile P2 Fig. 7. Computed axial load distribution along the length of piles P1 and P2 due to excavation G. Mobilized Shaft Resistance during Excavation To further understand the load transfer mechanism, the mobilized shaft resistance in pile P1 and P2 after the application of the working load (i.e. before the excavation) and during the excavation is shown in Figure 8. The depth below the ground surface (Z) is normalized by the pile length (Lp). The computed average mobilized unit shaft resistance f(Z) at various depths was calculated based on the following equation: ( ) ( ) ZF f Z s    (3) where ΔF is the difference between the computed axial loads at two consecutive depths, ΔZ is the vertical distance between the two consecutive depths, and s is the pile perimeter. As expected, after working load application, the positive shaft resistance (i.e. soil supporting the pile) is mobilized along the pile entire length. The mobilized positive shaft resistance increased along pile depth with maximum shaft resistance of 26kPa near the pile toe. On the completion of the excavation, the mobilized shaft resistance at the upper portion of the both piles (Z/Lp<0.5) decreased to negative and increased below formation level (0.5