Theoretical t-z Curves for Axially Loaded Piles | Journal of Geotechnical and Geoenvironmental Engineering | Vol 148, No 7

AbstractEstimation of nonlinear pile settlement can be simplified using one-dimensional “t-z” curves that conveniently divide the soil into multiple horizontal “slices.” This simplification reduces the continuum analysis to a two-point boundary-value problem of the Winkler type, which can be tackled by standard numerical procedures. Theoretical “t-z” curves can be established using the “shearing-of-concentric-cylinders” theory of Cooke and Randolph-Wroth, which involves two main elements: (1) a constitutive model cast in flexibility form, γ=γ(τ); and (2) an attenuation function of shear stress with radial distance from the pile, τ=τ(r). Soil settlement can then be determined by integrating shear strains over the radial coordinate, which often leads to closed-form solutions. Despite the simplicity and physical appeal of the method, only a few theoretical “t-z” curves are available in the literature. This paper introduces three novel attenuation functions for shear stresses, inspired by continuum solutions, which are employed in conjunction with eight soil constitutive models leading to a set of 32 “t-z” curves. Illustrative examples of pile settlement calculation in two soil types are presented to demonstrate application of the method.IntroductionNonlinear settlement estimation of axially loaded piles requires continuum analysis using constitutive models in multiple dimensions, which might be difficult to implement in practice and are unavailable for many real soils. The problem can be simplified to a single dimension using the Winkler model of soil reaction through the introduction of pertinent stress-displacement curves (“t-z” curves) along the pile shaft (Seed and Reese 1957; Coyle and Reese 1966) and a stress-displacement curve (“q-z” curve) at the pile tip (Frank and Zhao 1982; Chow 1986; Armaleh and Desai 1987; API/ANSI 2011). The pile settlement can then be estimated by discretizing the pile, which is modeled as an elastic-compressible rod, into a number of elements separated by nodes. A “t-z” curve for each node and a “q-z” curve at the base are selected. A numerical integration approach of a two-point nonlinear boundary value problem is then applied to compute the settlement at the pile head due to a given load [e.g., Poulos and Davis (1980); Salgado (2008); Guo (2012); Poulos (2017)]. Relevant algorithms are easy to implement and available in different formats, including finite-differences (Ensoft TZPile, Reese et al. 2014), finite-elements (Oasys Pile, Oasys 2017), and alternative formulations (Scott 1981; Motta 1994; Seo et al. 2009; Crispin et al. 2018; Psaroudakis et al. 2019). The stress-displacement response at the pile tip (i.e., “q-z” curves) and group effects lie beyond the scope of this paper; therefore, the Illustrative Examples presented in the following work involve only end-bearing solitary piles.Some empirical “t-z” curves are available in the literature [e.g., Coyle and Reese (1966); Coyle and Sulaiman (1967); Reese et al. (1969); Vijayvergiya (1977); Frank and Zhao (1982); Abchir et al. (2016)] but are limited to specific test configurations. In addition, the American Petroleum Institute and American National Standards Institute (API/ANSI) (2011) provides general recommendations to develop idealized “t-z” curves for axially loaded piles based on a database of field tests. Analogous (“p-y”) curves are available for laterally loaded piles [e.g., McClelland and Focht (1956); Matlock (1970); Lam and Martin (1986); Reese and Van Impe (2011); Khalili-Tehrani et al. (2014)].Instead, theoretical “t-z” curves have been derived from the consideration of a (plane-strain) horizontal soil “slice.” Early solutions (Cooke 1974; Randolph and Wroth 1978; Baguelin and Frank 1979) assumed linear-elastic soil and employed the concentric cylinder model to describe the attenuation of shear stress with a radial distance from the pile. These idealizations have been used to establish linearly elastic springs for axially loaded piles which are widely used in practice [Scott (1981), Guo and Randolph (1998), Mylonakis and Gazetas (1998), Salgado (2008), Fleming et al. (2009), Viggiani et al. (2012), Guo (2012), Crispin et al. (2018)]. Following the same approach, nonlinear “t-z” curves were obtained by Vardanega et al. (2012), Kraft et al. (1981) and Chang and Zhu (1998) using power law, hyperbolic, and modified hyperbolic constitutive models, respectively. These curves have been employed in various studies by Chow (1986), Zhu and Chang (2002), Randolph (1994), Guo (2012), and Voyagaki et al. (2021).This method is versatile because it can incorporate variable soil properties with radial distance from the pile, installation effects, interface action between soil and pile, variable soil properties with depth, negative skin friction, and dynamic phenomena (Akiyoshi 1982; El Naggar and Novak 1994; Michaelides et al. 1998; Crispin and Mylonakis 2021). In addition, the method provides information about the stress and displacement fields around a pile that can be used for developing solutions for pile groups, such as those by Randolph and Wroth (1979) and Mylonakis and Gazetas (1998). Moreover, the formulation can accommodate both total and effective stresses (by using appropriate constitutive relations); hence, it can handle both fine- and coarse-grained soil materials under drained and undrained conditions. Nevertheless, despite the elegance and practical appeal of this method, only a handful of theoretical “t-z” curves are reported in the literature.This paper derives an extended set of “t-z” curves pertaining to different conditions based on the same horizontal soil slice model as the aforementioned theoretical solutions. A “toolbox” of “t-z” curves is provided to allow a designer to select an appropriate curve from a selected soil model and an easily fitted attenuation function. The unloading/reloading case is not investigated in this work (yet, no relevant restrictions are imposed by the derived “t-z” curves). Two illustrative calculation examples involving end-bearing piles in natural Pisa and remolded kaolinite are provided to demonstrate the simplicity of calibrating the “t-z” curves for application to a design problem. In these examples, the soil constitutive models and attenuation functions are fitted to available triaxial test data and numerical data, respectively. A preliminary version of some of the work presented in this paper is available in Bateman (2019).Horizontal Soil Slice ModelFig. 1 indicates the soil strip considered in the horizontal soil slice model, which is based on the following main assumptions. 1.The soil around the pile is divided into a series of independent horizontal “slices” of infinitesimal thickness, Fig. 1, subjected to shearing on horizontal and vertical planes.2.Since the continuity of the soil medium with depth is not considered, the slices provide resistance only to shearing, Fig. 1. Additional resistance due to pile-induced normal stresses acting on the upper and lower faces of the slices is neglected.3.Shear strains depend solely on vertical soil displacements (i.e., γ≅∂u/∂r). The effect of radial displacements on shear strains is small and can be neglected.4.The attenuation of shear stresses with radial distance from the pile, τ=τ(r), can be established (often taken as proportional to 1/r by equilibrium considerations, as discussed in the following).5.Deformations are infinitesimal.Note that the second assumption can be relaxed using an extended set of attenuation functions, τ(r), to be discussed in the next section. To obtain explicit solutions, two additional assumptions (6 and 7) are needed:6.Soil stress-strain response in shear is cast in flexibility form, that is, γ=γ(τ), where τ and γ are the shear stress and shear strain acting on a soil element, respectively.7.The composite function γ=γ(τ(r)) is integrable over the radial coordinate, r.Using these assumptions, shear strain can be integrated over the radial distance to provide the soil settlement at the pile-soil interface due to an applied shear stress, τ0, at any given depth (Fig. 1) (1) u0=∫d/2∞γ(r)dr=∫d/2∞γ[τ(r)]drwhere u0 = soil settlement at the pile circumference; r = radial distance from the pile centerline; and d = pile diameter. Since volumetric strains do not explicitly appear in Eq. (1), the analysis highlights how immediate settlement develops in nearly incompressible media, such as saturated clay. Regarding the upper integration limit in Eq. (1), this paper shows that it is not always feasible to extend to infinity under certain conditions.After dividing the total side friction per unit pile length (πdτ0) by Eq. (1), the nonlinear secant Winkler modulus, k(τ0), of the soil slice due to the imposed shear stress τ0 can be readily obtained as (2) k(τ0)≡τ0πdu0=τ0πd∫d/2∞γ[τ(r)]drThis stiffness parameter naturally carries units of stress and represents a nonlinear extension of the classical modulus of subgrade reaction for axially loaded piles. It is noted that no assumptions have been made, or restrictions imposed, in Eqs. (1) and (2) regarding the spatial distribution of soil properties, including the low strain shear modulus and strength, which might vary with the radial distance from the pile [e.g., observed in field tests by Kalinski et al. (2001) and O’Neill (2001)]. However, the explicit solutions obtained later in this paper assume radially homogenous soil. Corrections to the solutions provided to account for radial inhomogeneity are available [e.g., Kraft et al. (1981), Bateman and Crispin (2020)].Eq. (1) indicates that two functions are required to obtain a theoretical “t-z” curve from the horizontal soil slice model describing, first, the attenuation of shear stresses, τ(r), and, second, the soil constitutive relationship (in flexibility form), γ(τ). Previous authors have utilized the concentric cylinder model to obtain an attenuation function, leading to a 1/r dependence (discussed later). This paper goes on to introduce three novel attenuation functions inspired by available continuum solutions [Mylonakis (2001a)], resulting in a combination of generalized power law and exponential decay functions. These attenuation functions are employed in conjunction with eight simplified soil constitutive models from the literature, leading to 32 families of theoretical “t-z” curves, the vast majority (28) of which are, to the authors’ knowledge, novel.Attenuation FunctionsTo derive attenuation functions of shear stresses, the cylindrical pile depicted in Fig. 2 is considered. Taking the vertical equilibrium of an infinitesimal soil element in cylindrical coordinates and focusing on pile-induced stresses yields [Randolph and Wroth (1978), Mylonakis (2001b), Anoyatis and Mylonakis (2012)] (3) ∂(rτ)∂r+∂σ∂zr=0where σ and τ = vertical normal and shear stresses, respectively; and z = depth below ground level.A simplified solution to this equation can be derived by assuming that variations in the vertical stress with a depth due to pile loading are negligible; accordingly, setting ∂σ/∂z=0 and integrating over r yields the elementary solution [Cooke (1974), Randolph and Wroth (1978); Baguelin and Frank (1979)] (4) This result is commonly known as the concentric cylinder (or plane strain) model, in which the soil is treated as a series of concentric cylinders with respective shear forces and no resistance by means of vertical normal stresses. This model provides a simple attenuation function that has been extensively employed in the literature (Scott 1981; Fleming et al. 2009; Guo 2012; Viggiani et al. 2012).A difficulty associated with the use of this equation lies in the singular nature of the associated displacement field. Indeed, integrating shear stresses over the radial coordinate leads to a logarithmic solution for displacements that diverges with increasing r. To correct this problem, Randolph and Wroth (1978) suggested an empirical radius, rm, beyond which the vertical soil settlement can be assumed negligible; rm is usually taken as proportional to the pile length, L, and specific values have been suggested by Randolph and Wroth (1978), Guo and Randolph (1998), Fleming et al. (2009), and Guo (2012).A more rigorous yet simple approach for tackling this problem is possible by first casting the governing equation [Eq. (3)] in displacement form through the introduction of the approximate stress-displacement relations τ≅−Gs∂u/∂r and σ≅−Ms∂u/∂z [as, for instance, done by Nogami and Novak (1976), Mylonakis (2001b), and Anoyatis et al. (2019)] (5) ∂∂r[r∂u(r,z)∂r]+η2r∂2u(r,z)∂z2=0where u = vertical soil settlement; η2=Ms/Gs is a dimensionless compressibility coefficient; Gs = elastic shear modulus; Ms = compressibility modulus of the soil material; and η = function of the Poisson’s ratio of the soil, νs. Suitable values of η (and, thus, Ms) are discussed in Appendix I and in Mylonakis (2001a, b).An approximate solution to Eq. (5) can be obtained using a technique analogous to the one employed for spread footings by Vlasov and Leontiev (1966), previously applied to this problem by Vallabhan and Mustafa (1996), Lee and Xiao (1999), and Mylonakis (2000, 2001b) (similar approaches have been employed for laterally loaded piles, including Guo and Lee (2001), Basu and Salgado (2008), Shadlou and Bhattacharya (2014), and Bhattacharya (2019), among others). By separating the displacement function in Eq. (5), u(r,z), into radial, u(r), and vertical, χ(z), components [u(r,z)=u(r)χ(z)], multiplying by a virtual displacement χ(z) and integrating over the vertical coordinate, z, these authors obtained an alternative form of Eq. (5) that is now independent of depth [modified from Mylonakis (2000)] (6) d2u(r)dr2+(1r)du(r)dr−(2qd)2u(r)=0where q = compressibility constant that results from the integration with depth and is proportional to η. The specific form of q is chosen to enable it to be employed as a dimensionless fitting parameter to simplify the problem.Contrary to the elementary concentric cylinder model, Eq. (6) duly accounts for the continuity of the soil in the vertical direction. Eq. (6) is of the Bessel type and admits the following solutions for displacements and stresses (Mylonakis 2000; Olver et al. 2010) (7a) (7b) where K0 and K1 = modified Bessel functions of the second kind and order 0 and 1, respectively. Interestingly, these solutions are identical to those of the related dynamic plane strain problem for axially loaded piles pioneered by Baranov (1967) and Novak (1974), and tend to zero with increasing r (thus rendering empirical corrections such as rm unnecessary). A discussion of this remarkable similarity is provided in Mylonakis (2001a).Due to its complexity, this expression does not result in closed-form “t-z” curves as desired. Concentrating on the solution for stresses in Eq. (7b), this can be simplified considering the asymptotic expansions of the associated Bessel function in different regimes. For instance, for small radial distances from the pile centerline, the asymptotic form of the modified Bessel function is the power law function K1(x)∼(1/2)Γ(v)(2/x)1, where x = independent variable and Γ(v) = Gamma function [Olver et al. (2010)]. Substituting this expression into Eq. (7b) yields Eq. (4), which implies that the first term in Eq. (3) governs the behavior of the solution in the vicinity of the pile, as assumed in the classical concentric cylinder model, regardless of soil compressibility.For large radial distances from the pile centerline, the pertinent asymptotic expression is K1(x)∼(π/2x)1/2e−x [Olver et al. (2010)], which yields the approximate solution (referred in the ensuing as “power-exponential decay”) (8) τ(r)=τ0(2rd)−1/2exp[−q(2rd−1)]Eqs. (4) and (8) are based on sound theoretical arguments that include some necessary simplifying assumptions. Therefore, to provide flexibility in the ensuing analyses, these can be generalized in the form shown in Eqs. (9) and (10). By applying the additional fitting parameters (m and n), a better fit can be achieved (9) (10) τ(r)=τ0(2rd)−nexp[−q(2rd−1)]where exp(x) = exponential function. The first [Eq. (9)] of these modified solutions is referred in the ensuing as “generalized concentric cylinder,” whereas the second [Eq. (10)] is referred to as “generalized power-exponential decay” model. The last model encompasses all others by using appropriate values for parameters n and q. In addition, these equations have previously been suggested as approximate attenuation functions for displacements and associated time derivatives in relevant dynamic problems involving frequency and damping (Mylonakis 1995; Gazetas et al. 1998; Su et al. 2019); here, they are used as semiempirical static solutions for stresses.Eqs. (9) and (10) enable theoretical “t-z” curves to be analytically derived without reliance on the empirical radius rm, which is employed in the concentric cylinder model [provided m is greater than 1 in Eq. (9)]. Instead, these equations require fitting the parameter m [Eq. (9)] or the parameters q and n [Eq. (10)] to field test data, numerical continuum solutions or more rigorous analytical models, similar to how rm has been fitted by previous authors (Randolph and Wroth 1978; Guo 2012). An example of fitting these parameters to a more rigorous analytical solution (that is unsuitable for direct use in developing “t-z” curves and discussed in Appendix I) is detailed subsequently in the paper. The approximate range of the parameters from these examples, given in Table 2, are 0.15d), which agrees with the observations from the asymptotic solutions in Eqs. (4) and (8). These results suggest that the attenuation of shear stresses in the vicinity of the pile is governed by equilibrium, not stress-strain behavior. Soil compressibility starts becoming important at higher distances, where shear stresses and strains are generally lower than the corresponding normal stresses and volumetric strains. These traits justify using some elasticity arguments [used to derive Eqs. (8) and (10)] for tackling the stress distribution in the nonlinear problem at hand. Similar patterns are observed in other relevant problems, such as the Boussinesq and Mindlin solutions [e.g., Davis and Selvadurai (1996)].Soil Constitutive ModelsEight soil constitutive models are considered in this paper. To obtain an explicit solution using Eq. (1), the selected constitutive relations must be cast in flexibility form, that is, γ=γ(τ). Each soil constitutive model (illustrated in Fig. 4) is given as follows in both the (a) stiffness and (b) flexibility forms. A table summarizing the models is provided in Appendix III. Each model contains a number of parameters that should be fitted to site-specific laboratory data and/or in situ testing. An example of fitting these constitutive models to the soil test data is subsequently discussed.The simplest constitutive relation considered is the linear stress-strain model defined by a single shear modulus value, G, given in stiffness and flexibility form by (for τ<τmax) (11a) (11b) The bilinear model provides a simple adaption by splitting the linear model into two regions, the first having a stiffness G1 and the second having a reduced stiffness G2, which is applicable after a predefined yield stress, τ1, is reached. Additional stiffness changes could be introduced to more closely model nonlinear soil behavior. However, doing so would be at the expense of an increasingly complex function involving many fitting parameters. The bilinear model is given by (for τ<τmax) (12a) τ={γG1γ≤τ1/G1G2(γ−τ1G1)+τ1γ≥τ1/G1(12b) γ={τG1τ≤τ1τ−τ1G2+τ1G1τ≥τ1Alternatively, a simple power law relationship can be employed, in this case defined with the shear strain when 50% of the soil shear strength is mobilized, γ50, and a soil nonlinearity exponent, b. This model was employed by Vardanega et al. (2012) to generate simple “t-z” curves [see also Williamson (2014), Vardanega (2015), Crispin et al. (2019)]. The power law relationship has the advantage that rm is not required due to an infinite initial stiffness (b<1). The power law model is given by (for τ<τmax) (13a) (13b) However, since infinite initial stiffness is not physically realizable, a linear model can be employed at low stress (τ≤τi), with an initial stiffness Gi, resulting in the linear-power law soil constitutive model (for τ<τmax) (14a) τ={γGiγ≤τi/Giτmax2(γγ50)bγ≥τi/Gi(14b) γ={τGiτ≤τiγ50(2ττmax)1bτ≥τiwhere (14c) τi=(τmax2)(2Giγ50τmax)bb−1A combination of the linear and power law models was suggested by Ramberg and Osgood (1943), which is defined here with a reference shear strain, γr, and two fitting constants, c1 and c2 (Ishihara 1996). This model is available in only the flexibility form (for τ<τmax) (15) γ=γr[ττmax+(c1ττmax)c2]The five constitutive relationships previously discussed do not approach a limiting stress; therefore, a shear strength cap, τmax, is introduced. For undrained conditions, this ultimate strength should be taken as equal to the undrained shear strength, τmax=cu, which, in turn, depends on the shear mode of the soil test [e.g., Mayne (1985), Beesley and Vardanega (2020)]. Nevertheless, a generic parameter τmax is employed because the “t-z” curves are not linked to a particular test or limited to undrained conditions. In practice, the shear strength of the resultant “t-z” curves might be capped at a lower value to represent slip at the pile-soil interface. To this end, empirical factors such as the familiar adhesion parameter α (Tomlinson 1957; Skempton 1959; Meyerhof 1976; Chakraborty et al. 2013) are available for undrained problems, which can be used to limit the maximum shear stress. Note that in the interest of space, simple names such as “Bilinear model” are adopted here over more accurate alternatives such as “capped Bilinear model” or “Bilinear-fully plastic model.”Instead, the hyperbolic model, which approaches a single value at large strain, can be employed. Different forms of this model are available in the literature [e.g., Kondner (1963)], and it was utilized by Kraft et al. (1981) to develop a “t-z” curve using the concentric cylinder model. The form of the hyperbolic model adopted by Kraft et al. (1981) is employed with a fitting parameter, Rf, which acts as a factor to τmax that alters the location of the asymptote. This enables a better fit of this model over a desired stress range, although this results in the model performing poorly at high stresses. In the classical model, Rf is expected to be less than 1 and requires a shear strength cap similar to the previous five models. However, in this paper, it has been employed as a fitting parameter and allowed to vary outside this range. The hyperbolic model is given by (16a) τ=τmaxRf[τmaxRfGiγ+1]−1(16b) γ=τmaxRfGi[τmaxRfτ−1]−1In addition, a modified hyperbolic function is employed with an exponent, c3, to give more control over the stiffness degradation, as discussed by Fahey and Carter (1993). This model has previously been utilized to develop a “t-z” curve in conjunction with the concentric cylinder model by Chang and Zhu (1998). The modified hyperbolic model is given only in the flexibility form (17) γ=τGi[1−(Rfττmax)c3]−1A similar asymptotic property can be obtained using an exponential relationship resulting in the exponential model indicated with a similar fitting parameter, Rf, employed as before (18a) τ=τmaxRf{1−exp[−γ(RfGiτmax)]}(18b) γ=−τmaxRfGiln(1−Rfττmax)Note that some of the formulations include an initial (low-strain) shear modulus, Gi, which is employed as a model parameter to enable a better fit of the constitutive relation to the available data over the strain range of interest. Alternatively, when small strains are important, this modulus might be set equal to the maximum shear modulus, Gmax, from high-quality experimental testing (e.g., seismic cone penetration testing).“t-z” CurvesBy substituting the four attenuation functions [Eqs. (4), (8), (9), and (10)] into each constitutive model [Eqs. (11)–(18)], theoretical “t-z” curves can be derived analytically using Eq. (1). As an example, substituting the concentric cylinder model [Eq. (4)] into the linear soil constitutive model Eq. (11b) results in a function, γ(r). By inputting this directly into Eq. (1), infinite settlement would be predicted when integrating over an infinite distance; thus, the upper integration limit is replaced with rm, and the following equation is derived. This is valid until slip occurs at the pile-soil interface (Randolph and Wroth 1978; Fleming et al. 2009) (19) u0d=τ02G∫d/2rm1rdr=τ02Gln(2rmd)Theoretical “t-z” curves derived using the concentric cylinder model Eq. (4) are given as follows and are summarized in Appendix IV.Bilinear [Eq. (12)] (20) u0d={τ02G1ln(2rmd)τ0≤τ1τ0(G1−G2)2G1G2[G2G1−G2ln(2rmdτ1τ0)+G1G1−G2ln(τ0τ1)+(τ1τ0−1)]τ0≥τ1Power Law, for b<1 (Vardanega et al. 2012; Williamson 2014; Vardanega 2015; Crispin et al. 2019) [Eq. (13)] (21) u0d=γ50b2(1−b)(2τ0τmax)1bLinear-Power Law [Eq. (14)] (22a) u0d={τ02Giln(2rmd)τ0≤τiγ50b2(1−b)(2τ0τmax)1b[1−(2rid)b−1b]+τ02Giln(rmri)τ0≥τiwhere (22b) 2rid=2τ0τmax(2γ50Giτmax)b1−bRamberg-Osgood [Eq. (15)] (23) u0d=τ0γr2τmaxln(2rmd)+γr2(c2−1)(c1τ0τmax)c2[1−(2rmd)1−c2]Hyperbolic [Kraft et al. (1981)] [Eq. (16)] (24) u0d=τ02Gi[ln(2rmd−Rfτ0τmax)−ln(1−Rfτ0τmax)]Modified Hyperbolic [Chang and Zhu (1998)] [Eq. (17)] (25) u0d=τ02Gic3[ln((2rmd)c3−(Rfτ0τmax)c3)−ln(1−(Rfτ0τmax)c3)]Exponential [Eq. (18)] (26) u0d=τmax2RfGi[ln(1−Rfτ0τmax)−2rmdln (1−Rfτ0τmax(d2rm))+Rfτ0τmax[ln(2rmd−Rfτ0τmax)−ln(1−Rfτ0τmax)]]Theoretical “t-z” curves derived using the generalized concentric cylinder model [Eq. (9)] are given in Eqs. (27)–(34) and are summarized in Appendix IV. The rm value has been included in these solutions; however, when m>1, the relevant terms vanish and the solutions can be further simplified.Linear (27) u0d=τ02G(m−1)[1−(2rmd)1−m]Bilinear (28) u0d={τ02G1(m−1)[1−(2rmd)1−m]τ0≤τ1τ0(G1−G2)2G1G2{1m−1[G1G1−G2−G2G1−G2(2rmd)1−m−(τ0τ1)1−mm]+[τ1τ0−(τ0τ1)1−mm]}τ0≥τ1Power Law, for b

Source link