ANALYSIS OF BRIDGE-PILE FOUNDATION SYSTEM IN MULTI-LAYERED NON-LINEAR SOIL STRATA USING ENERGY-BASED METHOD

Abstract

Disclosed is a method for analyzing displacement responses of piles subjected to at least one loading condition selected from an axial load, a lateral load, and a bending moment. Also disclosed is a non-transitory machine-readable storage medium including instructions executable by a processing resource of a computing device to cause the processing resource to perform the inventive method.

Description

BACKGROUND

1. Field of the Invention

This invention relates to structural/geotechnical engineering.

2. Description of Related Art

Structures such as high-rise buildings, bridges, offshore platforms, etc., which rest on piles embedded in multi-soil strata (see, e.g., FIG. 1), are acted upon by various horizonal forces (e.g., waves, currents and seismic events), vertical forces and moments. Although extensive research has been carried out in the area of soil pile interaction, analysis and design of a pile embedded in multi-soil strata still poses a challenge.

Numerous studies have been reported in published literature on analytical solutions for a pile foundation system installed in a homogeneous single soil layer. However, piles are rarely installed in an ideal homogeneous single soil layer.

Studies of pile foundation systems installed in multilayered soil have also been reported. See, e.g., Basu et al. “Analysis of laterally loaded piles in multilayered soil deposits.” Joint Transportation Research Program (2008): 330, and Hashem-Ali. “Analytical methods for predicting load-displacement behaviour of piles.” PhD diss., Durham University, 2014.

Despite the foregoing developments, it is desired to provide a method and a device for further advancing the existing analytical techniques and contribute to better understanding of the behavior of deep foundation systems in non-linear soil strata.

All references cited herein are incorporated herein by reference in their entireties.

SUMMARY

A first aspect comprises a method for analyzing displacement responses of piles subjected to at least one loading condition selected from an axial load, a lateral load, and a bending moment, includes the following steps: (a) inputting values r₀, L, λ, G₀, E_p, P_t, and ν; (b) calculating γ1, old and γ2, old; (c) determining integrations of soil parameters C, k and m; (d) calculating soil dimensionless function ϕ; (e) calculating pile displacement ν; (f) calculating γ1, new and γ2, new; (g) determining whether a condition (γ_old−γ_new)/γ_old<0.001 is true; and (h) repeating steps (a)-(g) until the condition is true.

A second aspect comprises a non-transitory machine-readable storage medium including instructions executable by a processing resource of a computing device to cause the processing resource to perform the above inventive method.

Another aspect comprises a method for analyzing displacement responses of a pile subjected to a lateral loading condition. One embodiment of the method comprises inputting values r₀, L, λ, G₀, E_p, and I_p, and applying force and moments, Q₀ and M₀; choosing initial values of y₁, y₂, y₃, y₄, y₅, y₆, y₇ and y₈, and calculating ϕ_r and ϕ_θ; calculating deviatoric strain ε_q, and using the deviatoric strain to calculate soil stiffness, G, where

$G={G_0\left(\frac{{\epsilon }_q}{{\epsilon }_{\mathrm{qo}}}\right)}^n;$

calculating n_s1, n_s2, m_s1, m_s2 and m_s3; calculating pile displacement u; calculating calculate y₁^new y₂^new y₃^new y₄^new y₅^new y₆^new y₇^new and y₈^new; determining whether a condition y_old−y_new)/y_old<0.001 is true; and repeating the above steps until the condition is true.

And yet another aspect comprises a non-transitory machine-readable storage medium including instructions executable by a processing resource of a computing device to cause the processing resource to perform the above method for analyzing displacement responses of a pile subjected to a lateral loading condition.

Other systems, methods, features and/or advantages will be or may become apparent to one with skill in the art upon examination of the following drawings and detailed description. It is intended that all such additional systems, methods, features and/or advantages be included within this description and be protected by the accompanying claims.

BRIEF DESCRIPTION OF SEVERAL VIEWS OF THE DRAWINGS

The invention will be described in conjunction with the following drawings in which like reference numerals designate like elements and wherein:

FIG. 1 is a schematic view of a high-rise building and its foundation resting on multilayered soil strata.

FIG. 2 is a screenshot of an exemplary graphical user interface of an embodiment of software according to the invention.

FIG. 3 is a graph of depth against pile deflection calculated by the method of the invention and by several comparative methods.

FIG. 4 illustrates an axially loaded pile in an isotropic non-linear elastic medium.

FIGS. 5A and 5B illustrate stresses and displacements within a soil continuum where FIG. 5A illustrates stresses within a soil continuum and FIG. 5B illustrates displacements within a soil continuum.

FIGS. 6A and 6B illustrate discretization of the soil where FIG. 6A illustrates a top view and FIG. 6B illustrates a front view.

FIG. 7 illustrates the variation of normalizing shear secant with logarithmic strain ε_q or normalized displacement.

FIG. 8 illustrates degradation of tangent with deviatoric strain.

FIG. 9 illustrates logarithmic scale of degradation of tangent stiffness with strain level.

FIG. 10 is an exemplary flowchart depicting the iterative solution procedure for calculating pile head displacements when axially (vertically) loaded.

FIG. 11 illustrates a laterally loaded pile in layered non-linear elastic medium.

FIG. 12 is an exemplary flowchart depicting the iterative solution procedure for calculating pile displacements when laterally loaded.

FIGS. 13, 14 and 15 are exemplary annotated screenshots of graphical user interfaces of embodiments of software according to the invention.

DETAILED DESCRIPTION

Before the present methods and systems are disclosed and described, it is to be understood that the methods and systems are not limited to specific synthetic methods, specific components, or to particular compositions. It is also to be understood that the terminology used herein is for the purpose of describing particular embodiments only and is not intended to be limiting.

As used in the specification and the appended claims, the singular forms “a,” “an” and “the” include plural referents unless the context clearly dictates otherwise. Ranges may be expressed herein as from “about” one particular value, and/or to “about” another particular value. When such a range is expressed, another embodiment includes—from the one particular value and/or to the other particular value. Similarly, when values are expressed as approximations, by use of the antecedent “about,” it will be understood that the particular value forms another embodiment. It will be further understood that the endpoints of each of the ranges are significant both in relation to the other endpoint, and independently of the other endpoint.

“Optional” or “optionally” means that the subsequently described event or circumstance may or may not occur, and that the description includes instances where said event or circumstance occurs and instances where it does not.

Throughout the description and claims of this specification, the word “comprise” and variations of the word, such as “comprising” and “comprises,” means “including but not limited to,” and is not intended to exclude, for example, other additives, components, integers or steps. “Exemplary” means “an example of” and is not intended to convey an indication of a preferred or ideal embodiment. “Such as” is not used in a restrictive sense, but for explanatory purposes.

Disclosed are components that can be used to perform the disclosed methods and systems. These and other components are disclosed herein, and it is understood that when combinations, subsets, interactions, groups, etc. of these components are disclosed that while specific reference of each various individual and collective combinations and permutation of these may not be explicitly disclosed, each is specifically contemplated and described herein, for all methods and systems. This applies to all aspects of this application including, but not limited to, steps in disclosed methods. Thus, if there are a variety of additional steps that can be performed it is understood that each of these additional steps can be performed with any specific embodiment or combination of embodiments of the disclosed methods.

The present methods and systems may be understood more readily by reference to the following detailed description of embodiments and to the Figures and their previous and following description.

As will be appreciated by one skilled in the art, the methods and systems may take the form of an entirely hardware embodiment, an entirely software embodiment, or an embodiment combining software and hardware aspects. Furthermore, the methods and systems may take the form of a computer program product on a computer-readable storage medium having computer-readable program instructions (e.g., computer software) embodied in the storage medium. More particularly, the present methods and systems may take the form of web-implemented computer software. Any suitable computer-readable storage medium may be utilized including hard disks, CD-ROMs, optical storage devices, or magnetic storage devices.

Embodiments of the methods and systems are described below with reference to block diagrams and flowchart illustrations of methods, systems, apparatuses and computer program products. It will be understood that each block of the block diagrams and flowchart illustrations, and combinations of blocks in the block diagrams and flowchart illustrations, respectively, can be implemented by computer program instructions. These computer program instructions may be loaded onto a general-purpose computer, special purpose computer, or other programmable data processing apparatus to produce a machine, such that the instructions which execute on the computer or other programmable data processing apparatus create a means for implementing the functions specified in the flowchart block or blocks.

These computer program instructions may also be stored in a computer-readable memory that can direct a computer or other programmable data processing apparatus to function in a particular manner, such that the instructions stored in the computer-readable memory produce an article of manufacture including computer-readable instructions for implementing the function specified in the flowchart block or blocks. The computer program instructions may also be loaded onto a computer or other programmable data processing apparatus to cause a series of operational steps to be performed on the computer or other programmable apparatus to produce a computer-implemented process such that the instructions that execute on the computer or other programmable apparatus provide steps for implementing the functions specified in the flowchart block or blocks.

Accordingly, blocks of the block diagrams and flowchart illustrations support combinations of means for performing the specified functions, combinations of steps for performing the specified functions and program instruction means for performing the specified functions. It will also be understood that each block of the block diagrams and flowchart illustrations, and combinations of blocks in the block diagrams and flowchart illustrations, can be implemented by special purpose hardware-based computer systems, including for example cloud-based systems, that perform the specified functions or steps, or combinations of special purpose hardware and computer instructions.

The increasing demand for adopting pile foundations in bridges has pointed towards the need to constantly improve the existing analytical techniques for better understanding of the behavior of such foundation systems. Disclosed and described herein are embodiments of an energy-based approach (simplified continuum model) to assess the displacement responses to general loading conditions (e.g., axial (e.g., vertical) load, lateral load and bending moment) of piles embedded in multilayered soil strata considering the soil non-linear behavior. The governing differential equations and the boundary conditions for a bridge pile embedded in multi-layered soil strata subjected to the general loading conditions are obtained using the Hamilton's principle employing variational principles and minimization of energies.

Soil discretization (non-linearity) has been incorporated through simple constitutive relationships that account for degradation of soil moduli with increasing strain values. The soil moduli (λ, and G) are assumed to vary in radial, circumference and depth directions according to the strain and stress levels. A simple power law, which was developed where the soil is assumed to be nonlinear elastic and perfectly plastic, is used. A Tresca yield surface is assumed to develop the soil stiffness variation with different strain levels that defines the non-linearity of the soil strata. This numerical technique can be applied to a pile foundation in multi-layered soil strata for a pier supporting the bridge and solved by customizing the software using, for example, MATLAB™ (The MathWorks, Inc., Natick, MA) R2019a (see FIG. 2). The inventive method is preferably based on the Finite Difference Method since it yields results faster when compared to the Finite Element Method.

The analysis yields the bridge pile displacements at any depth along the length of the pile. The results from the energy-based method are compared with those from the field test data as well as the Finite Element Analysis based on, for example, the software ANSYS Workbench 2021 R1 (ANSYS, Inc., Canonsburg, PA). As shown in FIG. 3, the results of the analysis of embodiments described herein are in good agreement with the published field data of McClelland and Focht (McClelland et al. “Soil modulus for laterally loaded piles.” Transactions of the American Society of Civil Engineers 123, no. 1 (1958): 1049-1063) and the three-dimensional finite element analysis results performed using the software ANSYS 2019R3. The results obtained are in good agreement with the field data when compared to the linear elastic solution that does not consider the soil non-linearity.

The methodology can be extended to study the response of the multi-strata soil supporting group piles underneath bridge piers.

The inventive method includes obtaining load-displacement responses on a single pile embedded in multi-layered non-linear soil strata subjected to the general loading conditions (axial load, lateral load and a moment). The method has been developed in the following three steps: (i) a single pile embedded in multi-layered non-linear soil strata is assessed for displacements when subjected to axial loads; (ii) a single pile embedded in multi-layered non-linear soil strata is assessed for displacements when subjected to lateral loads and moments; and (iii) a single pile embedded in multi-layered non-linear soil strata is assessed for displacements when subjected to axial loads, lateral loads and moments (generally loaded pile). This is an analytical solution to assess the load-displacement responses given any field data and is effective in capturing the responses quite well when compared to the existing field test data and the finite element analyses. By further validation, it was proven to be sufficient in analyzing the responses without undue experimentation. Software programs have been developed which embody certain aspects of the invention. GEOS-ALPILE has been developed for axial loads, GEOS-LLPILE has been developed for lateral loads and moments and GEOS-GLPILE has been developed for general loads. The developed software works on an iterative procedure and are very user-friendly as inputs of basic pile and soil properties provide outputs within seconds. It can be used by practitioners, researchers and for academia purposes (e.g., by students to understand the responses in an easy way).

I. Energy-Based Approach: Analysis of a Vertically Loaded Pile in Multi-Layered Non-Linear Soil Strata

Referring to Aryan, Prakash & Arockiasamy, Madasamy. (March, 2021). Energy-Based Approach for the Analysis of a Vertically Loaded Pile in Multi-layered Non-linear Soil Strata. 10.21203/rs.3.rs-337774/v1, which is fully incorporated by reference and made a part hereof, an analytical model has been developed based on the energy-based approach to predict the load-displacement responses of the piles embedded in multi-layered soil strata subjected to static axial loads.

Disclosed is an energy-based approach (simplified continuum model) to provide a displacement profile of an axially loaded pile in multi-layered soil strata considering the soil's non-linear elastic behavior. A mathematical framework is disclosed that considers the soil as a 3D elastic continuum and the pile is modeled following the elastic Euler-Bernoulli beam theory. The differential equations governing the soil and pile displacements have been obtained utilizing the principle of minimum potential energy assuming rational soil displacement functions. An iterative algorithm is adopted that solves the differential equations analytically and numerically. The study investigates the effect of explicit incorporation of soil properties and layering in order to understand the importance of predicting appropriate pile displacement responses in linear elastic soil system. The responses indicated that the effect of soil layers and their thicknesses, pile properties and the variation in soil moduli have a direct impact on the displacements of piles subjected to axial loading. Hence, proper emphasis has to be given to account for the soil non-linearity.

An axially loaded pile in an isotropic non-linear elastic multi-layered soil medium is shown in FIG. 4. This disclosure considers a pile of length L with circular cross section of radius r₀. The pile is embedded in n horizontal soil layers and is subjected to an axial (vertical) load P_t. The horizontal soil layers extend to infinity in the radial direction and the bottom n^th layer extends to infinity in the vertical direction.

In FIG. 4, the terms H₁, H₂, H₃ . . . H_n−1 denote the vertical height from the ground surface to the bottom of any layer i. Therefore, the thickness of any layer i is H_i−H_i−1 with H₀=0. Due to the axisymmetric problem behavior, a system of cylindrical coordinates (r-θ-z) is chosen with the origin coinciding with the center of the pile cross section at the pile head and the z axis coinciding with the pile axis. The pile head is considered to be free, and the tip of the pile is clamped. Another important assumption to be noted is that there is no slippage or separation between the pile and the surrounding soil and between soil layers. The stresses and the displacement within a soil continuum are shown below in FIGS. 5A and 5B.

At each point in the soil domain, soil elastic parameters G and λ have been calculated (see FIGS. 6A and 6B). The disclosed methods obtain pile head deflection caused by the action of the vertical load at the pile head.

For an axially loaded pile, the horizontal and tangential displacements can be neglected as they are accompanied by very small strains. For the case of a pile with circular cross-section, there are two functions to be considered: ν(z), which will represent the vertical displacement at depth z and the dimensionless functions ϕ(r), describing the variation of soil displacements in the radial direction.

The vertical displacement at any point of the soil is represented as a function in (r, z):

ν_r=0

ν_θ=0

ν_z(r, z)=ϕ(r)ν(z)

For a given uniform cross-sectional area of the pile along the length, ϕ(r)=1 when r=0 to r₀, whereas ϕ(r)=0 when r→∞. This explains the decay of the function ϕ(r) with an increase in the radial direction.

The pile and its surrounding elastic medium are subjected to a vertical displacement of the pile soil system when it is acted upon by a vertical load. The total potential energy of the pile and the soil is a summation of internal potential energy and external potential energy, which is given as:

$\begin{array}{cc}\Pi =\frac{1}{2}{E}_p{A}_p{\int }_0^L{\left(\varphi \frac{dv}{dz}\right)}^2+\frac{1}{2}{\int }_0^L{\int }_0^{2\pi }{\int }_{r_0}^{\infty }{\sigma }_{ij}{\epsilon }_{ij}rdrd\theta +\frac{1}{2}{\int }_L^{\infty }{\int }_0^{2\pi }{\int }_0^{\infty }{\sigma }_{ij}{\epsilon }_{ij}\mathrm{rdr}-{\mathrm{vP}}_{t\mathrm{at}z=0}& \left(1\right)\end{array}$

where E_p denotes the elastic Young modulus of the pile, A_p denotes the cross-section of the pile, ν represents the vertical displacement, P_t is the vertical load and σ_ij, ε_ij are stress and strain components, respectively. The first term of the equation represents potential pile energy, the second and third terms are potential energy from the surrounding soil and the soil below the pile, respectively.

The soil non-linearity is considered by varying the soil elastic parameters (G and λ) at each discretized nodal point in the soil domain. The stress-strain and strain-displacement relationships at any given nodal point in the soil medium are idealized by the following relationships. The stress-strain relationship is expressed as:

$\begin{array}{cc}\left[\begin{array}{c}{\sigma }_r\\ {\sigma }_{\theta }\\ {\sigma }_z\\ {\tau }_{r\theta }\\ {\tau }_{\mathrm{rz}}\\ {\tau }_{\theta z}\end{array}\right]=\left[\begin{array}{cccccc}\lambda +2G& \lambda & \lambda & 0& 0& 0\\ \lambda & \lambda +2G& \lambda & 0& 0& 0\\ \lambda & \lambda & \lambda +2G& 0& 0& 0\\ 0& 0& 0& G& 0& 0\\ 0& 0& 0& 0& G& 0\\ 0& 0& 0& 0& 0& G\end{array}\right]\left[\begin{array}{c}{\epsilon }_r\\ {\epsilon }_{\theta }\\ {\epsilon }_z\\ {\gamma }_{r\theta }\\ {\gamma }_{\mathrm{rz}}\\ {\gamma }_{\theta z}\end{array}\right]& \left(2\right)\end{array}$

where G and λ are the elastic constants of the soil. The strain-displacement relationship is given by:

$\begin{array}{cc}\left[\begin{array}{c}{\epsilon }_r\\ {\epsilon }_{\theta }\\ {\epsilon }_z\\ {\gamma }_{r\theta }\\ {\gamma }_{\mathrm{rz}}\\ {\gamma }_{\theta z}\end{array}\right]=\left[\begin{array}{c}-\frac{\partial u_r}{\partial r}\\ -\frac{u_r}{r}-\frac{1}{r}\frac{\partial u_{\theta }}{\partial \theta }\\ -\frac{\partial u_z}{\partial z}\\ -\frac{1}{r}\frac{\partial u_r}{\partial \theta }-\frac{\partial u_{\theta }}{\partial r}+\frac{u_{\theta }}{r}\\ -\frac{\partial u_z}{\partial r}-\frac{\partial u_r}{\partial z}\\ -\frac{1}{r}\frac{\partial u_z}{\partial \theta }-\frac{\partial u_{\theta }}{\partial z}\end{array}\right]=\left[\begin{array}{c}0\\ 0\\ -\varphi \left(r\right)\frac{\mathrm{dv}\left(z\right)}{\mathrm{dz}}\\ 0\\ -v\left(z\right)\frac{d\varphi \left(r\right)}{\mathrm{dr}}\\ 0\end{array}\right]& \left(3\right)\end{array}$

By substituting Equation (3) into Equation (2), the strain energy density function,

$W=\frac{{\sigma }_{\mathrm{ij}}{\epsilon }_{\mathrm{ij}}}{2}$

is obtained, where the summation implies the repetition of the indices i and j as required in indicial notation:

$\begin{array}{cc}\frac{1}{2}{\sigma }_{ij}{\epsilon }_{ij}=\frac{1}{2}\left[\left(\lambda +2G\right){\left(\varphi \frac{dv}{dz}\right)}^2+{G\left(v\frac{d\varphi }{\mathrm{dr}}\right)}^2\right]& \left(4\right)\end{array}$

Substituting Equation (4) into Equation (1) and integrating with respect to θ the potential energy equation becomes:

$\begin{array}{cc}\Pi =\frac{1}{2}{E}_p{A}_p{\int }_0^L{\left(\varphi \frac{dv}{dz}\right)}^2+\mathrm{dz}+\pi {\int }_0^L{\int }_{r_0}^{\infty }\left(\left(\lambda +2G\right){\left(\varphi \frac{dv}{dz}\right)}^2+{G\left(v\frac{d\varphi }{\mathrm{dr}}\right)}^2\right)\mathrm{rdrdz}+\pi {\int }_L^{\infty }{\int }_0^{\infty }\left(\left(\lambda +2G\right){\left(\varphi \frac{dv}{dz}\right)}^2+{G\left(v\frac{d\varphi }{\mathrm{dr}}\right)}^2\right)\mathrm{rdrdz}-{\mathrm{vP}}_{t\mathrm{at}z=0}& \left(5\right)\end{array}$

The variational principle has been used to calculate the potential energy δU and the external energy δW. As a result, the governing equations of the pile-soil system are obtained by minimizing the potential energy of soil and pile. The expression of potential energy contains different functions, such as ν(z), ϕ(r), dν(z)/dz and dϕ(r)/dr, and so minimizing the potential energy gives:

$\begin{array}{cc}\mathrm{\delta \Pi }=\left[A\left(v\right)\delta v+B\left(v\right)\delta \left(\frac{dv}{dz}\right)\right]+\left[C\left(\varphi \right)\mathrm{\delta \varphi }\right]& \left(v\right)\end{array}$

where A, B and C are the terms associated with variations δν, δ(dν/dz) and δϕ.

The variation of Equation (5) becomes

$\begin{array}{cc}\delta \prod =\frac{1}{2}{E}_p{A}_p{\int }_0^L\varphi \frac{dv}{dz}\delta \varphi \frac{dv}{dz}dz+\pi {\int }_0^L{\int }_{r_0}^{\infty }\left(\left(\lambda +2G\right)\left(\varphi \frac{dv}{dz}\delta \varphi \frac{dv}{dz}\right)+G\left(v\frac{d\varphi }{\mathrm{dr}}\delta v\frac{d\varphi }{\mathrm{dr}}\right)\right)\mathrm{rdrdz}+\pi {\int }_L^{\infty }{\int }_0^{\infty }\left(\left(\lambda +2G\right)\left(\varphi \frac{dv}{dz}\delta \varphi \frac{dv}{dz}\right)+G\left(v\frac{d\varphi }{\mathrm{dr}}\delta v\frac{d\varphi }{\mathrm{dr}}\right)\right)\mathrm{rdrdz}-{\mathrm{vP}}_{t\mathrm{at}z=0}& \left(6\right)\end{array}$

The governing equation of the pile is obtained for 0<z<L by collecting terms associated with δνdz, δ(dν/dz)dz and its derivative, δν and δ(dν/dz)dz≠0. The governing equation is obtained as follows:

$\begin{array}{cc}k\frac{d^{2}v}{d{z}^2}+C\frac{dv}{dz}+mv=0& \left(7\right)\end{array}$ $\mathrm{where}$ $\begin{array}{cc}C=2\pi {\int }_{r_0}^{\infty }r\left[\left(\lambda +2G\right){\varphi }^2\right]\mathrm{dr}& \left(8\right)\end{array}$ $\begin{array}{cc}k=E_p{I}_p+2\pi {\int }_{r_0}^{\infty }r\left[\left(\lambda +2G\right){\varphi }^2\right]\mathrm{dr}& \left(9\right)\end{array}$ $\begin{array}{cc}m=-2\pi {\int }_{r_0}^{\infty }r{G\left(\frac{d\varphi }{\mathrm{dr}}\right)}^2\mathrm{dr}& \left(10\right)\end{array}$

The tip of a pile is assumed to be clamped, which means that the displacement and the curvature are equal to zero at the base of the pile. The boundary conditions are obtained by collecting δν and δ(dν/dz). At the head of the pile (z=0):

$\begin{array}{cc}-E_p{A}_p-2\pi {\int }_{r_0}^{\infty }\left[\lambda +G\left(2+{\varphi }^2\right)\right]dr\frac{dv}{dz}+P_t=0& \left(11\right)\end{array}$

The displacement at the tip of the pile (z=L):

ν=0   (12)

The second order differential Equation (7) can be solved using a central finite difference scheme. Equation (7) becomes

$\begin{array}{cc}k\left(\frac{v_{i-1}-2v_i+v_{i+1}}{\Delta z^2}\right)+C\left(\frac{v_{i-1}+v_{i+1}}{\Delta z}\right)+m{v}_i=0& \left(13\right)\end{array}$

where i denotes the i^th node in z direction, and Δz is the distance between two nodes. This discretized analysis is then solved using software such as MATLAB R2019a or specifically-programmed software.

The governing equation of the soil is obtained for r₀≤r≤∞ by collecting terms associated with δϕdr:

$\begin{array}{cc}r\frac{d^2\varphi }{{\mathrm{dr}}^2}+{\gamma }_1\frac{d\varphi }{\mathrm{dr}}+{\gamma }_2\varphi =0& \left(14\right)\end{array}$ $\mathrm{where}$ $\begin{array}{cc}{\gamma }_1=\frac{{\int }_0^{L}G+\mathrm{rG}}{{\int }_0^{L}G}& \left(15\right)\end{array}$ $\begin{array}{cc}{\gamma }_2=\frac{-{\int }_0^L\left(\lambda +\mathrm{rG}\right){\left(\frac{\mathrm{dv}}{dz}\right)}^2}{{\int }_0^{L}G{v}^2}& \left(16\right)\end{array}$

Note that the dimensionless parameter γ, defined in the above equation, describes the function ϕ. Similar to the solution of Equation (7) by the central finite difference scheme, the governing differential Equation (14) is also solved using the FDM (using, for example, MATLAB R2019a software).

The variation of the shear stress with strain can be described using two parameters A and n that have been obtained experimentally using a pressuremeter test as shown in the equation below:

q=A(ε_q)ⁿ   (17)

where q represents the equivalent shear stress and ε_q is the deviator shear strain. It has been shown that the decay curves of the soil stiffness with strain can be divided into three regions as shown in FIG. 7. The first region in FIG. 7 represents the very small strain where the stiffness, G₀ is constant, the second region comprising small strains starts from ε₀ until ε=0.1% and the third region exceeds ε=0.1%, which is indicative of large strains.

In the second region, the stiffness decays rapidly and in the third region with large strain levels, the stiffness is the smallest, which concludes that the soil stiffness is high at the small strain and decreases with the large strain. FIG. 8 shows the degradation of soil stiffness with increasing strains for different clay types.

This disclosure assumes decay of soil stiffness with strain using a power law to describe the stress-strain behavior of soil:

$\begin{array}{cc}G={G_0\left(\frac{{\epsilon }_q}{{\epsilon }_{q0}}\right)}^n& \left(18\right)\end{array}$ $\begin{array}{cc}G=a{\epsilon }_q^n& \left(19\right)\end{array}$

where

$a=\frac{G_0}{{\epsilon }_{q0}^n}$

is a constant determined empirically; n describes soil nonlinearity which is equal to (−0.5) according to the experimental data analyzed by Osman et al. [Osman, A. S.; White, D. J.; Britto, A. M.; Bolton, M. D. Simple prediction of the undrained displacement of a circular surface foundation on non-linear soil. Géotechnique 2007, 57, 729-737, which is fully incorporated by reference], see FIG. 9.

${\epsilon }_q=\sqrt{\frac{2}{9}\left[{\left({\epsilon }_{\mathrm{rr}}-{\epsilon }_{\theta \theta }\right)}^2+{\left({\epsilon }_{\theta \theta }-{\epsilon }_{zz}\right)}^2+{\left({\epsilon }_{zz}-{\epsilon }_{\mathrm{rr}}\right)}^2\right]+\frac{4}{3}\left({\epsilon }_{r\theta }^2+{\epsilon }_{\theta z}^2+{\epsilon }_{\mathrm{zr}}^2\right)}$

represents the deviatoric strain; and ε_q0 is the maximum deviatoric strain with linear elastic behavior which is equal to 10⁻⁵. The soil stiffness G is estimated by calculating the strain at each location followed by the power law.

The pile deflection equation can be solved when the soil and geometry related parameters k, C and m are known; however, these parameters depend on the unknown dimensionless soil function ϕ, which can be estimated by calculating γ₁ and γ₂. Soil displacement is obtained when the initial numbers of these values γ₁ and γ₂, are inserted into Equation (14), from which the parameters k, C and m are obtained as a result of the pile displacement. New values of γ₁ and γ₂ (Equations (15) and (16)) are determined and then inserted into Equation (19) to evaluate ϕ then ν, therefore, an iteration technique is needed to obtain the condition

$\frac{{\gamma }_{old}-{\gamma }_{new}}{{\gamma }_{old}}<0.001.$

This iterative solution methodology (FIG. 10) is programmed in software (e.g., MATLAB R2019a) to obtain the pile head displacements.

FIG. 10 is an exemplary flowchart depicting the iterative solution procedure for calculating pile head displacements. FIG. 10 illustrates a method for analyzing displacement responses of piles subjected to at least one loading condition selected from the group consisting of an axial load, a lateral load, and a bending moment. The disclosed method comprises 1002, inputting values r₀, L, λ, G₀, E_p, P_t, and ν, and either assuming values for or calculating γ_{1, old} and γ_2,old. At 1004, integrations of soil parameters C, k and m are determined. At 1006, the soil dimensionless function ϕ is calculated. At 1008, the pile displacement ν is calculated. At 1010, γ_1,new and γ_2,new are calculated. At 1012, it is determined whether a condition (γ_old−γ_new)/γ_old<0.001 is true. At 1014, steps 1002-1012 are repeated until the condition ((γ_old−γ_new)/γ_old<0.001) is true.

II. Energy-Based Approach: Analysis of a Laterally Loaded Pile in Multi-Layered Non-Linear Elastic Soil Strata

Referring to Prakash, Aryan Ankitha & Arockiasamy, Madasamy. (July, 2022). Energy-Based Approach: Analysis of a Laterally Loaded Pile in Multi-Layered Non-Linear Elastic Soil Strata. Geotechnics. 2. 570-598. 10.3390/geotechnics2030028, which is fully incorporated by reference and made a part hereof, a pile deflection equation is developed.

A pile subjected to lateral loading can be modeled as a Euler-Bernoulli beam, especially for long and slender piles in which the pile shear deformation can be neglected for large slenderness ratios of L/D>10. A laterally loaded pile in an isotropic non-linear elastic multi-layered soil medium is shown in FIG. 11. This disclosure considers a pile of length L with circular cross section of radius r_p. The pile is embedded in n horizontal soil layers, and the pile head is subjected to a horizontal force F₀ accompanied with a bending moment M₀. The horizontal soil layers extend to infinity in the radial direction, and the bottom n^th layer extends to infinity in the vertical direction.

The terms H₁, H₂, H₃ . . . H_n−1 denote the vertical height from the ground surface to the bottom of any layer i. Therefore, the thickness of any layer i is H_i−H_i-1 with H₀=0. Due to the axisymmetric problem behavior, a system of cylindrical coordinates (r-θ-z) is chosen, with the origin coinciding with the center of the pile cross section at the pile head, and the z axis coinciding with the pile axis. The pile head is considered to be free, and the tip of the pile is clamped. Another important assumption to be noted is that there is no slippage or separation between the pile and the surrounding soil and between soil layers. The stresses and the displacement within a soil continuum are shown in FIGS. 5A and 5B. At each point in the soil domain, G and λ were calculated (see FIGS. 6A and 6B).

Soil behaves as a linear elastic material at an extremely low range of strain, and the shear modulus starts degrading at a strain as low as 10⁻⁵. Since the soil particles constantly change their position during the application of a load, the resistance offered by the soil mass against deformation also changes; this results in the change in the value of soil modulus with an increase in strain. A non-linear stress-strain curve can be incorporated in an elastic analysis by properly estimating the secant modulus for a given level of strain (or stress). Soil non-linearity is taken into the account by estimating the decay of soil stiffness with strain.

The variation of the shear stress with strain can be described using two parameters, A and n, which were obtained experimentally using a pressuremeter test, as shown in the equation below:

q=A(ε_q)ⁿ   (20)

where q represents the equivalent shear stress, ε_q is the deviator shear strain. The decay curves of soil stiffness with strain can be divided into three regions, as shown in FIG. 7. The first region in FIG. 7 represents the very small strain where the stiffness G₀ is constant; the second region comprising small strains starts from ε₀ until ε=0.1%, and the third region exceeds ε=0.1%, indicative of large strains.

In the second region, the stiffness decays rapidly, and in the third region, with large strain levels, the stiffness is the smallest, which concludes that soil stiffness is high at small strain and decreases with large strain. FIG. 8 shows the degradation of soil stiffness with increasing strains for different clay types.

The present disclosure is based on a non-linear elastic model developed by Osman et al. [see above], which assumes the decay of soil stiffness with strain using a power law to describe the stress-strain behavior of soil:

$\begin{array}{cc}G={G_0\left(\frac{{\epsilon }_q}{{\epsilon }_{q0}}\right)}^n& \left(21\right)\end{array}$ $\begin{array}{cc}G=a{\epsilon }_q^n& \left(22\right)\end{array}$

where

$a=\frac{G_0}{{\epsilon }_{q0}^n}$

is a constant determined empirically; n describes soil non-linearity, which is equal to (−0.5), according to the experimental data analyzed by Osman et al. [86] (see FIG. 9).

${\epsilon }_q=\sqrt{\frac{2}{9}\left[{\left({\epsilon }_{\mathrm{rr}}-{\epsilon }_{\theta \theta }\right)}^2+{\left({\epsilon }_{\theta \theta }-{\epsilon }_{zz}\right)}^2+{\left({\epsilon }_{zz}-{\epsilon }_{\mathrm{rr}}\right)}^2\right]+\frac{4}{3}\left({\epsilon }_{r\theta }^2+{\epsilon }_{\theta z}^2+{\epsilon }_{\mathrm{zr}}^2\right)}$

represents the deviatoric strain; and ε_q0 is the maximum deviatoric strain with linear elastic behavior, which is equal to 10⁻⁵. Soil stiffness G is estimated by calculating the strain at each location followed by the power law.

For a laterally loaded pile, the displacement at any point within the continuum (FIGS. 5A and 5B) can be expressed as a product of three separable functions, each accounting for one of the dimensions.

u _r =u(z)ϕ_r(r) cos θ  (23a)

u _θ =−u(z)ϕ_θ(r) sin θ  (23b)

u_z=0   (23c)

where u(z) is a displacement function (with a dimension of length) varying with depth z, representing the deflection of the pile axis; u_r, u_θand u_z are the soil displacements in the direction r, θ and z; ϕ_r(r) and ϕ_θ(r) are dimensionless soil displacement functions varying with the radial coordinate r, and θ is measured from a vertical reference section (r=r₀) that contains the applied force vector F₀ (FIGS. 5A and 5B). The above functions are set as ϕ_r(r)=1 and ϕ_θ(r)=1 at r=r₀ (i.e., at the pile-soil interface) and ϕ_r(r)=0 and ϕ_θ(r)=0 at ϕ_θ=∞. Such an assumption ensures the decrease in the displacement (owing to pile deflection) within the soil mass with increasing radial distance from the pile. Thus, ϕ_r and ϕ_θvary between 1 at the pile-soil interface to 0 at an infinite radial distance from the pile. It is also reasonable to assume that the vertical displacement of the pile owing to the lateral load and moment is negligible, and this justifies Equation (23c).

In the derivation of the governing differential equations that can capture the non-linear soil response, the soil within each layer is assumed to be elastic and isotropic but heterogeneous (with respect to r and θ but not with respect to z), with no sliding or separation between the soil layers or between the pile and the soil. By solving the equations (valid for elastic, heterogeneous soil) for different magnitudes of load (with appropriate values and variations of soil modulus), the analysis can trace the non-linear progression of pile deflection (due to soil non-linearity) with increasing applied load.

The total potential energy of the pile-soil system, including both the internal and external potential energies, is given by

$\begin{array}{cc}\prod =\frac{1}{2}{E}_p{I}_p{\int }_0^L{\left(\frac{d^{2}u}{d{z}^2}\right)}^{2}dz+{\int }_0^{\infty }{\int }_0^{2\pi }{\int }_{r_p}^{\infty }\frac{1}{2}{\sigma }_{\mathrm{ij}}{\epsilon }_{\mathrm{ij}}rdrd\theta \mathrm{dz}+{\int }_L^{\infty }{\int }_0^{2\pi }{\int }_0^{r_p}\frac{1}{2}{\sigma }_{\mathrm{ij}}{\epsilon }_{\mathrm{ij}}rdrd\theta dz-F_{o}u{❘}_{z=0}+M_o\frac{du}{dz}{❘}_{z=0}& \left(24\right)\end{array}$

where u is the lateral pile deflection, and σ_ij, ε_ij are the stress and strain tensors in the soil (FIGS. 5A and 5B). The first integral represents the internal potential energy of the pile. The second and third integrals represent the internal potential energy of the continuum. The remaining two terms represent the external potential energy. Soil non-linearity is considered by varying the soil elastic parameters (G and λ) at each discretized nodal point in the soil domain. The stress-strain and strain-displacement relationships at any given nodal point in the soil medium are idealized by the following relationships. The stress-strain relationship is expressed as

$\begin{array}{cc}\left[\begin{array}{c}{\sigma }_{\mathrm{rr}}\\ {\sigma }_{\theta \theta }\\ {\sigma }_{\mathrm{zz}}\\ {\tau }_{r\theta }\\ {\tau }_{\mathrm{rz}}\\ {\tau }_{\theta z}\end{array}\right]=\left[\begin{array}{cccccc}\lambda +2G& \lambda & \lambda & 0& 0& 0\\ \lambda & \lambda +2G& \lambda & 0& 0& 0\\ \lambda & \lambda & \lambda +2G& 0& 0& 0\\ 0& 0& 0& G& 0& 0\\ 0& 0& 0& 0& G& 0\\ 0& 0& 0& 0& 0& G\end{array}\right]\left[\begin{array}{c}{\epsilon }_{\mathrm{rr}}\\ {\epsilon }_{\theta \theta }\\ {\epsilon }_{\mathrm{zz}}\\ {\gamma }_{r\theta }\\ {\gamma }_{\mathrm{rz}}\\ {\gamma }_{\theta z}\end{array}\right]& \left(25\right)\end{array}$

where G and λ are the elastic constants of the soil. The strain-displacement relationship is given by

$\begin{array}{cc}\left[\begin{array}{c}{\epsilon }_{\mathrm{rr}}\\ {\epsilon }_{\theta \theta }\\ {\epsilon }_{\mathrm{zz}}\\ {\gamma }_{r\theta }\\ {\gamma }_{\mathrm{rz}}\\ {\gamma }_{\theta z}\end{array}\right]=\left[\begin{array}{c}-\frac{\partial u_r}{\partial r}\\ -\frac{u_r}{r}-\frac{1}{r}\frac{\partial u_{\theta }}{\partial \theta }\\ -\frac{\partial u_z}{\partial z}\\ -\frac{1}{r}\frac{\partial u_r}{\partial \theta }-\frac{\partial u_{\theta }}{\partial r}+\frac{u_{\theta }}{r}\\ -\frac{\partial u_z}{\partial r}-\frac{\partial u_r}{\partial z}\\ -\frac{1}{r}\frac{\partial u_z}{\partial \theta }-\frac{\partial u_{\theta }}{\partial z}\end{array}\right]\left[\begin{array}{c}-u\left(z\right)\frac{d{\varphi }_r\left(r\right)}{\mathrm{dr}}\cos \theta \\ -u\left(z\right)\frac{{\varphi }_r\left(r\right)-{\varphi }_{\theta }\left(r\right)}{r}\cos \theta \\ 0\\ u\left(z\right)\left\{\frac{{\varphi }_r\left(r\right)-{\varphi }_{\theta }\left(r\right)}{r}+\frac{d{\varphi }_{\theta }\left(r\right)}{\mathrm{dr}}\right\}\sin \theta \\ -\frac{\mathrm{du}\left(z\right)}{\mathrm{dz}}{\varphi }_r\left(r\right)\cos \theta \\ \frac{\mathrm{du}\left(z\right)}{\mathrm{dz}}{\varphi }_{\theta }\left(r\right)\sin \theta \end{array}\right]& \left(26\right)\end{array}$

Combining Equations (25) and (26) gives the strain energy density within any layer:

$\begin{array}{cc}\frac{1}{2}{\sigma }_{\mathrm{ij}}{\epsilon }_{\mathrm{ij}}=\frac{1}{2}\left[\left(\lambda +2G\right){u^2\left(\frac{d{\varphi }_r}{\mathrm{dr}}\right)}^2{\cos }^2\theta +2\lambda u^2\frac{d{\varphi }_r}{\mathrm{dr}}\frac{{\varphi }_r-{\varphi }_{\theta }}{r}{\cos }^2\theta +\left(\lambda +2G\right)u^2\frac{{\left({\varphi }_r-{\varphi }_{\theta }\right)}^2}{r^2}{\cos }^2\theta +G{u}^2\frac{{\left({\varphi }_r-{\varphi }_{\theta }\right)}^2}{r^2}+{\mathrm{Gu}}^2{\left(\frac{d{\varphi }_{\theta }}{\mathrm{dr}}\right)}^2{\sin }^2\theta +2G{u}^2\frac{{\varphi }_r-{\varphi }_{\theta }}{r}\frac{d{\varphi }_{\theta }}{\mathrm{dr}}{\sin }^2\theta +G{\left(\frac{\mathrm{du}}{\mathrm{dz}}\right)}^2{\varphi }_r^2{\cos }^2\theta +{G\left(\frac{\mathrm{du}}{\mathrm{dz}}\right)}^2{\varphi }_{\theta }^2{\sin }^2\theta \right]& \left(27\right)\end{array}$

The equilibrium equations are obtained using the principle of minimum potential energy, according to which δΠ=0. Substituting Equation (27) into Equation (24) and applying δΠ=0, we obtain

$\begin{array}{cc}\mathrm{\delta \Pi }={\int }_0^L{E}_p{I}_p\frac{d^{2}u}{d{z}^2}\delta \left(\frac{d^{2}u}{d{z}^2}\right)\mathrm{dz}+{\int }_0^{\infty }{\int }_{r_p}^{\infty }{\int }_0^{2\pi }\left[\left(\lambda +2G\right)u\delta {u\left(\frac{d{\varphi }_r}{\mathrm{dr}}\right)}^2{\cos }^2\theta +\left(\lambda +2G\right)u^2{\cos }^2\theta \left(\frac{d{\varphi }_r}{\mathrm{dr}}\right)\delta \left(\frac{d{\varphi }_r}{\mathrm{dr}}\right)+2\lambda u\delta u{\cos }^2\theta \frac{d{\varphi }_r}{\mathrm{dr}}\frac{\left({\varphi }_r-{\varphi }_{\theta }\right)}{r}+\lambda u^2{\cos }^2\theta \frac{1}{r}\left({\varphi }_r-{\varphi }_{\theta }\right)\delta \left(\frac{d{\varphi }_r}{\mathrm{dr}}\right)+\lambda u^2{\cos }^2\theta \frac{1}{r}\frac{d{\varphi }_r}{\mathrm{dr}}\delta {\varphi }_r-\lambda u^2{\cos }^2\theta \frac{1}{r}\frac{d{\varphi }_r}{\mathrm{dr}}\delta {\varphi }_{\theta }+\left(\lambda +2G\right)u\delta u{\cos }^2\theta \frac{1}{r^2}{\left({\varphi }_r-{\varphi }_{\theta }\right)}^2+\left(\lambda +2G\right)u^2{\cos }^2\theta \frac{1}{r^2}\left({\varphi }_r-{\varphi }_{\theta }\right)\delta {\varphi }_r-\left(\lambda +2G\right)u^2{\cos }^2\theta \frac{1}{r^2}\left({\varphi }_r-{\varphi }_{\theta }\right)\delta {\varphi }_{\theta }+\mathrm{Gu}\delta u{\sin }^2\theta \frac{1}{r^2}{\left({\varphi }_r-{\varphi }_{\theta }\right)}^2+2\mathrm{Gu}\delta u{\sin }^2\theta \frac{1}{r}\frac{d{\varphi }_{\theta }}{\mathrm{dr}}\left({\varphi }_r-{\varphi }_{\theta }\right)+{\mathrm{Gu}}^2{\sin }^2\theta \frac{1}{r^2}\left({\varphi }_r-{\varphi }_{\theta }\right)\delta {\varphi }_r-G{u}^2{\sin }^2\theta \frac{1}{r^2}\left({\varphi }_r-{\varphi }_{\theta }\right)\delta {\varphi }_{\theta }+{\mathrm{Gu}}^2{\sin }^2\theta \frac{1}{r}\frac{d{\varphi }_{\theta }}{\mathrm{dr}}\delta {\varphi }_r-{\mathrm{Gu}}^2{\sin }^2\theta \frac{1}{r}\frac{d{\varphi }_{\theta }}{\mathrm{dr}}\delta {\varphi }_{\theta }+\mathrm{Gu}\delta u{\sin }^2{\theta \left(\frac{d{\varphi }_{\theta }}{\mathrm{dr}}\right)}^2+{\mathrm{Gu}}^2{\sin }^2\theta \frac{d{\varphi }_{\theta }}{\mathrm{dr}}\delta \left(\frac{d{\varphi }_{\theta }}{\mathrm{dr}}\right)+G{u}^2{\sin }^2\theta \frac{1}{r}\left({\varphi }_r-{\varphi }_{\theta }\right)\delta \left(\frac{d{\varphi }_{\theta }}{\mathrm{dr}}\right)+G\frac{du}{dz}\delta \left(\frac{du}{dz}\right){\varphi }_r^2{\cos }^2\theta +{G\left(\frac{du}{dz}\right)}^2{\cos }^2{\mathrm{\theta \varphi }}_r\delta {\varphi }_r+G\frac{du}{dz}\delta \left(\frac{du}{dz}\right){\varphi }_{\theta }^2{\sin }^2\theta +{G\left(\frac{du}{dz}\right)}^2{\sin }^2{\mathrm{\theta \varphi }}_{\theta }\delta {\varphi }_{\theta }\right]rd\theta d\mathrm{rdz}+\pi r_p^2{\int }_L^{\infty }G\frac{du}{dz}\delta \left(\frac{du}{dz}\right)\mathrm{dz}-F_o\delta u{❘}_{z=0}+M_o\delta \left(\frac{du}{dz}\right){❘}_{z=0}=0& \left(28\right)\end{array}$

The variations δu,δϕ_r, and δϕ_θ are functions of u(z),ϕ_r (r) and ϕ_θ (r) and are independent. To obtain the pile governing equation, from Equation (28), all terms associated with δu, δ(du/dz) will be collected. Furthermore, the terms that are related to δϕ_r and δϕ_θ will be collected to obtain the governing equation of the soil.

The governing equation of piles can be obtained by collecting the terms associated with δu and δ(du/dz) for (0≤z≤L) from Equation (28) and then equating the summation to zero.

$\begin{array}{cc}{\int }_0^L{E}_p{I}_p\frac{d^{2}u}{d{z}^2}\delta \left(\frac{d^{2}u}{d{z}^2}\right)dz+\pi r_p^2{\int }_L^{\infty }G\frac{du}{dz}\delta \left(\frac{du}{dz}\right)dz+{\int }_0^{\infty }{\int }_{r_p}^{\infty }{\int }_0^{2\pi }\left[\left(\lambda +2G\right)u\delta {u\left(\frac{d{\varphi }_r}{\mathrm{dr}}\right)}^2{\cos }^2\theta +2\lambda u\delta u{\cos }^2\theta \frac{d{\varphi }_r}{\mathrm{dr}}\frac{\left({\varphi }_r-{\varphi }_{\theta }\right)}{r}+\mathrm{Gu}\delta u{\sin }^2{\theta \left(\frac{d{\varphi }_{\theta }}{\mathrm{dr}}\right)}^2+G\frac{du}{dz}\delta \left(\frac{du}{dz}\right){\varphi }_r^2{\cos }^2\theta +G\frac{du}{dz}\delta \left(\frac{du}{dz}\right){\varphi }_{\theta }^2{\sin }^2\theta +\mathrm{Gu}\delta u{\sin }^2\theta \frac{1}{r^2}{\left({\varphi }_r-{\varphi }_{\theta }\right)}^2+\left(\lambda +2G\right)u\delta u{\cos }^2\theta \frac{1}{r^2}{\left({\varphi }_r-{\varphi }_{\theta }\right)}^2\right]rd\theta d\mathrm{rdz}-F_o\delta u{❘}_{z=0}+M_o\delta \left(\frac{du}{dz}\right){❘}_{z=0}=0& \left(29\right)\end{array}$

The governing equation becomes

$E_p{I}_p\frac{d^{4}u}{d{z}^4}-C\frac{d^{2}u}{d{z}^2}+ku=0$ $\mathrm{where}$ $C={\int }_{r_p}^{\infty }{\int }_0^{2\pi }{G}_{si}\left({\varphi }_r^2{\cos }^2\theta +{\varphi }_{\theta }^2{\sin }^2\theta \right)\mathrm{rd}\theta$

where i denotes the layer number from the surface to the length of the pile

C _s=∫_rp^∞∫₀^2πG_sj(ϕ_r² cos² θ+ϕ_θ² sin² θ)rdθ

and where j represents the layer number from pile length to infinity.

$k={\int }_{r_p}^{\infty }{\int }_0^{2\pi }\left(\left(\lambda +G\right){\left(\frac{d{\varphi }_r}{dr}\right)}^2{\cos }^2\theta +2\lambda \frac{1}{r}\frac{d{\varphi }_r}{dr}\left({\varphi }_r-{\varphi }_{\theta }\right){\cos }^2\theta +\left(\lambda +G\right)\frac{{\left({\varphi }_r-{\varphi }_{\theta }\right)}^2}{r^2}{\cos }^2\theta +G\frac{{\left({\varphi }_r-{\varphi }_{\theta }\right)}^2}{r^2}{\sin }^2\theta +2G\frac{\left({\varphi }_r-{\varphi }_{\theta }\right)}{r}\frac{d{\varphi }_{\theta }}{dr}{\sin }^2\theta +{G\left(\frac{d{\varphi }_{\theta }}{dr}\right)}^2{\sin }^2\theta \right)\mathrm{rd}\theta \mathrm{dr}$

The fourth-order differential Equation (30) can be solved using a central finite difference scheme represented as:

$E_p{I}_p\left(\frac{u_{i-2}-4u_{i-1}+6u_i-4u_{i+1}+u_{i+2}}{\Delta z^4}\right)-C\left(\frac{-u_{i-1}+u_{i+1}}{\Delta z}\right)+k{u}_i=0$

where i denotes the ith node in z direction, and Δz is the distance between two nodes. At each point in the soil domain, λ and G were calculated (see FIGS. 6A and 6B).

Similar to the pile governing equation, the soil governing equation can be formulated by considering the variations of ϕ_r(r) and ϕ_θ(r). Collecting the terms associated with ϕ_r and ϕ_θ over the domain r₀≤r≤∞ and equating the summation to zero yields

$\begin{array}{cc}{\int }_0^{\infty }{\int }_{r_p}^{\infty }{\int }_0^{2\pi }\left[\left(\lambda +2G\right)u^2{\cos }^2\theta \left(\frac{d{\varphi }_r}{\mathrm{dr}}\right)\delta \left(\frac{d{\varphi }_r}{\mathrm{dr}}\right)+{\mathrm{Gu}}^2{\sin }^2\theta \frac{1}{r}\left(\frac{d{\varphi }_{\theta }}{\mathrm{dr}}\right)\delta {\varphi }_r+\lambda u^2{\cos }^2\theta \frac{\left({\varphi }_r-{\varphi }_{\theta }\right)}{r}\delta \left(\frac{d{\varphi }_r}{\mathrm{dr}}\right)+\lambda u^2{\cos }^2\theta \frac{1}{r}\left(\frac{d{\varphi }_r}{\mathrm{dr}}\right)\delta {\varphi }_r+{\mathrm{Gu}}^2{\sin }^2\theta \frac{1}{r^2}\left({\varphi }_r-{\varphi }_{\theta }\right)\delta {\varphi }_r+\left(\lambda +2G\right)u^2{\cos }^2\theta \frac{1}{r^2}\left({\varphi }_r-{\varphi }_{\theta }\right)\delta {\varphi }_r+G{\cos }^2{\theta \left(\frac{du}{dz}\right)}^2{\varphi }_r\delta {\varphi }_r\right]rdrd\theta dz=0& \left(31\right)\end{array}$ $\begin{array}{cc}{\int }_0^{\infty }{\int }_{r_p}^{\infty }{\int }_0^{2\pi }\left[-G{u}^2{\sin }^2\theta \left(\frac{d{\varphi }_{\theta }}{\mathrm{dr}}\right)\delta {\varphi }_{\theta }-\lambda u^2{\cos }^2\theta \frac{1}{r}\left(\frac{d{\varphi }_r}{\mathrm{dr}}\right)\delta {\varphi }_{\theta }-G{u}^2{\sin }^2\theta \frac{1}{r^2}\left({\varphi }_r-{\varphi }_{\theta }\right)\delta {\varphi }_{\theta }-\left(\lambda +2G\right)u^2{\cos }^2\theta \frac{1}{r^2}\left({\varphi }_r-{\varphi }_{\theta }\right)\delta {\varphi }_{\theta }+G{u}^2{\sin }^2\theta \left(\frac{d{\varphi }_{\theta }}{\mathrm{dr}}\right)\delta \left(\frac{d{\varphi }_{\theta }}{\mathrm{dr}}\right)+{\mathrm{Gu}}^2{\sin }^2\theta \frac{\left({\varphi }_r-{\varphi }_{\theta }\right)}{r}\delta \left(\frac{d{\varphi }_{\theta }}{\mathrm{dr}}\right)+G{\sin }^2{\theta \left(\frac{\mathrm{du}}{\mathrm{dz}}\right)}^2{\varphi }_{\theta }\delta {\varphi }_{\theta }\right]rdrd\theta dz=0& \left(32\right)\end{array}$

Simplifying Equations (31) and (32) gives

$\begin{array}{cc}{\int }_{r_p}^{\infty }\left(m_{s1}\left(\frac{d{\varphi }_r}{\mathrm{dr}}\right)\delta \left(\frac{d{\varphi }_r}{\mathrm{dr}}\right)+m_{s3}\frac{1}{r}\left(\frac{d{\varphi }_r}{\mathrm{dr}}\right)\delta {\varphi }_r+m_{s3}\frac{\left({\varphi }_r-{\varphi }_{\theta }\right)}{r}\delta \left(\frac{d{\varphi }_r}{\mathrm{dr}}\right)+m_{s1}\frac{1}{r^2}\left({\varphi }_r-{\varphi }_{\theta }\right)\delta {\varphi }_r+m_{s2}\frac{1}{r^2}\left({\varphi }_r-{\varphi }_{\theta }\right)\delta {\varphi }_r+m_{s2}\frac{1}{r}\left(\frac{d{\varphi }_{\theta }}{\mathrm{dr}}\right)\delta {\varphi }_r+n_{s1}{\varphi }_r\delta {\varphi }_r\right)dr=0& \left(33\right)\end{array}$ $\begin{array}{cc}{\int }_p^{\infty }\left(-m_{s2}\left(\frac{d{\varphi }_{\theta }}{\mathrm{dr}}\right)\delta \left(\frac{d{\varphi }_{\theta }}{\mathrm{dr}}\right)-m_{s3}\frac{1}{r}\left(\frac{d{\varphi }_r}{\mathrm{dr}}\right)\delta {\varphi }_{\theta }-m_{s2}\frac{1}{r^2}\left({\varphi }_r-{\varphi }_{\theta }\right)\delta {\varphi }_{\theta }-m_{s1}\frac{1}{r^2}\left({\varphi }_r-{\varphi }_{\theta }\right)\delta {\varphi }_{\theta }+m_{s2}\left(\frac{d{\varphi }_{\theta }}{\mathrm{dr}}\right)\delta \left(\frac{d{\varphi }_{\theta }}{\mathrm{dr}}\right)+m_{s2}\frac{\left({\varphi }_r-{\varphi }_{\theta }\right)}{r}\delta \left(\frac{d{\varphi }_{\theta }}{\mathrm{dr}}\right)n_{s2}{\varphi }_{\theta }\delta {\varphi }_{\theta }\right)dr=0& \left(34\right)\end{array}$

where n_s1, n_s2, m_s1, m_s2 and m_s3 for homogeneous and non-homogeneous soil comprising multi-layer soil summations with i representing the number of layers and n being the last layer number are given below

$\begin{array}{cc}{m}_{s1}\left(r\right)={\int }_0^{\infty }{\int }_0^{2\pi }\left({\lambda }_s+2G_s\right)u^2{\cos }^2\theta \mathrm{rd}\theta \mathrm{dz}={\sum }_{i=n}^n{\int }_{H_{i-1}}^{H_i}{\int }_0^{2\pi }\left({\lambda }_{si}+2G_{si}\right)u_i^2{\cos }^2\theta \mathrm{rd}\theta \mathrm{dz}& \left(35a\right)\end{array}$ $\begin{array}{cc}{m}_{s2}\left(r\right)={\int }_0^{\infty }{\int }_0^{2\pi }{G}_s{u}^2{\sin }^2\theta rd\theta dz={\sum }_{i=n}^n{\int }_{H_{i-1}}^{H_i}{\int }_0^{2\pi }{G}_{si}{u}_i^2{\sin }^2\theta \mathrm{rd}\theta \mathrm{dz}& \left(35b\right)\end{array}$ $\begin{array}{cc}{m}_{s3}\left(r\right)={\int }_0^{\infty }{\int }_0^{2\pi }{\lambda }_s{u}^2{\cos }^2\theta \mathrm{rd}\theta \mathrm{dz}={\sum }_{i=n}^n{\int }_{H_{i-1}}^{H_i}{\int }_0^{2\pi }{\lambda }_{si}{u}_i^2{\cos }^2\theta \mathrm{rd}\theta \mathrm{dz}& \left(35c\right)\end{array}$ $\begin{array}{cc}{n}_{s1}\left(r\right)={\int }_0^{\infty }{\int }_0^{2\pi }{G}_s\left(\frac{\mathrm{du}}{\mathrm{dz}}\right){\cos }^2\theta rd\theta dz={\sum }_{i=n}^n{\int }_{H_{i-1}}^{H_i}{\int }_0^{2\pi }{G_{si}\left(\frac{d{u}_i}{dz}\right)}^2{\cos }^2\theta \mathrm{rd}\theta \mathrm{dz}& \left(35d\right)\end{array}$ $\begin{array}{cc}{n}_{s1}\left(r\right)={\int }_0^{\infty }{\int }_0^{2\pi }{G}_s\left(\frac{\mathrm{du}}{\mathrm{dz}}\right){\sin }^2\theta \mathrm{rd}\theta \mathrm{dz}={\sum }_{i=n}^n{\int }_{H_{i-1}}^{H_i}{\int }_0^{2\pi }{G_{si}\left(\frac{d{u}_i}{dz}\right)}^2{\sin }^2\theta \mathrm{rd}\theta \mathrm{dz}& \left(35e\right)\end{array}$

Equations (33) and (34) can be simplified further using integration by parts of the terms that consist of δ((dϕ_r)/dr) and δ((dϕ_θ)/dr), respectively. The governing equation is then obtained by dividing the resulting equations of (33) by (−m_s1) and (34) by (−m_m2).

$\begin{array}{cc}\frac{d^2{\varphi }_r}{{\mathrm{dr}}^2}+\frac{1}{m_{S1}}\frac{d{m}_{S1}}{\mathrm{dr}}\frac{d{\varphi }_r}{\mathrm{dr}}-\left\{\frac{1}{r^2}\frac{m_{S1}+m_{S2}+m_{S3}}{m_{S1}}-\frac{1}{r}\frac{1}{m_{S1}}\frac{d{m}_{S3}}{\mathrm{dr}}+\frac{n_{S1}}{m_{S1}}\right\}{\varphi }_r=\frac{m_{S2}+m_{S3}}{m_{S1}}\frac{1}{r}\frac{d{\varphi }_{\theta }}{\mathrm{dr}}-\left\{\frac{1}{r^2}\frac{m_{S1}+m_{S2}+m_{S3}}{m_{S1}}-\frac{1}{r}\frac{1}{m_{S1}}\frac{d{m}_{S3}}{\mathrm{dr}}\right\}{\varphi }_{\theta }& \left(36\right)\end{array}$ $\begin{array}{cc}\frac{d^2{\varphi }_{\theta }}{{\mathrm{dr}}^2}+\frac{1}{m_{S2}}\frac{d{m}_{S2}}{\mathrm{dr}}\frac{d{\varphi }_{\theta }}{\mathrm{dr}}-\left\{\frac{1}{r^2}\frac{m_{S1}}{m_{S2}}+\frac{1}{r}\frac{1}{m_{S2}}\frac{d{m}_{S2}}{\mathrm{dr}}+\frac{n_{S2}}{m_{S2}}\right\}{\varphi }_{\theta }=\frac{m_{S2}+m_{S3}}{m_{S2}}\frac{1}{r}\frac{d{\varphi }_r}{\mathrm{dr}}-\left\{\frac{1}{r^2}\frac{m_{S1}}{m_{S2}}+\frac{1}{r}\frac{1}{m_{S2}}\frac{d{m}_{S2}}{\mathrm{dr}}\right\}{\varphi }_r& \left(37\right)\end{array}$

Considering the soil boundary conditions ϕ_r=1 at r=r₀; ϕ_r=0 at r→∞; ϕ_θ=1 at r=r₀; and ϕ_θ=0 at r→∞, the governing differential equations of soil can be rewritten as

$\begin{array}{cc}\frac{d^2{\varphi }_r}{{\mathrm{dr}}^2}+y_1\frac{d{\varphi }_r}{\mathrm{dr}}-y_2{\varphi }_r=y_3\frac{d{\varphi }_{\theta }}{\mathrm{dr}}-y_4{\varphi }_{\theta }& \left(38\right)\end{array}$ $\begin{array}{cc}\frac{d^2{\varphi }_{\theta }}{{\mathrm{dr}}^2}+y_5\frac{d{\varphi }_{\theta }}{\mathrm{dr}}-y_6{\varphi }_{\theta }=-y_7\frac{d{\varphi }_r}{\mathrm{dr}}-y_8{\varphi }_r& \left(39\right)\end{array}$ $\mathrm{where}$ $\begin{array}{cc}{y}_1=\frac{1}{m_{S1}}\frac{d{m}_{S1}}{\mathrm{dr}}& \left(40a\right)\end{array}$ $\begin{array}{cc}{y}_1=\frac{1}{m_{S1}}\frac{d{m}_{S1}}{\mathrm{dr}}& \left(40b\right)\end{array}$ $\begin{array}{cc}{y}_3=\frac{1}{r}\frac{m_{S2}+m_{S3}}{m_{S1}}& \left(40c\right)\end{array}$ $\begin{array}{cc}{y}_4=\frac{1}{r^2}\frac{m_{S1}+m_{S2}+m_{S3}}{m_{S1}}-\frac{1}{r}\frac{1}{m_{S1}}\frac{d{m}_{S3}}{\mathrm{dr}}& \left(40d\right)\end{array}$ $\begin{array}{cc}{y}_5=\frac{1}{m_{S2}}\frac{d{m}_{S2}}{\mathrm{dr}}& \left(40e\right)\end{array}$ $\begin{array}{cc}{y}_6=\frac{1}{r^2}\frac{m_{S1}}{m_{S2}}+\frac{1}{r}\frac{1}{m_{S2}}\frac{d{m}_{S2}}{\mathrm{dr}}+\frac{n_{S2}}{m_{S2}}& \left(40f\right)\end{array}$ $\begin{array}{cc}{y}_7=\frac{1}{r}\frac{m_{S2}+m_{S3}}{m_{S2}}& \left(40g\right)\end{array}$ $\begin{array}{cc}{y}_8=\frac{1}{r^2}\frac{m_{S1}}{m_{S2}}+\frac{1}{r}\frac{1}{m_{S2}}\frac{d{m}_{S2}}{\mathrm{dr}}& \left(40h\right)\end{array}$

The disclosed embodiments operate on an iterative basis. The pile deflection equation (Equation (30)) can be solved when the soil and geometry related parameters k, C and m are known; however, these parameters are dependent on the unknown dimensionless soil functions, ϕ_r and ϕ_θ, which can be estimated by calculating y₁, y₂, y₃, y₄, y₅, y₆, y₇ and y₈. Soil deformation is calculated in radial and circumference directions, then calculated in the depth direction. To obtain the soil displacement, the initial numbers of these values, y₁, y₂, y₃, y₄, y₅, y₆, y₇ and y₈, must be assumed. Then, they are inserted into Equations (38) and (39), from which the soil parameters k and C are obtained as a result of the pile displacement. New values of y₁, y₂, y₃, y₄, y₅, y₆, y₇ and y₈ can then be inserted into Equations (38) and (39) to evaluate ϕ_r and ϕ_θwhich are then inserted into Equation (30) to obtain the pile displacement u, so an iteration technique is needed to satisfy the condition

$\frac{{\gamma }_{old}-{\gamma }_{new}}{{\gamma }_{old}}<0.001.$

FIG. 12 shows the flow chart of the iterative solution procedure for a laterally-loaded pile. The disclosed method comprises 1202, inputting values r₀, L, λ, G₀, E_p, I_p, and applying force and moments, Q₀ and M₀. At 1204, choosing initial values of y₁, y₂, y₃, y₄, y₅, y₆, y₇ and y₈, and calculating ϕ_r and ϕ_θ. At 1206, calculating deviatoric strain ε_q, ad using the deviatoric strain to calculate soil stiffness, G, where

$G={G_0\left(\frac{{\epsilon }_q}{{\epsilon }_{\mathrm{qo}}}\right)}^n.$

At 1208, calculate n_s1, n_s2, m_s1, m_s2 and m_s3, and then calculate the pile displacement u. At 1210, use the pile displacement u to calculate y₁^new y₂^new y₃^new y₄^new y₅^new y₆^new y₇^new and y₈^new. At 1212, it is determined whether a condition ((y_old−y_new)/y_old<0.001)is true. At 1214, steps 1202-1212 are repeated until the condition ((y_old−y_new)/y_old<0.001) is true.

In examples, a computer system may be included and/or operated within which a set of instructions for implementing the disclosed embodiments of the methods are execute to form a machine and for causing the machine (e.g., computer, cytometer, smartphone) to perform any one or more of the methodologies discussed herein, may be executed. In alternative embodiments, the machine may be connected (e.g., networked) to other machines in a local area network (LAN), an intranet, an extranet, or the Internet. The machine may operate in the capacity of a server or a client machine in a client-server network environment, or as a peer machine in a peer-to-peer (or distributed) network environment. It is understood that the machine may be a personal computer (PC), a tablet PC, a set-top box (STB), a personal digital assistant (PDA), a cellular telephone, a web appliance, a server, a network router, switch or bridge, or any machine capable of executing the computer application or another other set of instructions (sequential or otherwise) that specify actions to be taken by that machine. Further, the term “machine” shall also be taken to include any collection of machines that individually or jointly execute a set (or multiple sets) of instructions to perform any one or more of the methodologies discussed herein.

The exemplary machine may include a processing device, a main memory (e.g., read-only memory (ROM), flash memory, dynamic random access memory (DRAM) (such as synchronous DRAM (SDRAM) or Rambus DRAM (RDRAM), etc.), a static memory (e.g., flash memory, static random access memory (SRAM), etc.), and a data storage device, which communicate with each other via a bus.

Processing device represents one or more general-purpose processing devices such as a microprocessor, central processing unit, or the like. More particularly, the processing device may be complex instruction set computing (CISC) microprocessor, reduced instruction set computer (RISC) microprocessor, very long instruction word (VLIW) microprocessor, or a processor implementing other instruction sets, or processors implementing a combination of instruction sets. Processing device may also be one or more special-purpose processing devices such as an application specific integrated circuit (ASIC), a field programmable gate array (FPGA), a digital signal processor (DSP), network processor, or the like. Processing device is configured to execute listings manager logic for performing the operations and steps discussed herein.

Computer system may further include a network interface device (e.g., GUI). Computer system also may include a video display unit (e.g., a liquid crystal display (LCD) or a cathode ray tube (CRT)), an alphanumeric input device (e.g., keyboard, keypad), a cursor control device (e.g., mouse. touchpad), and a signal generation device (e.g., a speaker). For example, FIGS. 13, 14 and 15 are exemplary annotated screenshots of graphical user interfaces of embodiments of software according to the invention.

Data storage device may include a machine-readable storage medium (or more specifically a computer-readable storage medium) having one or more sets of instructions (e.g., reference generation module) embodying any one or more of the methodologies of functions described herein. The reference generation module may also reside, completely or at least partially, within main memory and/or within processing device during execution thereof by computer system; main memory and processing device also constituting machine-readable storage media. The reference generation module may further be transmitted or received over a network via network interface device.

Machine-readable storage medium may also be used to store the device queue manager logic persistently. While a non-transitory machine-readable storage medium is shown in an exemplary embodiment to be a single medium, the term “machine-readable storage medium” should be taken to include a single medium or multiple media (e.g., a centralized or distributed database, and/or associated caches and servers) that store the one or more sets of instructions. The term “machine-readable storage medium” shall also be taken to include any medium that is capable of storing or encoding a set of instruction for execution by the machine and that causes the machine to perform any one or more of the methodologies of the present invention. The term “machine-readable storage medium” shall accordingly be taken to include, but not be limited to, solid-state memories, and optical and magnetic media.

The components and other features described herein can be implemented as discrete hardware components or integrated in the functionality of hardware components such as ASICs, FPGAs, DSPs or similar devices. In addition, these components can be implemented as firmware or functional circuitry within hardware devices. Further, these components can be implemented in any combination of hardware devices and software components.

Some portions of the detailed descriptions are presented in terms of algorithms and symbolic representations of operations on data bits within a computer memory. These algorithmic descriptions and representations are the means used by those skilled in the data processing arts to most effectively convey the substance of their work to others skilled in the art. An algorithm is here, and generally, conceived to be a self-consistent sequence of steps leading to a desired result. The steps are those requiring physical manipulations of physical quantities. Usually, though not necessarily, these quantities take the form of electrical or magnetic signals capable of being stored, transferred, combined, compared, and otherwise manipulated. It has proven convenient at times, principally for reasons of common usage, to refer to these signals as bits, values, elements, symbols, characters, terms, numbers, or the like.

The instructions may include, for example, computer-executable instructions and data which cause a general purpose computer, special purpose computer, or special purpose processing device to perform a certain function or group of functions. Computer-executable instructions also include program modules that are executed by computers in stand-alone or network environments. Generally, program modules include routines, programs, objects, components, and data structures, and the like that perform particular tasks or implement particular abstract data types. Computer-executable instructions, associated data structures, and program modules represent examples of the program code means for executing steps of the methods disclosed herein. The particular sequence of such executable instructions or associated data structures represents examples of corresponding acts for implementing the functions described therein.

In the aforementioned description, numerous details are set forth. It will be apparent, however, to one skilled in the art, that the disclosure may be practiced without these specific details. In some instances, well-known structures and devices are shown in block diagram form, rather than in detail, in order to avoid obscuring the disclosure.

The disclosure is related to an apparatus for performing the operations herein. This apparatus may be specially constructed for the required purposes or it may comprise a general purpose computing device selectively activated or reconfigured by a computer program stored therein. Such a computer program may be stored in a non-transitory computer readable storage medium, such as, but not limited to, any type of disk including floppy disks, optical disks, CD-ROMs and magnetic-optical disks, read-only memories (ROMs), random access memories (RAMs), EPROMs, EEPROMs, magnetic or optical cards, flash memory devices including universal serial bus (USB) storage devices (e.g., USB key devices) or any type of media suitable for storing electronic instructions, each of which may be coupled to a computer system bus.

While the invention has been described in detail and with reference to specific examples thereof, it will be apparent to one skilled in the art that various changes and modifications can be made therein without departing from the spirit and scope thereof.

Claims

1. A method for analyzing displacement responses of a pile subjected to at least one loading condition selected from the group consisting of an axial load, a lateral load, and a bending moment, said method comprising: (a) inputting values r0, L, λ, G0, Ep, Pt, and ν;(b) calculating γ1,old and y2,old;(c) determining integrations of soil parameters C, k and m;(d) calculating soil dimensionless function ϕ;(e) calculating pile displacement ν;(f) calculating γ1,new and γ2,new;(g) determining whether a condition (γold−γnew)/γold<0.001 is true; and(h) repeating steps (a)-(g) until the condition is true.

2. The method of claim 1, wherein the pile is installed in nonlinear soil.

3. The method of claim 2, wherein the soil is assumed to be nonlinear elastic and perfectly plastic.

4. The method of claim 2, wherein soil moduli (λ and G0) are assumed to vary in radial, circumference and depth directions according to strain and stress levels.

5. The method of claim 1, wherein the pile comprises a plurality of piles.

6. The method of claim 1, wherein the method is coded into executable instructions, sored in a memory, and executed by a processor.

7. The method of claim 5, wherein the executable instructions are executed using a finite difference method (FDM).

8. The method of claim 1, wherein soil moduli (λ and G0) are assumed to vary in radial, circumference and depth directions according to strain and stress levels.

9. The method of claim 1, wherein the loading condition is the axial (vertical) load and a head of the pile (pile head) is assumed to be free, and a tip of the pile is assumed to be clamped.

10. The method of claim 8, wherein calculating pile displacement ν comprises determining pile head deflection caused by the action of the vertical load at the pile head.

11. A non-transitory machine-readable storage medium comprising instructions executable by a processing resource of a computing device to cause the processing resource to perform a method of analyzing displacement responses of a pile subjected to an axial loading condition, said method comprising: (a) inputting values r0, L, λ, G0, Ep, Pt, and ν;(b) calculating γ1,old and γ2,old;(c) determining integrations of soil parameters C, k and m;(d) calculating soil dimensionless function ϕ;(e) calculating pile displacement ν;(f) calculating γ1,new and γ2,new;(g) determining whether a condition (γold−γnew)/γold<0.001 is true; and(h) repeating steps (a)-(g) until the condition is true.

12. The non-transitory machine-readable storage medium of claim 11, wherein the pile is installed in nonlinear soil.

13. The non-transitory machine-readable storage medium of claim 12, wherein the soil is assumed to be nonlinear elastic and perfectly plastic.

14. The non-transitory machine-readable storage medium of claim 12, wherein soil moduli (λ and G0) are assumed to vary in radial, circumference and depth directions according to strain and stress levels.

15. The non-transitory machine-readable storage medium of claim 11, wherein the pile comprises a plurality of piles.

16. The non-transitory machine-readable storage medium of claim 15, wherein the instructions executable by a processing resource are executed using a finite difference method (FDM).

17. The non-transitory machine-readable storage medium of claim 11, wherein soil moduli (λ and G0) are assumed to vary in radial, circumference and depth directions according to strain and stress levels.

18. The non-transitory machine-readable storage medium of claim 11, wherein the loading condition is the axial (vertical) load and a head of the pile (pile head) is assumed to be free, and a tip of the pile is assumed to be clamped.

19. The non-transitory machine-readable storage medium of claim 18, wherein calculating pile displacement v comprises determining pile head deflection caused by the action of the vertical load at the pile head.

20. A method for analyzing displacement responses of a pile subjected to a lateral loading condition, said method comprising the following steps: (a) inputting values r0, L, λ, G0, Ep, and Ip, and applying force and moments, Q0 and M0;(b) choosing initial values of y1, y2, y3, y4, y5, y6, y7 and y8, and calculating ϕr and ϕθ;(c) calculating deviatoric strain εq, and using the deviatoric strain to calculate soil stiffness, G, where

21. A non-transitory machine-readable storage medium comprising instructions executable by a processing resource of a computing device to cause the processing resource to perform a method of analyzing displacement responses of a pile subjected to a lateral loading condition, said method comprising: (a) inputting values r0, L, λ, G0, Ep, and Ip, and applying force and moments, Q0 and M0;(b) choosing initial values of y1, y2, y3, y4, y5, y6, y7 and y8, and calculating ϕr and ϕθ;(c) calculating deviatoric strain εq, and using the deviatoric strain to calculate soil stiffness, G, where