Acessibilidade / Reportar erro

Static Limit Analysis of Reinforced Soil Structures by a Simple Finite Element and Second-Order Cone Programming

Abstract

To discretize reinforced soil structures in plane strain and predict their collapse load, a simple three-node triangular finite element is formulated based on the static theorem of the limit analysis. The element satisfies the equilibrium equations and the mechanical boundary conditions in a weak sense. A modified Mohr-Coulomb yield surface is adopted to describe the reinforced soil behavior from a macromechanics point of view. It is also taken into account the possibility of tension failure of the reinforcement and failure of the reinforcement interface. The stated nonlinear convex optimization problem is cast as second-order cone programming. Numerical examples illustrate the predictive accuracy of the above scheme as well as the efficiency and speed of an interior-point method to reach optimal solutions.

Keywords:
reinforced soil; static theorem; finite element; second-order cone programming

1 INTRODUCTION

In the design of reinforced soil structures, the maximum load to be resisted at the impending collapse must be evaluated. The ability of the methods to accurately estimate ultimate limit states depends on the fulfillment of theoretical requirements derived from continuum mechanics concerning the equilibrium, strain-displacement relations, constitutive behavior and boundary conditions. Under certain restrictions, the theory of plasticity allows the prediction of the collapse load by means of the limit analysis based on the bound theorems (Drucker et al., 1952Drucker, D.C., Prager, W., Greenberg, H.J. (1952). Extended limit design theorems for continuous media. Quarterly Journal of Applied Mathematics 9(4):381-389.). Numerical solutions are most facilitated in this case coupling such theorems with the finite element method, as originally proposed by Lysmer (1970Lysmer, J. (1970). Limit analysis of plane problems in soil mechanics. Journal of the Soil Mechanics and Foundations Division 96(SM4):1311-1334.) and Bottero et al. (1980Bottero, A., Negre, R., Pastor, J., Turgeman, S. (1980). Finite element method and limit analysis for soil mechanics problems. Computer methods in Applied Mechanics and Engineering 22(1):131-149.). This is a straightforward approach for computing lower and upper bounds on the load at incipient collapse which have gained increasing attention (Lyamin and Sloan, 2002Lyamin, A.V., Sloan S.W. (2002). Lower bound limit analysis using non-linear programming. International Journal for Numerical Methods in Engineering 55(5):573-611.; Krabbenhoft and Damkilde, 2003Krabbenhoft, K., Damkilde, L. (2003). A general non-linear optimization algorithm for lower bound limit analysis. International Journal for Numerical Methods in Engineering 56(2):165-184.; Makrodimopoulos and Martin, 2006Makrodimopoulos, A., Martin, C.M. (2006). Lower bound limit analysis of cohesive-frictional materials using second-order cone programming. International Journal for Numerical Methods in Engineering 66(4):604-634.; Yu and Tin-Loi, 2006Yu, X., Tin-Loi, F. (2006). A simple mixed finite element for static limit analysis. Computers & Structures 84(29-30):1906-1917.; Ciria et al., 2008Ciria, H., Peraire, J., Bonet, J. (2008). Mesh adaptive computation of upper and lower bounds in limit analysis. International Journal for Numerical Methods in Engineering 75(8):899-944.; Munõz et al., 2009; Le et al., 2010Le, C.V., Nguyen-Xuan, H., Nguyen-Dang, H. (2010). Upper and lower bound limit analysis of plates using FEM and second-order cone programming. Computers & Structures 88(1-2):65-73.; Bleyer and de Buhan, 2014Bleyer, J., de Buhan, P. (2014). Lower bound static approach for the yield design of thick plates. International Journal for Numerical Methods in Engineering 100(11):814-833.; Nguyen-Thoi et al., 2015Nguyen-Thoi, T., Phung-Van, P., Nguyen-Thoi, M.H., Dang-Trung, H. (2015). An upper-bound limit analysis of Mindlin plates using CS-DSG3 method and second-order cone programming. Journal of Computational and Applied Mathematics 281:32-48.). This is partly due to the development of robust optimization methods and solvers on which they strongly rely (Drud, 1996Drud, A.S. (1996). CONOPT: A system for large scale nonlinear optimization, Reference manual for CONOPT subroutine library. ARKI Consulting and Development A/S (Bagsvaerd).; Sturm, 1999Sturm, J.F. (1999). Using SeDuMi 1.02, a MATLAB toolbox for optimization over symmetric cones. Optimization Methods and Software 11-12:625-653.; Lyamin and Sloan, 2002Lyamin, A.V., Sloan S.W. (2002). Lower bound limit analysis using non-linear programming. International Journal for Numerical Methods in Engineering 55(5):573-611.; Murtagh and Saunders, 2003Murtagh, B.A., Saunders, M.A. (2003). MINOS 5.51 user’s guide. Technical Report Sol 83-20R, Department of Management Science and Engineering, Stanford University (Stanford).; Tütüncü et al., 2003Tütüncü, R.H., Toh, K.C., Todd, M.J. (2003). Solving semidefinite-quadratic-linear programs using SDPT3. Mathematical Programming, Series B, 95:189-217.; Andersen et al., 2003Andersen, E.D., Roos, C., Terlaky, T. (2003). On implementing a primal-dual interior-point method for conic quadratic optimization. Mathematical Programming, Series B, 95(2):249-277.; Krabbenhoft and Damkilde, 2003Krabbenhoft, K., Damkilde, L. (2003). A general non-linear optimization algorithm for lower bound limit analysis. International Journal for Numerical Methods in Engineering 56(2):165-184.).

The identification of failure modes through experiments has had an essential impact on the development of reinforced soil mechanics. It has been realized that, in the case of sufficiently good bonding between soil and reinforcement, reinforced soil behaves, in the macroscopic scale, like a homogeneous material. This important observation was confirmed by laboratory investigations reported in Long et al. (1972Long, N.T., Guégan, Y., Legeay, G. (1972). Étude de la terre armée à l’appareil triaxial. Rapport de Recherches 17, Laboratoire Central des Ponts et Chaussées (Paris).) and Chapuis (1980Chapuis, R.P. (1980). La compression triaxiale des sols armés. Canadian Geotechnical Journal 17(2):153-164.). Reinforced soil can then be treated as a composite material formed by the association of frictional soil and tension-resistant fiber elements, such as geosynthetics, with behavior presumed homogeneous and anisotropic from a macromechanics point of view. The material strength is estimated from the strength characteristics of its components and the interaction between them. Limit analysis procedures derived from the referred macromechanics approach have been successfully applied to predict the observed response of reinforced superficial foundations and earth walls (Sawicki, 1983Sawicki, A. (1983). Plastic limit behavior of reinforced earth. Journal of Geotechnical Engineering 109(7):1000-1005.; Kulczykowski, 1985Kulczykowski, M. (1985). Analysis of bearing capacity of reinforced subsoil loaded by a footing, Ph.D. Thesis (in Polish), Institute of Hydro-Engineering, Poland.; Sawicki and Lesniewska, 1987Sawicki, A., Lesniewska, D. (1987). Failure modes and bearing capacity of reinforced soil retaining walls. Geotextiles and Geomembranes 5(1):29-44., 1989; de Buhan et al., 1989de Buhan, P., Siad, L. (1989). Influence of a soil-strip interface failure condition on the yield strength of reinforced earth. Computers and Geotechnics 7(1-2):3-18.; de Buhan and Siad, 1989; Yu and Sloan, 1997Yu, H.S., Sloan, S.W. (1997). Finite element limit analysis of reinforced soils. Computers & Structures 63(3):567-577.).

In this paper, specific attention is focused on the static approach of the limit analysis coupled with the finite element method to provide limit load solutions for reinforced soil structures in plane strain. In the description of the soil behavior, the isotropic Mohr-Coulomb yield surface is modified to include the effect of anisotropy caused by the presence of reinforcement. It is also taken into account the possibility of tension failure of the reinforcement and failure of the soil-reinforcement interface (de Buhan and Siad, 1989de Buhan, P., Siad, L. (1989). Influence of a soil-strip interface failure condition on the yield strength of reinforced earth. Computers and Geotechnics 7(1-2):3-18.). The finite element proposed for limit analysis, which has its roots in the classic paper by Anderheggen and Knöpfel (1972Anderheggen, E., Knöpfel, H. (1972). Finite element limit analysis using linear programming. International Journal of Solids and Structures 8(12):1413-1431.), satisfies the equilibrium equations and the mechanical boundary conditions on average. In this sense, a solution obtained with this weak equilibrium model is not a strict lower bound. The formulation of the static theorem is then written within the framework of a nonlinear convex optimization technique, known as second-order cone programming (Lobo et al., 1998Lobo, M.S., Vandenberghe, L., Boyd, S., Lebret, H. (1998). Applications of second-order cone programming. Linear Algebra and its Applications 284(1-3):193-228.).

There is no analytical approach for the solution of general convex optimization problems, but there are very effective methods for solving them (Boyd and Vandenberghe, 2009Boyd, S., Vandenberghe, L. (2009). Convex Optimization. Cambridge University Press.). The interior-point method employed in this paper works extremely well in practice. It was developed by Andersen et al. (2003Andersen, E.D., Roos, C., Terlaky, T. (2003). On implementing a primal-dual interior-point method for conic quadratic optimization. Mathematical Programming, Series B, 95(2):249-277.) and implemented in the solver MOSEK (MOSEK ApS, 2016MOSEK ApS (2016). Modeling cookbook. Available from https://mosek.com/resources/doc/.
https://mosek.com/resources/doc/...
) to address nonlinear convex optimization problems, such as second-order cone programs. Predictions of the collapse load for some well-known stability problems of soil mechanics, namely the bearing capacity of footings and the stability of slopes, are computed after the soil is reinforced and the discrete problem is formulated as second-order cone programming (SOCP). For comparison, the examples are also treated as linear programming (LP) and the results with both schemes compared with exact or approximate solutions available in the literature.

2 FAILURE CONDITIONS

The reinforced soil in this paper will be treated as a homogeneous composite material with anisotropic properties. The reinforcement is assumed to be unidirectional with thickness d very small compared to the space h between two reinforcements (d/h ≪ 1), as shown in Figure 1. Three stress tensors are defined at every point in the homogenized continuum: tensor ⎣σt σt τtn⎦T for macrostresses, and tensors σs=σtsσnsτtnsT and σr=σtrσnrτtnrT for microstresses which act on the soil and reinforcement, respectively. The tensor components are related by

σ t = σ t s + d h σ t r = σ t s + σ r σ n = σ n s = σ n r τ t n = τ t n s = τ t n r (1)

where σr=d/hσtr (Yu and Sloan, 1997Yu, H.S., Sloan, S.W. (1997). Finite element limit analysis of reinforced soils. Computers & Structures 63(3):567-577.).

Figure 1:
Stresses on reinforced soil.

The relations between the macrostress components in the orthogonal Cartesian coordinate systems xy e tn are

σ x = σ t cos 2 θ + σ n sin 2 θ 2 τ t n sin θ cos θ σ y = σ t sin 2 θ + σ n cos 2 θ + 2 τ t n sin θ cos θ τ x y = ( σ t σ n ) sin θ cos θ + τ t n ( cos 2 θ sin 2 θ ) (2)

where θ represents the angle between the horizontal axis x and the reinforcement direction t, measured counterclockwise. Substitution of (1) and relations

σ t s = σ x s cos 2 θ + σ y s sin 2 θ + 2 τ x y s sin θ cos θ σ n s = σ x s sin 2 θ + σ y s cos 2 θ 2 τ x y s sin θ cos θ τ t n s = ( σ x s σ y s ) sin θ cos θ + τ x y s ( cos 2 θ sin 2 θ ) (3)

into (2) yields

σ x s = σ x σ r cos 2 θ σ y s = σ y σ r sin 2 θ τ x y s = τ x y σ r sin θ cos θ . (4)

The reinforcement inside the soil is supposed to act merely as tensile load carrying elements (de Buhan et al., 1989de Buhan, P., Mangiavacchi, R., Nova, R., Pellegrini, G., Saleçon, J. (1989). Yield design of reinforced earth walls by a homogenization method. Géotechnique 39(2):189-201.). The bounds

0 σ r σ o (5)

are thereby imposed on σr , where σo = (d/h)σyield is the tensile yield strength σyield of the reinforcement times its volumetric fraction.

The soil mass with cohesion c and angle of internal friction ϕ is initially assumed to obey the Mohr-Coulomb yield surface in plane strain (Muñoz et al., 2009Muñoz, J.J., Bonet, J., Huerta, A., Peraire, J. (2009). Upper and lower bounds in limit analysis: adaptive meshing strategies and discontinuous loading. International Journal for Numerical Methods in Engineering 77(4):471-501.):

F ¯ s = ( σ x s σ y s ) 2 + ( 2 τ x y s ) 2 [ 2 c cos ϕ ( σ x s + σ y s ) sin ϕ ] (6)

This expression is then modified to include the effect of anisotropy caused by the presence of reinforcement using relations (4):

F s = ( σ x σ y σ r cos 2 θ ) 2 + ( 2 τ x y σ r sin 2 θ ) 2 [ 2 c cos ϕ ( σ x + σ y σ r ) sin ϕ ] (7)

The soil-reinforcement interface failure surface is assumed to be described by

F i = | τ t n | c i + σ n tan ϕ i (8)

where ci and ϕi denote the interface cohesion and angle of internal friction, respectively (Yu and Sloan, 1997Yu, H.S., Sloan, S.W. (1997). Finite element limit analysis of reinforced soils. Computers & Structures 63(3):567-577.). In view of (2), the failure surface (8) takes the form

F i = 1 2 | ( σ y σ x ) sin 2 θ + 2 τ x y cos 2 θ | c i + ( σ x sin 2 θ + σ y cos 2 θ τ x y sin 2 θ ) tan ϕ i = 0 (9)

in the system xy.

3 STATIC THEOREM

The limit analysis relies on the assumption of an elastic-perfectly plastic material with a flow rule associated to a convex yield surface, and also on small displacement gradients so that the body does not undergo large deformation at collapse (problem readily stated on the undeformed body geometry). The static approach of the limit analysis requires that the assumed stress field must satisfy the equilibrium equations, the mechanical boundary conditions and the yield criterion everywhere. Under these idealized conditions, the computed limit load is a lower bound on the true collapse load (Chen, 1975Chen, W.F. (1975). Limit Analysis and Soil Plasticity, Elsevier (Amsterdam).).

3.1 Equilibrium Equations

Let Ω be the region occupied by the reinforced soil structure subjected to the body force b = ⌊bx by⌋T expressed in the system xy. The macrostress field σ = ⌊σx σy σxy⌋T must satisfy the equilibrium equations

D σ + b = 0 (10)

throughout the domain Ω, where the differential operator

D = [ x 0 y 0 x x ] (11)

3.2 Mechanical Boundary Conditions

Let Γ be the boundary of the domain Ω, with unit outward normal vector denoted by n = ⌊nx ny⌋T . In addition to (10), the macrostress field must also satisfy the mechanical boundary conditions

t = N σ = t ¯ (12)

on the portion Γt of Γ on which the traction (stress vector) t = ⌊tx ty⌋T is prescribed as being t¯. Matrix

N = [ n x 0 n y 0 n y n x ] (13)

contains the components of n.

3.3 Yield Criterion

To assure that the yield criterion is satisfied it is necessary to impose simultaneously

F s 0 F i 0 0 σ r σ o (14)

As the function Fi defined in (9) has a term in modulus, the second of inequalities (14) splits into two other conditions given by

F i = A j σ x + B j σ y + C j τ x y c i 0 j = 1,2 (15)

where

A 1 = sin 2 θ tan ϕ i 1 2 sin 2 θ B 1 = cos 2 θ tan ϕ i + 1 2 sin 2 θ C 1 = sin 2 θ tan ϕ i + cos 2 θ A 2 = sin 2 θ tan ϕ i + 1 2 sin 2 θ B 2 = cos 2 θ tan ϕ i 1 2 sin 2 θ C 2 = sin 2 θ tan ϕ i cos 2 θ (16)

4 FINITE ELEMENT FORMULATION

Suppose that the reinforced soil is divided into a number of triangular elements and treated as an assembly of them. To apply the above equations to a typical element Ωє shown in Figure 2, the definition of the boundary should be extended to include the traction continuity on the interelement portion Γi:

( t ) + + ( t ) = 0 (17)

where the superscripts “+” and “-” denote the two sides of Γi. Expressions (10), (12) and (17) may now be viewed as the equilibrium within the element and on its boundary Γt and Γi, respectively (Pian and Wu, 2006Pian, T.H.H., Wu, C.C. (2006). Hybrid and Incompatible Finite Element Methods, Chapman (Boca Raton).; Washizu, 1982Washizu, K. (1982). Variational Methods in Elasticity and Plasticity, 3rd edn., Pergamon Press (Oxford).).

Figure 2:
Finite element Ωє with the unit normal vector n on its boundary Γє.

Equations (10), (12) and (17) can be enforced to be satisfied on average over the element and on its boundary by means of

Ω e w T ( D σ + b ) d x d y Γ t w T ( t t ¯ ) d s Γ i w T t d s = 0 (18)

where w = ⎣wx wy⎦T is an arbitrary weight function that is continuous on the element domain and across its interfaces. The last integral when considered jointly with those of the neighborhood elements enforces (17), since it represents one of the terms in ΓiwT[(t)++(t)] ds=0.

In view of the divergence theorem, one writes

Ω e w T ( D σ ) d x d y = Γ e w T t d s Ω e ( D T w ) T σ d x d y (19)

where w is additionally assumed to be differentiable once with respect to x and y. Since Γe= Γu ⋃ Γt⋃ Γi, expression (18)Lysmer, J. (1970). Limit analysis of plane problems in soil mechanics. Journal of the Soil Mechanics and Foundations Division 96(SM4):1311-1334. is then simplified to

Γ u w T t d s + Γ t w T t ¯ d s + Ω e w T b d x d y Ω e ( D T w ) T σ d x d y = 0 (20)

where the portion Γu of the element boundary Γe falls on the reinforced soil boundary with prescribed displacement.

The macrostress field is linearly approximated over the element by

σ = N 1 σ 1 + N 2 σ 2 + N 3 σ 3 (21)

where

σ i = σ x i σ y i τ x y i T N i = 1 2 A ( α i + β i x + γ i y ) i = 1,2,3 (22)

are the macrostress nodal values and the shape functions, respectively, with A standing for the triangle area. To evaluate

α i = x j y k x k y j β i = y j y k γ i = x k x j (23)

from the node coordinates, the indices i, j, k should be permuted in a natural order (i ≠ j ≠ k). The weight function is also taken to vary linearly,

w = N 1 w 1 + N 2 w 2 + N 3 w 3 (24)

with nodal values

w i = w x i w y i T i = 1,2,3 (25)

Substitution of (21) and (24) into (20) yields

{ w 1 w 2 w 3 } T ( F [ G G G ] { σ 1 σ 2 σ 3 } ) = 0 (26)

where

G = 1 6 [ y 2 y 3 0 x 3 x 2 0 x 3 x 2 y 2 y 3 y 3 y 1 0 x 1 x 3 0 x 1 x 3 y 3 y 1 y 1 y 2 0 x 2 x 1 0 x 2 x 1 y 1 y 2 ] F = Γ t [ N 1 N 2 N 3 ] T t ¯ d s + Ω e [ N 1 N 2 N 3 ] T b d x d y (27)

and

N i = [ N i 0 0 N i ] (28)

Because (26) holds for any arbitrary weight function, it follows that

[ G G G ] { σ 1 σ 2 σ 3 } = F (29)

Elimination of the integral over Γu by choosing w null over there introduces zero components of w1, w2 or w3 into (26) which must be accounted for in (29) by removing the respective equations. Thus, the element enforces the equilibrium equations (10) and the mechanical boundary conditions (12) by the discrete equation (29) on average according to (20). Note that the traction continuity (17) is rigorously enforced because the assumed stress field (21) is continuous across the element interfaces. The element development has its roots in the classic paper by Anderheggen and Knöpfel (1972Anderheggen, E., Knöpfel, H. (1972). Finite element limit analysis using linear programming. International Journal of Solids and Structures 8(12):1413-1431.).

5 THE SOCP PROBLEM

In the static limit analysis, the structure equilibrium and the yield criterion, expressed in terms of nodal stresses, are constraints of an optimization problem for the applied load maximization. The optimal solution

λ * = { max λ | l σ ¯ = λ f 1 + f 2 , g 1 ( σ ¯ , σ ¯ r ) 0 , g 2 ( σ ¯ ) 0 , 0 σ ¯ r σ o } (30)

identifies the collapse load, where the applied load has been split into two parts: ⋋f1 which is adjusted during the optimization by means of the load factor ⋋, and f2 which is kept constant. Vectors σ¯ and σ¯r collect the nodal values of σ and σr respectively.

Under the continuity provided by the approximations (21) and (24), the equality constraint

l σ ¯ = λ f 1 + f 2 (31)

which arises from the assembly of (29), represents the discrete equilibrium of the whole reinforced soil structure. The inequalities constraints

g 1 ( σ ¯ , σ ¯ r ) 0 g 2 ( σ ¯ ) 0 0 σ ¯ r σ o (32)

stem from the evaluation at each node of the first, second and third of inequalities (14), respectively, since it is sufficient to enforce them at each node in order that they are satisfied throughout the mesh.

The nonlinear convex optimization problem (30) is not a SOCP problem. To formulate it as such (Lobo et al., 1998Lobo, M.S., Vandenberghe, L., Boyd, S., Lebret, H. (1998). Applications of second-order cone programming. Linear Algebra and its Applications 284(1-3):193-228.), we introduce the auxiliary variables

v = { v 1 v 2 v 3 } = [ 1 1 0 cos 2 θ 0 0 0 2 sin 2 θ 0 sin ϕ sin ϕ 0 sin ϕ 2 c cos ϕ ] { σ x σ y τ x y σ r 1 } (33)

into (7) to state Fs ≤ 0 as the three-dimensional second-order cone

v 1 2 + v 2 2 v 3 v 3 0 (34)

sketched in Figure 3. Now, the problem can be treated as SOCP by just replacing g1(σ¯, σ¯r)0 with

h 1 ( σ ¯ , σ ¯ r , v ¯ ) = 0 h 2 ( v ¯ ) 0 v ¯ 3 0 (35)

where the vectors v¯ and v¯3 collect the nodal values of v and v3. The new constraints (35) stem from the evaluation of (33) and (34) at each node, and all the problem nonlinearity concentrates in the constraint h2 ≤ 0 related to the cone definition.

Figure 3:
Second-order cone (34).

The optimization problem solution is carried out following the steps: (a) set up of the SOCP problem with YALMIP (Löfberg, 2004Löfberg, J. (2004). YALMIP: A toolbox for modeling and optimization in MATLAB. In Proceedings of the IEEE CSS International Symposium on Computer Aided Control System Design:284-289 (Taipei).) in the MATLAB environment: problem (30) with g1 ≤ 0 replaced by (35); (b) solution by MOSEK using the interior-point optimizer developed for nonlinear convex optimization by Andersen et al. (2003Andersen, E.D., Roos, C., Terlaky, T. (2003). On implementing a primal-dual interior-point method for conic quadratic optimization. Mathematical Programming, Series B, 95(2):249-277.).

6 RESULTS

Predictions of the collapse load for some well-known stability problems of soil mechanics, namely the bearing capacity of footings and the stability of slopes, are computed after the soil is reinforced. This paper focuses on the limit analysis based on the static theorem as well as on the efficiency and speed with which the problems are solved as SOCP. For comparison, the piecewise linearization of Fs adopted by Yu and Sloan (1997Yu, H.S., Sloan, S.W. (1997). Finite element limit analysis of reinforced soils. Computers & Structures 63(3):567-577.) is accomplished, using 48 sides in the linearized yield polygon presented therein, so that the problems can also be solved as LP by performing a sequence of optimization started with the Newton barrier method and concluded with the simplex method at FICO® Xpress Optimization Suite. The role played by the simplex method in such a procedure is to guarantee that the optimal solution is achieved. The computations are performed on a Sony VAIOTM All-in-One machine (i5 dual core 2.50 GHz CPU, 12 GB RAM), running a 64-bit Windows 7. The reported CPU times refer to the time spent on the optimization iterations.

6.1 Strip Footing on Cohesionless Reinforced Soil

For the strip footing of width B indicated in Figure 4a, only half of the full geometry needs to be considered due to symmetry conditions. To deal with the semi-infinite domain, a rectangular portion is discretized as in Figure 4b extending horizontally from the symmetry plane to the right over a length 3B and vertically downwards to a Depth 2B. This domain is sufficiently extended so that conditions at the far boundary do not have significant effect on the calculated limit load.

Figure 4:
Strip footing: (a) half of the adopted domain with horizontal reinforcement; (b) unstructured mesh with 2051 elements.

In this example, f1 present in (30), denotes all the loads applied on Γt and f2 = 0 also present in (30), accounts for the null body forces. The adopted mechanical boundary conditions are

t s = 0 on the left ( symmetric ) edge t n = 0 t s = 0 on the upper edge without strip footing t n = q on the upper edge with strip footing (36)

where q is the force per unit area exerted by the footing, tn and ts are the traction components normal and tangential to the edges.

The weight function components wn and ws are imposed to be null at nodes on the edges where the traction components tn and ts are unknown, respectively:

w n = 0 on the left ( symmetric ) edge w s = 0 on the upper edge with strip footing (37)

In addition, wn and ws are made to be null on the far right and bottom edges because tn and ts are taken to be unknown over there.

It is well-known that the bearing capacity of strip footing is null when the soil is purely frictional, weightless and subjected to no surcharge (Terzaghi, 1943Terzaghi, K. (1943). Theoretical Soil Mechanics, John Wiley Sons (New York).; Davis and Booker, 1971Davis, E.H., Booker, J.R. (1971). The bearing capacity of strip footings from the standpoint of plasticity theory. In Proceedings of the first Australian-New Zealand Conference on Geomechanics:276-282 (Melbourne).). The simple inclusion of horizontal reinforcement, however, rises the bearing capacity to

q e σ o = ( 1 + sin ϕ ) e ( π 2 + ϕ ) tan ϕ (38)

This is the exact solution, as demonstrated by Kulczykowski (1985Kulczykowski, M. (1985). Analysis of bearing capacity of reinforced subsoil loaded by a footing, Ph.D. Thesis (in Polish), Institute of Hydro-Engineering, Poland.), when the reinforced soil is treated as a homogeneous and anisotropic composite material with perfectly rough soil-reinforcement interface. Table 1 summarizes our results obtained by simulating the above property using ci = c = 0, ϕ i = ϕ. It is also depicted the number of iterations and the CPU time required to solve the SOCP and the LP problems. The small difference between the SOCP and LP predictions can be reduced increasing the number of sides in the linearized yield polygon. However, the computational cost of an LP solution may become prohibitive for a larger number of sides because of the amount of iterations and CPU time involved. The SOCP number of iterations agrees with already published numerical experiments, in the sense that an interior-point method demands typically between 5 and 50 iterations to solve any SOCP problem (Lobo et al., 1998Lobo, M.S., Vandenberghe, L., Boyd, S., Lebret, H. (1998). Applications of second-order cone programming. Linear Algebra and its Applications 284(1-3):193-228.). In the LP solutions, the simplex phase answers for 99.2% to 99.9% of the total number of iterations. For instance, out of 10103 iterations for ϕ = 10° the count of 10037 takes place in the simplex phase.

Table 1:
Bearing capacity of strip footing on cohesionless reinforced soil.

The exact solution (38) is graphically compared in Figure 5 to our SOCP predictions. The rigorous lower bounds provided by Yu and Sloan (1997Yu, H.S., Sloan, S.W. (1997). Finite element limit analysis of reinforced soils. Computers & Structures 63(3):567-577.) are also shown, which are based on the discretization of the whole domain into triangular elements devised by Lysmer (1970Lysmer, J. (1970). Limit analysis of plane problems in soil mechanics. Journal of the Soil Mechanics and Foundations Division 96(SM4):1311-1334.) along with extension elements developed by Pastor (1978Pastor, J. (1978). Analyse limite: détermination numérique de solutions statistiques complètes. Application au talus vertical. Journal de Mecánique Appliquée 2(2):167-196.). In that work, as stated previously, the optimization problem is conceived as linear programming. For a given mesh, the scheme adopted by Yu and Sloan leads to a much larger problem size because of the use of Lysmer element and the piecewise linearization of Fs. The gap to the exact solution displayed by our predictions is almost steady for all angles of internal friction considered, whereas it is noticeable the deterioration of Yu and Sloan lower bounds for higher values of ϕ.

Figure 5:
Strip footing on cohesionless reinforced soil: q/σo versus ϕ.

6.2 Strip Footing on Cohesive-Frictional Reinforced Soil

The exact bearing capacity of strip footing on weightless unreinforced soil is given by

q e = c [ e π tan ϕ tan 2 ( π 4 + ϕ 2 ) 1 ] cot ϕ (39)

as demonstrated by Prandtl (1920Prandtl, L. (1920). Über die härte plastischer körper. Nachr K Ges Wiss, Göttingen, Mathematish-Physikalische Klasse:74-85.). The advantage of having horizontal reinforcement can be estimated by analyzing numerically the change in (39) due to the reinforcement inclusion. The same mesh of Figure 4b is adopted and the mechanical boundary conditions are taken as

t s = 0 on the left ( symmetric ) edge t n = 0 t s = 0 on the upper edge without strip footing t n = q t s = 0 on the upper edge with strip footing. (40)

The weight function component wn is imposed to be null at nodes on the left edge, and wn and ws are made to be null on the far right and bottom edges as before.

Table 2 presents the ratio qr/qe between the bearing capacities of strip footing on reinforced and unreinforced soil for different σo/c and angles of internal friction. It is assumed perfectly rough soil-reinforcement interface (ci = c, , ϕi = ϕ). Although the SOCP and LP results differ from each other by a small amount, the difference of iterations and CPU time between both optimization techniques is remarkably large. Results for σo/c = 0 illustrate the accuracy of our predictions for the bearing capacity of strip footing on unreinforced soil.

Table 2:
Ratio qr/qe between the bearing capacities of strip footing on cohesive-frictional soil.

Figure 6 shows the bearing capacity ratio qr/qu obtained with SOCP, where qu refers to the unreinforced scenario. For a given ϕ, the advantage of the reinforcement increases with σo/c. Indeed, the ratio qr/qu rises almost linearly with σo/c, where σo can be augmented by either increasing the reinforcement yield strength σyield or the reinforcement volumetric fraction d/h. We can observe that the benefit from horizontal reinforcement is less for soils with higher angle of internal friction.

Figure 6:
Strip footing on cohesive-frictional soil: qr/qu versus σo/c for different angles of internal friction ϕ.

6.3 Cohesionless Reinforced Earth Wall

To deal with the reinforced earth wall, the dimensions depicted in Figure 7a are considered sufficiently large to represent adequately the semi-infinite domain and the mesh adopted in the analysis is shown in Figure 7b. In this example, f1, present in (30), represents the vertical body force and f2 = 0, also present in (30), accounts for the null loads acting on Γt. The mechanical boundary conditions are

t n = 0 t s = 0 on the three free edges. (41)

Figure 7:
Earth wall: (a) adopted domain with horizontal reinforcement; (b) unstructured mesh with 4147 elements.

On the far left, right and bottom edges, where tn and ts are supposed to be unknown, wn and ws are set equal to zero.

Table 3 condenses our results for horizontal reinforcement under perfectly rough soil-reinforcement interface (ci = c = 0, ϕi = ϕ). Note again the small difference between SOCP and LP predictions, but a much higher computational cost of an LP solution. Figure 8 shows that our results obtained with SOCP are bracketed by the lower bounds of Yu and Sloan (1997Yu, H.S., Sloan, S.W. (1997). Finite element limit analysis of reinforced soils. Computers & Structures 63(3):567-577.) and the upper bounds

γ H u p p e r σ o = 2 tan 2 ( π 4 + ϕ 2 ) (42)

of Sawicki and Lesniewska (1987Sawicki, A., Lesniewska, D. (1987). Failure modes and bearing capacity of reinforced soil retaining walls. Geotextiles and Geomembranes 5(1):29-44.).

Table 3:
Critical height of cohesionless reinforced earth wall.

Figure 8:
Cohesionless reinforced earth wall: γH/σo versus ϕ.

6.4 Cohesive-Frictional Reinforced Earth Wall

To estimate the advantage of having horizontal reinforcement in cohesive-frictional earth walls, analyses are performed for cases with and without reinforcement. The same mesh and boundary conditions of the previous example are considered. Assuming perfectly rough soil-reinforcement interface (ci = c, ϕi = ϕ), Table 4 presents the ratio Hr/Hc between the critical heights of reinforced and unreinforced earth wall for different σo/c and angles of internal friction. The values Hc are elaborate upper bounds computed by Chen (1975Chen, W.F. (1975). Limit Analysis and Soil Plasticity, Elsevier (Amsterdam).) for the unreinforced soil using a logarithmic spiral surface discontinuity mechanism passing through the vertical slope toe. Although the difference between the SOCP and LP solutions is small, the computational cost of one technique differs excessively from the other. Results for σo/c = 0 illustrate that our solutions and Chen upper bounds are in close agreement.

Table 4:
Ratio Hr/Hc between the critical heights of cohesive-frictional earth wall.

Figure 9 depicts the ratio Hr/Hu predicted with SOCP, where Hu refers to the critical height of unreinforced soil. The ratio increases almost linearly with σo/c for a given ϕ, similar to what is observed in Figure 6 for strip footing. However, the benefit of reinforcement inclusion in earth wall is greater for soils with higher angle of internal friction and shows to be more dependent on that angle.

Figure 9:
Cohesive-frictional earth wall: Hr/Hu versus σo/c for different angles of internal friction ϕ.

6.5 Cohesionless Reinforced Earth Wall Under Surcharge

Suppose the earth wall treated in the third example is now replaced by one which is weightless and subjected to a vertical load p (force per unit area) applied to the upper edge. The load is uniformly distributed over a length H measured from the vertical edge. Sawicki and Lesniewska (1987Sawicki, A., Lesniewska, D. (1987). Failure modes and bearing capacity of reinforced soil retaining walls. Geotextiles and Geomembranes 5(1):29-44.) managed to identify

p e σ o = tan 2 ( π 4 + ϕ 2 ) (43)

as the exact solution for the collapse load of this problem. Our analysis is carried out using the mesh of Figure 7b and the mechanical boundary conditions

t n = p t s = 0 on the upper edge with surcharge t n = 0 t s = 0 on the upper edge without surcharge and on the two free edges. (44)

The weight function components wn and ws are made to be null on the far left, right and bottom edges as before.

Table 5 outlines our results normalized to the exact solution (43), where excellent accuracy is observed for all angles of internal friction considered. Any attempt to solve the problem as LP has faced significant slowness in the simplex phase (>3600s).

Table 5:
Bearing capacity of cohesionless reinforced earth wall under surcharge.

7 CONCLUSIONS

The paper proposes a numerical approach for computing static limit loads of reinforced soil structures by combining discretization into simple triangular elements and SOCP. Examples illustrate that the above scheme provides excellent results and the fact of not satisfying the equilibrium equations and the mechanical boundary conditions rigorously is far from being a severe handicap for the developed element. For comparison, the examples are also treated as LP after a piecewise linearization of the reinforced soil yield condition. All the optimization problems idealized either as SOCP or LP are of large scale in the sense that they involve more than 10000 variables and constraints. The SOCP problems are solved by MOSEK and the LP problems are solved by FICO® Xpress Optimization Suite, both of which are state of the art optimization packages. It is clear from the results that the solution of large scale SOCP problems with the MOSEK interior-point optimizer is highly efficient and fast when compared with the combination of Newton barrier and simplex methods at FICO® Xpress Optimization Suite.

References

  • Anderheggen, E., Knöpfel, H. (1972). Finite element limit analysis using linear programming. International Journal of Solids and Structures 8(12):1413-1431.
  • Andersen, E.D., Roos, C., Terlaky, T. (2003). On implementing a primal-dual interior-point method for conic quadratic optimization. Mathematical Programming, Series B, 95(2):249-277.
  • Bleyer, J., de Buhan, P. (2014). Lower bound static approach for the yield design of thick plates. International Journal for Numerical Methods in Engineering 100(11):814-833.
  • Bottero, A., Negre, R., Pastor, J., Turgeman, S. (1980). Finite element method and limit analysis for soil mechanics problems. Computer methods in Applied Mechanics and Engineering 22(1):131-149.
  • Boyd, S., Vandenberghe, L. (2009). Convex Optimization. Cambridge University Press.
  • Chapuis, R.P. (1980). La compression triaxiale des sols armés. Canadian Geotechnical Journal 17(2):153-164.
  • Chen, W.F. (1975). Limit Analysis and Soil Plasticity, Elsevier (Amsterdam).
  • Ciria, H., Peraire, J., Bonet, J. (2008). Mesh adaptive computation of upper and lower bounds in limit analysis. International Journal for Numerical Methods in Engineering 75(8):899-944.
  • Davis, E.H., Booker, J.R. (1971). The bearing capacity of strip footings from the standpoint of plasticity theory. In Proceedings of the first Australian-New Zealand Conference on Geomechanics:276-282 (Melbourne).
  • de Buhan, P., Mangiavacchi, R., Nova, R., Pellegrini, G., Saleçon, J. (1989). Yield design of reinforced earth walls by a homogenization method. Géotechnique 39(2):189-201.
  • de Buhan, P., Siad, L. (1989). Influence of a soil-strip interface failure condition on the yield strength of reinforced earth. Computers and Geotechnics 7(1-2):3-18.
  • Drucker, D.C., Prager, W., Greenberg, H.J. (1952). Extended limit design theorems for continuous media. Quarterly Journal of Applied Mathematics 9(4):381-389.
  • Drud, A.S. (1996). CONOPT: A system for large scale nonlinear optimization, Reference manual for CONOPT subroutine library. ARKI Consulting and Development A/S (Bagsvaerd).
  • Fair Isaac Corporation (2013). FICO® Xpress-Optimizer reference manual.
  • Krabbenhoft, K., Damkilde, L. (2003). A general non-linear optimization algorithm for lower bound limit analysis. International Journal for Numerical Methods in Engineering 56(2):165-184.
  • Kulczykowski, M. (1985). Analysis of bearing capacity of reinforced subsoil loaded by a footing, Ph.D. Thesis (in Polish), Institute of Hydro-Engineering, Poland.
  • Le, C.V., Nguyen-Xuan, H., Nguyen-Dang, H. (2010). Upper and lower bound limit analysis of plates using FEM and second-order cone programming. Computers & Structures 88(1-2):65-73.
  • Lobo, M.S., Vandenberghe, L., Boyd, S., Lebret, H. (1998). Applications of second-order cone programming. Linear Algebra and its Applications 284(1-3):193-228.
  • Löfberg, J. (2004). YALMIP: A toolbox for modeling and optimization in MATLAB. In Proceedings of the IEEE CSS International Symposium on Computer Aided Control System Design:284-289 (Taipei).
  • Long, N.T., Guégan, Y., Legeay, G. (1972). Étude de la terre armée à l’appareil triaxial. Rapport de Recherches 17, Laboratoire Central des Ponts et Chaussées (Paris).
  • Lyamin, A.V., Sloan S.W. (2002). Lower bound limit analysis using non-linear programming. International Journal for Numerical Methods in Engineering 55(5):573-611.
  • Lysmer, J. (1970). Limit analysis of plane problems in soil mechanics. Journal of the Soil Mechanics and Foundations Division 96(SM4):1311-1334.
  • Makrodimopoulos, A., Martin, C.M. (2006). Lower bound limit analysis of cohesive-frictional materials using second-order cone programming. International Journal for Numerical Methods in Engineering 66(4):604-634.
  • MOSEK ApS (2016). Modeling cookbook. Available from https://mosek.com/resources/doc/
    » https://mosek.com/resources/doc/
  • Muñoz, J.J., Bonet, J., Huerta, A., Peraire, J. (2009). Upper and lower bounds in limit analysis: adaptive meshing strategies and discontinuous loading. International Journal for Numerical Methods in Engineering 77(4):471-501.
  • Murtagh, B.A., Saunders, M.A. (2003). MINOS 5.51 user’s guide. Technical Report Sol 83-20R, Department of Management Science and Engineering, Stanford University (Stanford).
  • Nguyen-Thoi, T., Phung-Van, P., Nguyen-Thoi, M.H., Dang-Trung, H. (2015). An upper-bound limit analysis of Mindlin plates using CS-DSG3 method and second-order cone programming. Journal of Computational and Applied Mathematics 281:32-48.
  • Pastor, J. (1978). Analyse limite: détermination numérique de solutions statistiques complètes. Application au talus vertical. Journal de Mecánique Appliquée 2(2):167-196.
  • Pian, T.H.H., Wu, C.C. (2006). Hybrid and Incompatible Finite Element Methods, Chapman (Boca Raton).
  • Prandtl, L. (1920). Über die härte plastischer körper. Nachr K Ges Wiss, Göttingen, Mathematish-Physikalische Klasse:74-85.
  • Sawicki, A. (1983). Plastic limit behavior of reinforced earth. Journal of Geotechnical Engineering 109(7):1000-1005.
  • Sawicki, A., Lesniewska, D. (1987). Failure modes and bearing capacity of reinforced soil retaining walls. Geotextiles and Geomembranes 5(1):29-44.
  • Sawicki, A., Lesniewska, D. (1989). Limit analysis of cohesive slopes reinforced with geotextiles. Computers and Geotechnics 7(1-2):53-66.
  • Sturm, J.F. (1999). Using SeDuMi 1.02, a MATLAB toolbox for optimization over symmetric cones. Optimization Methods and Software 11-12:625-653.
  • Terzaghi, K. (1943). Theoretical Soil Mechanics, John Wiley Sons (New York).
  • Tütüncü, R.H., Toh, K.C., Todd, M.J. (2003). Solving semidefinite-quadratic-linear programs using SDPT3. Mathematical Programming, Series B, 95:189-217.
  • Washizu, K. (1982). Variational Methods in Elasticity and Plasticity, 3rd edn., Pergamon Press (Oxford).
  • Yu, H.S., Sloan, S.W. (1997). Finite element limit analysis of reinforced soils. Computers & Structures 63(3):567-577.
  • Yu, X., Tin-Loi, F. (2006). A simple mixed finite element for static limit analysis. Computers & Structures 84(29-30):1906-1917.

Publication Dates

  • Publication in this collection
    Dec 2017

History

  • Received
    05 Feb 2017
  • Reviewed
    21 July 2017
  • Accepted
    19 Sept 2017
Individual owner www.lajss.org - São Paulo - SP - Brazil
E-mail: lajsssecretary@gmsie.usp.br