Boundary Value Problem for Analysis of Portal Double-Row Stabilizing Piles

Abstract

This paper presents a new numerical approach for computing the internal force and displacement of portal double-row piles used to stabilize potential landslide. First, the new differential equations governing the mechanical behaviour of the stabilizing pile are formulated and the boundary conditions are mathematically specified. Then, the problem is numerically solved by the high-accuracy Runge-Kutta finite difference method. A program package has been developed in MATLAB depending on the proposed algorithm. Illustrative examples are presented to demonstrate the validity of the developed program. In short, the proposed approach is a practical new idea for analyzing the portal double-row stabilizing pile as a useful supplement to traditional methods such as FEM.

(ProQuest: ... denotes non-US-ASCII text omitted.)

Academic Editor:Fayun Liang

State Key Laboratory of Geomechanics and Geotechnical Engineering, Institute of Rock and Soil Mechanics, Chinese Academy of Sciences, Wuhan 430071, China

Received 10 May 2013; Revised 26 July 2013; Accepted 29 July 2013

This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

1. Introduction

Double-row stabilizing piles have been widely used in slope reinforcement engineering and treatment of landslide geological disasters, which have some advantages such as larger rigidity, less displacement in the top of the piles, and large resisting force. Existing methods for the analysis of double-row stabilizing piles can be generally classified into the following two categories [1-6]: (1) coupled method (continuum analysis) that simultaneously solves pile response and slope stability [7]; (3) uncoupled method which deals with pile and slope separately. In the uncoupled method, pile-soil interaction is commonly represented by equivalent Winkler or p -y springs [8-13].

As for coupled method, the finite element method is certainly the most comprehensive approach to study pile-slope stability. However, its use generally requires high numerical costs and accurate measurements of material properties. This makes the use of this method rather unattractive for practical applications [6].

To date, in practical engineering applications, the uncoupled method is the most widely used approach to design the double-row reinforcing piles to increase slope stability due to its simplicity of use. First, the lateral force acting on the pile segment above the slip surface due to soil movement is evaluated usually by the limit equilibrium method. Second, the response of the double-row pile subjected to lateral loading is analysed by FEM modeling it as a beam resting on linear or nonlinear soil/rock spring supports. The FEM modeling is reasonably accurate but complicated and time consuming.

In this paper, a new uncoupled method to compute the response of portal double-row piles subjected to lateral earth pressure loading based on new boundary value problem approach is introduced. First, the new governing differential equations including six variables (three internal forces and three displacements) are formulated and the boundary condition is specified. Second, the high-accuracy Runge-Kutta differential method is used to solve the corresponding system of differential equations to obtain the pile's internal forces and displacements. A program for pile response analysis and graphics edit is developed. At last, the program was verified against the FEM analysis results in terms of pile deflection, bending moment, and shear force along the length of the pile.

The objective of this study is to provide an alternative method for the design of portal double-row pile used for slope stabilization or earth retaining.

2. Derivation of Governing Differential Equations

In order to solve the complicated engineering problem of the response of double-row stabilizing pile under laterally loading by using accurately mathematical model, we often need to define its boundary value problem which involves the governing differential equations and corresponding boundary conditions. Then, closed-form or numerical solutions for the engineering problem can be obtained by many appropriate mathematical methods.

Under the scheme of uncoupled analysis of the pile (as shown in Figure 1), the new governing differential equations for stabilizing piles embedded in slope will be developed in the following according to the general principles of solid and structural mechanics including static force equilibrium, deformation compatibility, and constitutive relationship.

Figure 1: The uncoupled analysis model for the double-row portal stabilizing piles.

[figure omitted; refer to PDF]

2.1. The Loading Condition

Below the slip surface, the sliding bed is assumed to be stable and cannot move. Before the active force induced by the sliding mass act upon the pile segment above the slip surface, the earth pressure acting at the front and back sides of the pile segment below the slip surface is in equilibrium. So it can be neglected. Only the active force induced by the movement of the sliding mass is considered in the calculation model, which are loaded on the pile segment above the slip surface. Their distribution along the pile shaft can be assumed as uniform, triangular, trapezoidal, and rectangular profiles. Then, the reaction force acting at the pile segment below the slip surface is calculated. These active forces and reaction forces considered are in equilibrium [14].

2.2. The Soil-Pile Interaction Model

Due to its simplicity and reasonable accuracy, the Winkler foundation model is adopted in current analysis to describe the pile-soil interaction behavior. The Winkler method assumed that the substratum is composed of independent horizontal springs. Under the Winkler hypothesis, the soil reaction pressures (p ) acting on the pile can be modeled by discrete independent linear or nonlinear springs in the form of the following equation: [figure omitted; refer to PDF] where k is the spring constant, also called the modulus of horizontal subgrade reaction (it has a unit of force/length³ ). The main difference between the different Winkler foundation models available is in the selection of the foundation stiffness coefficients. p is the horizontal soil reaction pressure (it has a unit of force/length² ). w is the horizontal displacement (it has a unit of length).

2.3. The New Equilibrium Differential Equations

Let us consider an isolated free portion of pile, as shown in Figure 2, having a infinitesimal length of ds and acted upon by external distributed normal load _qn and tangential load _qτ . The free segment can be imagined to be cut out of the pile, and the internal forces (M,N,Q ) in the original pile may become external forces on the isolated free portion.

Figure 2: The free-body diagram of infinitesimal isolated segment of pile.

[figure omitted; refer to PDF]

Figure 2 also shows the coordinate system used in this paper for introducing the governing equations.

We adopt sign conventions so that the six variables as shown in Figure 2 are positive. The sign convention adopted for forces is that positive sign indicates tensile axial force N , the positive shearing force Q should be directed so that they will tend to rotate the element counterclockwise, and the positive bending moment M will tend to make the element concave leftward. The sign convention adopted for displacements is that the positive normal displacement v points outward normal, the positive tangential displacement u points right when facing outward normal, and the positive [straight phi] is counterclockwise. The lateral pressure _qn is considered positive when applied from left to right. The lateral pressure _qτ is considered positive when applied from up to down.

Thus, considering the equilibrium of the above infinitesimal pile segment, as it bends under the action of the applied loads (shown in Figure 2), we arrive at two force equilibrium equations in the directions of u and v and one moment equilibrium equation: [figure omitted; refer to PDF]

Simplifying the above equations and neglecting the higher-order term and the term with the square of the differential, the equilibrium equations now take the following forms: [figure omitted; refer to PDF] where H={1,contact,0,noncontact, the "one" indicates that Winkler reaction force exists and the "zero" indicates that Winkler reaction force does not exist. S is the position coordinate; D is perimeter of the cross-section. b is the width of the cross-section, _ks is the Winkler modulus of vertical subgrade reaction, and _kn is the Winkler modulus of horizontal subgrade reaction.

For the sake of convenience of formula deducing, let [figure omitted; refer to PDF]

The set of equations of equilibrium (3) can be rewritten in the following matrix form: [figure omitted; refer to PDF] where B=[0000000-10] ; L=[_ks HD000_kn Hb0000] .

2.4. Geometric and Constitutive Equations

When the deformation (du- , dv- , d[straight phi]- ) of the differential element (shown in Figure 2) induced by the internal forces (N,Q,M ) is considered given, the corresponding strains can be expressed as [figure omitted; refer to PDF] According to the related theory of elastic beam, the internal forces (N,Q,M ) can be related to strains as in the following linear constitutive equation: [figure omitted; refer to PDF] where α is the constant related to the shape of pile cross-section (α=6/5 for rectangular cross-section; α=10/9 for circular cross-section); A is the cross-sectional area; EI is the flexural rigidity of the pile's cross-section.

Considering the deformation, dU- can be decomposed into two parts. One part is the dZ induced by the displacement on its direction and another part is the projection of other displacement onto this direction which takes the form B-Zds , where B- is the a undetermined third-order square matrix. Then, the deformation dU- can be expressed as [figure omitted; refer to PDF] Applying the principle of virtual work to the isolated differential element of pile (shown in Figure 2), B- can be determined. We suppose that each point of the body is given an infinitesimal virtual displacement δZ satisfying displacement boundary conditions where prescribed. The virtual deformation associated with the infinitesimal virtual displacement is δU- . The virtual work of the external surface forces is -∫^P-T (δZ)ds , where P-=P+LZ . The virtual work of the internal forces is ∫^XT d(δU-) . By equating the external work to the internal work, we have -∫^P-T (δZ)ds=∫^XT d(δU-) . Substituting (5) and (8) into the above equation and simplifying yields ∫^XT (^BT -B-)δZds=0 . Since this equation is satisfied for arbitrary δZ , the terms in the brackets in the integral must vanish at every point which means that B-=^BT . At last, we develop the following geometric and constitutive equations: [figure omitted; refer to PDF]

2.5. Matrix Form of the Governing Equations

For the sake of convenience of problem solving, combining the three equilibrium differential equations (3) and the three geometric and constitutive equations (9) leads to a system of six equations: [figure omitted; refer to PDF]

Then, we can use matrix notation to present the above equations (10) in the form: [figure omitted; refer to PDF] where K is the coefficients square matrix of six order and X and p represent two column matrices.

This system of six independent differential equations can be solved for six unknown functions (three independent forces and three independent displacements).

3. Boundary Conditions

As mentioned above, there are total six unknowns to be determined (N,Q,M,u,v,[straight phi] ). Therefore, six boundary conditions are needed for the problem solving.

The boundary conditions for (12) are determined according to the way in which the pile's head and base are supported or restrained. There are three conditions at the base point and three conditions at the head point. We use matrix notation to present these boundary conditions in the following form: [figure omitted; refer to PDF] where S=0 indicates the beginning point of calculation (the pile's base point), S=L indicates the end point of calculation (the pile's head point), C is the matrix of boundary condition on the beginning point, and D is the matrix of boundary condition on the end point. They are 3×6 matrices. In this study, the following possible pile end conditions were considered.

(1) Free head (allows both displacement and rotation): N=Q=0 , M=0 . The corresponding matrix of boundary condition is [figure omitted; refer to PDF]

(2) In case of bottom end hinged (allows rotation without displacement): u=v=0 , M=0 . The corresponding matrix of boundary condition is [figure omitted; refer to PDF]

(3) In case of bottom end fixed (allows neither displacement nor rotation): u=v=0 , [straight phi]=0 . The corresponding matrix of boundary condition is [figure omitted; refer to PDF]

(4) In case of bottom end partially hinged (allows rotation without vertical displacement): u=0 , Q=0 , and M=0 . The corresponding matrix of boundary condition is [figure omitted; refer to PDF]

(5) Elastic vertical support at the bottom end: N=_Kv ·u·Area , Q=0 , and M=0 . The corresponding matrix of boundary condition is [figure omitted; refer to PDF]

where _Kv is the modulus of vertical compressibility. Area is the cross-sectional area of the pile.

4. Definition of Boundary Value Problem

Next, we impose the boundary conditions (13) at the pile head and base upon the derived new governing differential equations (12) to define a boundary value problem of the following equations: [figure omitted; refer to PDF] To this end, the response of double-row portal stabilizing pile is mathematically idealised as the boundary value problem of (19).

Thus, many numerical methods to solve the ordinary differential equations can be adopted to solve the boundary value problem of (19).

It should be noted that the existence and uniqueness of solution for the boundary value problem of (19) should be mathematically proved. This matter is however outside the scope of the writer's major. According to the physical character of the problem, we can imagine that the solution exists and is unique. The solution can be validated through comparative studies.

5. Method of Solution

5.1. Uniformity Preprocessing

The orders of magnitude of the section internal forces (N,Q,M ) are so much higher than those of the displacements (u,v,[straight phi] ) that numerical solving of the equations may meet singularity difficulty. So for reasons of numerical stability, it is necessary to perform uniformity preprocessing for the order of magnitude of the element in the coefficient matrix K . We multiply the displacements (u,v,[straight phi] ) by E and substitute the original displacement variables by the expressions (Eu,Ev,E[straight phi] ). So, we redefine two variables as follows: [figure omitted; refer to PDF]

Finally, we obtain the following system of ordinary differential equations: [figure omitted; refer to PDF]

This system of six independent differential equations subjected to boundary conditions can be numerically solved for six unknown functions (three forces and three displacements).

5.2. The Runge-Kutta Finite Difference Algorithm

The Runge-Kutta algorithm is commonly used for the solution of the ordinary differential equation of the form dX/ds=F(X,s) . So, it is chosen to solve (21).

5.2.1. Derivation of the Recursion Formula

The following finite difference formula (22) is one format of the Runge-Kutta methods: [figure omitted; refer to PDF]

The Runge-Kutta algorithm of this type (22) is a numerical method of fifth order, where δ is the stepsize of difference and F(X,s)=K-X+p . For convenience of computation, this formulation may be rewritten in the following form: [figure omitted; refer to PDF]

As we can see, the value of function at point (n+1 ) can be determined from the value of function at point (n ), where n=0 represents the beginning point of calculation and n=m represents the end point of calculation.

In the above equation, _Gn and _Hn can be obtained from the following recursion formula: [figure omitted; refer to PDF] where I is the identity matrix.

5.2.2. Determination of the Initial Vector _X0

The initial value is the start point of the recursion formula. Now, we discuss in the following how to obtain the initial vector _X0 by using the recursion formula of (23) and imposing the boundary conditions at pile head and base.

Considering the recursion formula of (23), _Xn can be expressed in terms of _X0 as follows: [figure omitted; refer to PDF]

In the case of n=0 , we have ^D(0) =I, ^F(0) =0 and substitute it into the recursion formula of (23). We get [figure omitted; refer to PDF]

It can be rewritten as follows: [figure omitted; refer to PDF]

Comparing the above two equations, the recursion formula for ^D(n) and ^F(n) is obtained as follows: [figure omitted; refer to PDF]

Now considering the case of boundary point: _Xm =^D(m)_X0 +^F(m) , we substitute the boundary conditions at end point C_X0 =o , D_Xm =o into the above equation. This leads to the equation to solve for _X0 : [figure omitted; refer to PDF]

The above set of linear algebraic equations can be solved for _X0 by using the method of Gaussian elimination with pivot selection. Once _X0 is known, _Xn can be obtained in sequence using the recursion formula of (23).

5.2.3. Coordinate Transformation between Pile and Beam

When solving the problem of double-row portal piles using the recursion formula of (23), due to the direction change of the axes of the pile and the beam at the connection point (shown in Figure 3), first we need to distinguish K- and p in the recursion formula for the two connected segments. And then in order to satisfy the equilibrium of internal forces and maintain the continuity of displacements at the connection point, we need to introduce the so-called connection matrix to the corresponding formula when dealing with _Dm , _Fm , and _Xn , recursively.

Figure 3: Scheme of local coordinate system transformation between pile and connection beam.

[figure omitted; refer to PDF]

The centroidal axes' rotation from pile to beam at connection point means mathematically that the local coordinate system rotates clockwise by β degree at the connection point (shown in Figure 3).

According to the principle of equilibrium of internal forces, we get [figure omitted; refer to PDF] And the bending moment remains unchanged.

According to the principle of vector analysis, we get [figure omitted; refer to PDF] And the rotation displacement remains unchanged.

Combining the above equations, the corresponding transformation matrix for the local coordinate system rotating clockwise by β degree can be expressed as follows: [figure omitted; refer to PDF]

As we can see, the above connection matrix is the so-called orthogonal matrix whose inverse is its transposed matrix.

As shown in Figure 3, the node of difference (no. n ) is the connection node. First we need to distinguish K- and p in the recursion formula for the two connected segments and then insert the connection matrix for transformation from _Xn to _X-n when dealing with the _Dm , _Fm , and _Xn , recursively. The detailed procedure is as follows.

(1) For the calculation of _Xn , [figure omitted; refer to PDF] _X-n =_CnXn =_Cn (_Gn-1Xn-1 +_Hn-1 δ) , where _Cn is the transformation matrix for point n . Thus, _Xn+1 =_GnX-n +_Hn δ . Next, [figure omitted; refer to PDF] where m indicates the end point of calculation.

(2) For the calculation of _Dm : because _Xn =^D(n)_X0 +^F(n) , and [figure omitted; refer to PDF] So, the above procedure also applies to calculation of ^D(m) as follows: [figure omitted; refer to PDF]

(3) For the calculation of _Fm , [figure omitted; refer to PDF]

where m denotes the end point of calculation.

5.2.4. The Solution Flow Process

In short, the proposed solution procedure involves four main steps:

(1) calculating the value of _Gn and _Hn using the given Equation (25);

(2) calculating ^D(m) and ^F(m) using the given recursion formula of (29);

(3) calculating the vector _X0 by solving linear algebraic Equations (30);

(4) calculating _Xn using the given recursion formula of (23).

Because the equations and solution formula are all given in form of matrices, a simple computer program has been written on the platform of MATLAB to run this procedure. At last, we can get the shear, bending-moment, and deflection diagram along the pile.

6. Verification

The practical examples of portal double-row piles used to stabilize an potential landslide (shown in Figure 4) are considered herein to verify the developed numerical calculation techniques. Soil strength parameters used in the stability analysis are from laboratory shear testing on the undisturbed soil samples. The resisting (shear) force required to achieve the desired safety factor and transferred by the pile is estimated to be 2147 kN/m. Before the pile is installed, the slope is approaching limit state and the safety factor can be assumed to be one. The values of c and [straight phi] can be determined based on experiment data and satisfying this limit state condition. Then, the required resisting force 2147 kN/m can be obtained using back analysis. When this additional force is applied to the specified place of the landslide, the safety factor calculated using the limit equilibrium method can achieve the desired value 1.3. The manual digging discrete reinforced concrete piles were designed to be installed at a spacing of 4 m to increase the factor of safety of the whole slope to the required value of 1.3.

Figure 4: Scheme of portal double-row piles used to stabilize a landslide (units: m).

[figure omitted; refer to PDF]

The front row pile has a length L=16.0 m and sectional dimension 1.5×2.0 m . The length of the portion embedded into the sliding surface is 6.0 m. The back row pile has a length L=37.0 m and sectional dimension 1.5×2.8 m . The length of the portion embedded into the sliding surface is 18.0 m. The sectional dimension of the connection beam is 0.8 × 0.8 m. The pile is constructed using C30 concrete (assuming that concrete does not crack during working). A lateral force F=2147 × 4 kN is assumed to act upon the pile segment above the sliding surface. The pressure distribution is considered to have a rectangular shape, as proposed by the Chinese design code (Code for design on retaining structures of railway subgrade no. TB10025-2006).

The conceptual calculation model used to simulate the lateral response of the pile is shown in Figure 5. The boundary condition at the pile base is considered as free head which allows both lateral displacement and rotation.

Figure 5: The conceptual calculation model for portal double-row stabilizing pile.

[figure omitted; refer to PDF]

For simplicity of program editing, the friction force along the pile-soil interface can be neglected. The modulus of horizontal subgrade reaction for the mudstone below the slip surface is listed in Table 1 for the three different mudstone layers. The coefficient of subgrade reaction was determined according to the suggestions for the mudstone in the related Chinese design code and in situ tests. The related Chinese design code provides detailed experiment procedure to determine the modulus of horizontal subgrade reaction.

Table 1: Parameters of material properties.

	Young's modulus(GPa)	Shear modulus(GPa)	Modulus of subgrade reaction(MPa/m)
C30 concrete	30	12	/
Completely decomposed mudstone	/	/	80
Intensely weathered mudstone	/	/	110
Moderately weathered mudstone	/	/	150

Then, the developed method is applied to analyze the lateral response of the double-row pile which is also computed by the FEM program we developed. We employ an elastic beam column element to model the pile and horizontal spring element to represent the reactions of the surrounding soil in the FEM model. Comparisons of shear force, bending moment, and deflection of the pile between boundary value method (BVM) and FEM are presented in Figures 6 and 7. Complete agreement between them can be observed.

Figure 6: Comparison of internal force and deflection of front row pile between the boundary value method and FEM.

[figure omitted; refer to PDF]

Figure 7: Comparison of internal force and deflection of back row pile between the boundary value method and FEM.

[figure omitted; refer to PDF]

Through the above comparative studies, it has been found that the program we developed works very well and can replace the existing numerical methods that have been used to design the portal double-row stabilizing pile.

7. Summary and Conclusions

In this paper, a new numerical uncoupled method for calculating the response of portal double-row stabilizing piles is proposed. The theoretical background and a detailed derivation of the proposed numerical solution scheme are described. The feasibility of the method developed was verified using the comparative case study. The proposed method has more higher modeling and computing efficiency than the FEM and can be an alternative method for analyzing the behavior of portal double-row piles used for slope stabilization.

Acknowledgments

This study was jointly financially supported by grants from the National Natural Science Foundation of China (Grant no. 51008298) and Key Research Program of Chinese Academy of Sciences (no. KZZD-EW-05). The support is gratefully acknowledged. The author also thanks the anonymous referees whose comments helped to improve the work and the presentation of this paper.

References

[1] Y. K. Chow, "Analysis of piles used for slope stabilization," International Journal for Numerical and Analytical Methods in Geomechanics , vol. 20, no. 9, pp. 635-646, 1996.

[2] J. A. Smethurst, W. Powrie, "Monitoring and analysis of the bending behaviour of discrete piles used to stabilise a railway embankment," Geotechnique , vol. 57, no. 8, pp. 663-677, 2007.

[3] S. Hassiotis, J. L. Chameau, M. Gunaratne, "Design method for stabilization of slopes with piles," Journal of Geotechnical and Geoenvironmental Engineering , vol. 123, no. 4, pp. 314-323, 1997.

[4] H. G. Poulos, "Design of reinforcing piles to increase slope stability," Canadian Geotechnical Journal , vol. 32, no. 5, pp. 808-818, 1995.

[5] T. Ito, T. Matsui, W. P. Hong, "Design method for stabilizing piles against landslide--one row of piles," Soils and Foundations , vol. 21, no. 1, pp. 21-37, 1981.

[6] C.-C. Fan, J. H. Long, "Assessment of existing methods for predicting soil response of laterally loaded piles in sand," Computers and Geotechnics , vol. 32, no. 4, pp. 274-289, 2005.

[7] J.-R. Peng, M. Rouainia, B. G. Clarke, "Finite element analysis of laterally loaded fin piles," Computers and Structures , vol. 88, no. 21-22, pp. 1239-1247, 2010.

[8] S. Jeong, B. Kim, J. Won, J. Lee, "Uncoupled analysis of stabilizing piles in weathered slopes," Computers and Geotechnics , vol. 30, no. 8, pp. 671-682, 2003.

[9] H. Cheng, R. Wei-Zhong, "Boundary value method of ordinary differential equation system for analysing interaction between landslide and stabilizing pile," Chinese Journal of Computational Mechanics , vol. 29, no. 3, pp. 421-426, 2012.

[10] H.-W. Ying, Z.-H. Chu, B.-H. Li, X.-W. Liu, "Study on calculation method of retaining structure with double-row piles and its application," Rock and Soil Mechanics , vol. 28, no. 6, pp. 1145-1150, 2007.

[11] C.-Y. Zhou, Z.-Q. Liu, W. Shang, H. Chen, S.-R. Wen, "New mode for calculation of portal double row anti-sliding piles," Rock and Soil Mechanics , vol. 26, no. 3, pp. 441-449, 2005.

[12] M. Yamin, R. Y. Liang, "Limiting equilibrium method for slope/drilled shaft system," International Journal for Numerical and Analytical Methods in Geomechanics , vol. 34, no. 10, pp. 1063-1075, 2010.

[13] W. D. Guo, F. H. Lee, "Load transfer approach for laterally loaded piles," International Journal for Numerical and Analytical Methods in Geomechanics , vol. 25, no. 11, pp. 1101-1129, 2001.

[14] T. Ito, T. Matsui, "Methods to estimate lateral force acting on stabilizing piles," Soils and Foundations , vol. 15, no. 4, pp. 43-59, 1975.