A bipotential approach for plastic limit loads of strip footings with non-associated materials

International Journal of Non-Linear Mechanics 90 (2017) 1–10

Contents lists available at ScienceDirect

International Journal of Non–Linear Mechanics journal homepage: www.elsevier.com/locate/nlm

A bipotential approach for plastic limit loads of strip footings with nonassociated materials

MARK

⁎

Madani Hamlaoui, Abdelbacet Oueslati , Géry De Saxcé Univ. Lille, CNRS, Arts et Métiers Paris Tech, Centrale Lille, FRE 3723 – LML – Laboratoire de Mécanique de Lille, F-59000 Lille, France

A R T I C L E I N F O

A B S T R A C T

Keywords: Non-associated plasticity Bipotential Mohr–Coulomb criterion Collapse load

This paper is concerned with a bipotential approach for estimating the plastic collapse loads of a half-space made with a non-associated Mohr–Coulomb material and indented by a rigid punch. In geotechnics, this problem is called the bearing capacity of shallow strip footing for which the analytical solution is derived by Prandtl (1920) [46] and Hill (1950) [35] in the context of associated plasticity. However, when the plastic model is not associated, no analytical methods have yet been developed. Here we explore this issue in a rigorous mathematical framework coupling the bipotential concept and limit analysis. First, the method proposed makes use of the method of characteristics to build a statically and plastically admissible stress field that enables a lower estimate of the plastic limit loads. Next, the extended kinematic theorem of limit analysis to non-standard plasticity is applied to derive an upper quasi-bound of the collapse loads. For this aim, the internal rate of plastic dissipation is obtained thanks to the bipotential functional depending on both a trial stress field and a Prandtllike collapse mechanism. The analytic estimates are compared to the formulae and numerical results provided in literature.

1. Introduction Limit analysis [15,16,50,20] is a powerful method for the direct determination of the collapse load of structures subjected to proportional loadings and operating beyond the elastic limit. The constitutive laws are supposed rigid-perfectly plastic, modelled by a plastic domain and an associated plastic flow rule. Typically, the task is to predict the ultimate load factors using the lower and upper theorems related to static and kinematic approaches respectively. We recall that the static method is based on a trial stress distribution in internal equilibrium with the applied loading and satisfies the traction boundary conditions and the plastic criterion. The minimum multiple of this distribution that will cause the solid to collapse provides a lower bound of the ultimate load. The dual (kinematic) approach is based on a failure mechanism that involves a kinematically and plastically admissible velocity vector. The rate of change of internal energy balances the power of external loads which assess the upper bound of the plastic limit load. The exact collapse load is obtained when the static and kinematic approaches provide the same value of the limit load. For a long time, limit analysis has been used to evaluate the bearing capacity of foundations and the slop stability for geomaterials modelled by associated plastic laws. In particular, a problem of practical interest

⁎

for engineers is the bearing capacity of a semi-infinite soil foundation with the standard law of Mohr–Coulomb, a cohesion c and a friction angle φ. Prandtl–Hill analytical solutions provide the exact limit load [46,35,15]. Another development in computing the bearing capacity of foundations is achieved by using the finite element method and the finite difference method [29,30,40,52,27,47]. Note also that numerical limit analysis bounds involving linear/non-linear programming are proposed in the literature [34,49,36]. A notable advantage of these numerical methods is the study of three-dimensional problems with complex geometries and loadings. It is noteworthy that the classical limit analysis theorems are restricted to standard materials with associated flow rule (the plastic strain rate is normal to the yielding surface). However, many experimental observations showed that for geomaterials and polymers, the dilatancy angle is lower than the friction one and thus the plasticity is not associated. Decidedly this affects the failure mechanism and the plastic limit loads. The assessment of the closed-form expression of the bearing capacity of a strip footing remains open. Many numerical results are proposed in the literature [16,52,39,33]. Moreover, a simple but widely used formula for computing the limit load has been proposed by Drescher and Detournay [25]. The classical approach for modelling the non-associated constitutive

Corresponding author. E-mail address: [email protected] (A. Oueslati).

http://dx.doi.org/10.1016/j.ijnonlinmec.2016.12.001 Received 10 December 2015; Received in revised form 30 November 2016; Accepted 1 December 2016 Available online 25 December 2016 0020-7462/ © 2017 Elsevier Ltd. All rights reserved.

International Journal of Non-Linear Mechanics 90 (2017) 1–10

M. Hamlaoui et al.

rule is not relevant for many non-linear materials such as geomaterials and polymers [3,15], cyclic plasticity of metals [4], frictional contact [6] and plasticity with damage [38]. The standard way to model the non-associated plasticity is based on the use of the yield function f to define the elastic domain

laws relies on a dissipation plastic potential different from the plastic yield. An alternative theory is provided by the concept of bipotential introduced by de Saxcé [21] and defining the framework of the Implicit Standard Materials. The application of this concept permits the generalization of the classical variational principles to the boundary value problems with non-standards materials. This formulation combined with the incremental finite element method was applied in [1,2] to study the vertical bearing capacity of the shallow footing with the nonassociated Drucker–Prager criterion. Concerning the theory of limit analysis, the rigorous extension of the upper and lower theorems the Implicit Standard Materials has been proposed in [21,23]. The present work provides analytical estimates of the collapse loads of strip footing with the non-associated Mohr–Coulomb model. The paper is organized as follows. In Section 2, a brief review of the basic notions of the bipotential concept and the extension of the limit analysis tools to the non-associated plasticity are presented. In Section 3, we set out the non-standard Mohr–Coulomb and the plastic flow rule. We prove that this criterion belongs to the Implicit Standard Materials class and we build its bipotential functional. Next, in Section 4, we develop the analytic computations for estimating the collapse load by using the slip-line method and the bipotential tools. Section 5 presents a comparison between the analytic estimates of the bearing capacity of the strip against results and numerical computations provided in the literature. The paper concludes by listing the major contributions of this study and some perspectives for future work.

K = {π such that f (π ) ≤ 0} and on Melan's plastic potential g to define the non-associated flow rule:

∃ λ ≥ 0 such that κ˙ = λ

Unfortunately, although widespread in the literature, this formalism is not a convenient framework to develop variational approaches based on functionals. The first attempt to address this challenge was the hemivariational inequality approach proposed by Panagiatopoulos [45]. Later, de Saxcé [24,21] developed the bipotential theory based on an extension of Fenchel's inequality [28], and the generalization of Moreau's superpotential [43] for the non-associated models. By definition, a bipotential is a functional b : X × Y → ] − ∞, +∞], defined by the following properties: (a) b is convex and lower semicontinuous in each argument. (b) For any κ˙′ and π′ we have

(c) For κ˙ and π we have the equivalences:

π ∈ ∂b (·,π )(κ˙)⟺κ˙ ∈ ∂b (κ˙, ·)(π )⟺b (κ˙, π ) = κ˙·π

Let X be a topological vector space of velocities κ˙ , and Y be its dual space of like-stress variables π, put in duality through a dual pairing X × Y →  : (κ˙, π ) ↦ κ˙·π . In convex analysis, the sub-differential of a function φ in a point κ˙ is the (possibly empty) set:

κ˙ =

(1)

If φ is a smooth and convex function, then the law is uni-valued and we get

∂b (κ˙, π ), ∂π

π=

∂b (κ˙, π ) ∂κ˙

(8)

Relations (7) and (8) can be qualified by implicit normality laws in the sense that (by reference to the implicit function theorems) the unknown κ˙ (resp. π) belongs to both left and right hand members of the relation. For this reason, the dissipative materials admitting a bipotential are called Implicit Standard Materials. Generally speaking, these materials are non-Druckerian because the constitutive law is non-associated. This behaviour is unstable and softening may occur [7]. In the particular case of standard materials, it is easy to show that the bipotential is separated in two parts: the classical superpotential of dissipation φ and its conjugate function φ*:

∂φ (κ˙) = {Dφ (κ˙)} For more details on convex analysis, the reader is referred for instance to [42,48]. On this basis, the concept of potential can be extended in a weak form. We do not require more than the function φ to be convex and lower semi-continuous (with possible infinite values) and we consider multivalued laws generated by φ according to (2)

The function φ is called a superpotential. The converse law takes a similar form:

κ˙ ∈ ∂φ* (π )

(7)

From a mechanical point of view, the bipotential represents the plastic dissipation power (by volume unit) and the two former conditions in (7) are the constitutive law and its inverse one. The couples (κ˙, π ) satisfying the latter conditions in (7) are called extremal couples. If the bipotential is differentiable, the two fist relations in (7) reads

2.1. Definition and basic relations

π ∈ ∂φ (κ˙)

(6)

b (κ˙′, π ′) ≥ κ˙′·π ′

2. The bipotential concept

∂φ (κ˙) = {π ∈ Y | ∀ κ˙′ ∈ X , φ (κ˙′) − φ (κ˙) ≥ (κ˙′ − κ˙)· π}

∂g (π ) ∂π

b (κ˙, π ) = φ (κ˙) + φ* (π )

(9)

Thus the fundamental inequality (6) degenerates into Fenchel's one (4).

(3)

where φ* is Fenchel's transform (or conjugate) of φ

2.2. Limit analysis of implicit standard materials

φ* (π ) = sup(κ˙·π − φ (κ˙))

Classical theorems of limit analysis require that the structures be rigid-perfectly plastic obeying to the normal flow rule. Fortunately, the bipotential approach paves the way to variational formulation which leads to the extension of the static and kinematic theorems of limit analysis to non-associated materials [22,23,6,53,13,17] as we shall see hereafter. Consider a rigid-perfectly plastic solid occupying the volume V with a regular boundary S. This surface S is split into two disjoint parts St and Sv such that St ∪ Sv = S . The solid is subjected to volume forces F , a surface traction distribution t on St and a velocity field v on Sv. A velocity field v′ is said to be kinematically admissible (K.A.) if v′ = v on Sv and d (v′) = grads (v′) in V. A stress field σ′ is statically admissible (S.A.) if div σ ′ + F = 0 within V and t (σ ′) = σ ′·n = t on St. σ′ is said to be Plastically Admissible (P.A.) if σ ′ ∈ K where K is the

κ˙∈ X

As φ is convex and lower semicontinuous, φ* so is and φ** = φ . Consequently, the superpotential φ and its Fenchel's conjugate φ* satisfy Fenchel's inequality [28]

∀ κ˙′ ∈ X , π ′ ∈ Y ,

φ (κ˙′) + φ* (π ′) ≥ κ˙′·π ′

(4)

Moreover, (2) and (3) are equivalent to

φ (κ˙) + φ* (π ) = κ˙·π

(5)

For instance, the associated plasticity is obtained by taking φ* as the indicator function IK of the convex strength domain K (equal to zero in K and +∞ otherwise), and by considering the normal flow rule (3). Nevertheless, it has been recognized experimentally that the normality 2

International Journal of Non-Linear Mechanics 90 (2017) 1–10

M. Hamlaoui et al.

αL ≤ αk

strength convex domain. The boundary value problem (B.V.P.) corresponding to given external loadings (F , t , v ) is the following: find (v, σ ) such that v is K.A., σ is S.A. and d (v ), σ are related by the non-associated plastic law

This result constitutes the extended kinematic theorem. The proof can be found for instance in [21,23,13].

(10)

σ ∈ ∂b (·,σ )(d (v ))⟺d (v ) ∈ ∂b (d (v ), ·)(σ ),

(18)

3. The non-associated Mohr–Coulomb model

which is also equivalent to the extremality condition The behaviour of geomaterials and cohesive-frictional soils is usually modelled by the Mohr–Coulomb criterion with a cohesion c and an internal friction angle φ. The yielding condition is expressed by:

(11)

b (d (v ), σ ) = d (v ): σ

where the dual pairing is the double contracted product denoted by “:”. The variational formalism is typically used in limit analysis for numerical simulations by the finite element method. The basic approaches are based on one-field principles, namely Markov's velocity principle [41] and Hill's stress principle [35]. More accurate numerical approximations can be obtained using two-field principles. For instance, a hybrid-displacement principle was proposed by Nguyen [44] and a saddle-point problem was introduced by Christiansen [14] or Casciaro and co-workers [11,12]. A survey of variational principles in limit analysis is proposed in [50]. For the non-associated plasticity, the corresponding variational formulation is obtained thanks to the bifunctional which belongs to this large class of two-field functionals. The necessity to consider two fields stems here from the non-associativity of the constitutive law. This concept allows us to extend the calculus of variation for material admitting a bipotential [21,8]. The bifunctional is given by:

B (v , σ ) =

∫V b (d (v), σ ) dV − ∫V v·F

dV −

∫S v·t

dS −

t

∫S

τ + σ tan(φ) = c

(19)

where σ is the normal stress and τ is the shear strength on any surface element at failure. In the present study, the traction stress is taken as positive. In the (σ , τ ) space, the yield curve is composed of two edges of the cone Kσ defined by

Kσ = {(σ , τ ) such that τ + σ tan(φ) − c ≤ 0}

(20)

⎞ ⎛ c with vertex at ⎜ tan(φ) , 0⎟. For sake of simplicity, the vertex position is ⎠ ⎝ c denoted by H = tan(φ) in the sequel. Let ε˙ the normal strain velocity and γ˙ the shear angle velocity, corresponding respectively to σ and τ. At any regular point of the yield curve, the plastic flow rule is characterized by a constant plastic dilatancy angle ψ such that

t (σ )·v dS

v

ε˙ = ± γ˙ tan(ψ )

(21)

(12) where the plus (resp. minus) sign corresponds to the stress state on the upper (resp. lower) edge of the cone Kσ . The angle ψ belongs to the range from zero to φ. The event φ = ψ corresponds to the associated material. Otherwise, the material is said to be non-associated. The flow rule claims that the couple (ε˙, γ˙) belongs to the cone  ϵ defined as

It satisfies the following relations:

B (v′, σ ′) ≥ 0 for any v′ K.A. and any σ′ S.A. B (v, σ ) = 0 for the couple (v, σ ) solution of the (B.V.P).

(i) (ii)

 ε = {(ε˙, γ˙) such that γ˙ tan(ψ ) ≤ ε˙}

The couple (v, σ ) solution of the (B.V.P) is the solution of the following variational problems:

B (v, σ ) = min B (v′, σ ) = min B (v, σ′) = 0 v′ K.A.

and to the edges when (σ , τ ) is not the vertex of  σ as shown in Fig. 1. The flow rule may be written in a compact form as follows:

(13)

σ′ S.A.

(ε˙ + (tan(φ) − tan(ψ ))|γ˙|, γ˙) ∈ ∂IKσ

Let us now specify the variational versions of the limit analysis formalism [50] by considering a vanishing prescribed velocity field on Sv and proportional loading (αF 0, αt 0 ). The positive coefficient α is called the load factor and (F 0, t 0 ) represent the reference load distributions. For any solution (v, σ ) of the classical boundary value problem (B.V.P.) under the proportional loading, the variational solution of the mechanical problem satisfies:

t

⎧0 if x ∈ Kσ IKσ (x ) = ⎨ ⎩+ ∞ otherwise

∫V b (v, σ ) dV

(14)

(15)

and the external power 0 Pext (v ) =

∫V v·F 0 dV + ∫S v·t 0 dS t

(16)

The kinematical factor associated to a K.A. velocity field v′ satisfying 0 Pext (v ) > 0 is defined by:

αk =

Pint (v′, σ ) 0 Pext (v′)

(24)

Hereafter, it will be shown that the non-standard Mohr–Coulomb model belongs to the class of Implicit Standard Materials and admits a bipotential given by:

where αL is called a limit load factor. Let us introduce the internal power

Pint (v, σ ) =

(23)

where ∂ is the subgradient operator and IKσ is the indicator function of the convex set Kσ . Recall that the indicator function is defined by

⎞

⎛

∫V b (v, σ ) dV = ∫V d (v): σdV = α L ⎜⎝∫V v·F 0 dV + ∫S v·t 0 dS⎟⎠

(22)

(17)

It depends on the stress field σ, the exact solution of the (B.V.P.). In this case it bounds limit factors Fig. 1. The non-associated Mohr–Coulomb yield condition.

3

International Journal of Non-Linear Mechanics 90 (2017) 1–10

M. Hamlaoui et al.

⎧ Hε˙ + (tan(ψ ) − tan(φ))(σ − H ) |γ˙| ⎪ b ((ε˙, γ˙), (σ , τ )) = ⎨ if (ε˙, γ˙) ∈ Kϵ and (σ , τ ) ∈ Kσ ⎪ ⎩+ ∞ otherwise

σ≤H− (25)

|γ˙| ≤

Proof. To prove this proposition, it is enough to verify that, for any couples (ε˙, γ˙) ∈ Kϵ and (σ , τ ) ∈ Kσ , the following inequality holds: (26)

σε˙ ≤ Hε˙ −

γτ ˙ ≤ |γ˙||τ| ≤

(27)

(29)

Also, since σ ≤ H and using (29) we get

tan(ψ )(σ − H )|γ˙| ≥ (σ − H ) ε˙

|τ | ε˙ tan(φ)

(38)

|τ|ε˙ tan(ψ )

(39)

Addition of (38) and (39) yields to desired inequality (35).□ In most problems of geotechnics and metal forming, the failure involves the development of zones or bands of high shear strain. The deformation is therefore localized within this band of very small thickness h. Moreover, the tangential velocity change δu on a sheared zone is accompanied by the normal velocity change δv such that δv = δu tan(ψ ). This property is called the dilatancy of the material. A physical interpretation of this property can be provided by considering a block shearing on an horizontal plane. Indeed, the block dilates under shearing and causes the volume to change and the displacement (or velocity) to form an angle ψ with the slip band. In addition, since the shear band is thin, then the plastic strain rates can be stated as

(28)

On the other hand, since (γ˙, ε˙) ∈ Kϵ , we have

|γ˙| tan(ψ ) ≤ ε˙

(37)

Moreover, the following relations hold:

Therefore,

−tan(φ) (σ − H )|γ˙| ≥ |τ| |γ˙| ≥ τγ˙

ε˙ tan(ψ )

Subsequently,

Let us first assume that (σ , τ ) ∈ Kσ . Then, it is easy to check that

⎛ c ⎞ σ tan(φ) − c ≤ −|τ| ⇒ −tan(φ) ⎜σ − ⎟ ≥ |τ | tan(φ) ⎠ ⎝

(36)

Next, since (ε˙, γ˙) ∈ Kσ it follows that:

Proposition 1. The function (25) is a bipotential.

Hε˙ + (tan(ψ ) − tan(φ))(σ − H ) |γ˙| ≥ σε˙ + τγ˙

|τ | tan(φ)

(30)

Finally, adding (27) and (30) yields the desired result, (26) is obtained and the proof is achieved.□ Proposition 2. The flow rule (23) is satisfied by the extremal pairs of the bipotential (25) and conversely. Proof. Firstly, it can be easily shown that application of first relation in Eqs. (19) leads to the plastic flow rule (21). Conversely, consider a couple (ε˙, γ˙) satisfying by the flow rule (21) and prove that the bipotential satisfies the extremality condition (11) expressed by:

ε˙ = δv / h ,

(40)

γ˙ = δu / h

In terms of incremental variables, the relation (22) is recast as

Kε = {(δv, δu ) such that δu tan(ψ ) ≤ δv}

(41)

Before it should be pointed out that by performing simple algebraic computations, the finite part of the bipotential given by (25) may be written as follows:

The evaluation of the rate of dissipation of internal energy is crucial in limit analysis computations. For the associated Mohr–Coulomb model (φ = ψ ), it is well known that the rate dissipation along a discontinuity band between tow sliding blocks is independent of the thickness of the band and reads [15,16]:

b ((ε˙, γ˙), (σ , τ )) = H (ε − |γ˙| tan(ψ )) + σ |γ˙| tan(ψ ) + |γ˙|(H − σ )tan(φ)

D = c δu

Hε˙ + (tan(ψ ) − tan(φ))(σ − H ) |γ˙| = σ ε˙ + τγ˙

(31)

where c is the cohesion. For non-associated materials (ψ < φ), formula (42) cannot be longer applied and specific tools are required.

(32) First, it is obvious that, if the plastic strains rate (ε˙, γ˙) vanishes, then (32) is satisfied. Now, consider any regular point of the yield curve Kσ . Then we have

τ + σ tan(φ) = c

Proposition 4. The finite value of the bipotential (25) along a velocity surface discontinuity is given by:

ε˙ = |γ˙| tan(ψ )

⎞ ⎛ sin(2φ) b (p , δu ) = ⎜K c + (1 − K ) (p + H ) ⎟ |δu| ⎠ ⎝ 2

It follows that

b ((ε˙, γ˙), (σ , τ )) = εσ ˙ + |τ||γ˙|

(33)

Proof. Let us denote by R the radius of the Mohr circle which is tangent to the cone Kσ and by p its centre (see Fig. 5). It is easy to establish

Proposition 3. The bipotential for the non-associated Mohr–Coulomb constitutive law is not unique and can be represented by the alternative expression:

R = (p + H )sin(φ)

tan(φ) − tan(ψ ) ε˙ ≥ σε˙ + τγ˙ tan(ψ )tan(φ)

and

σ + p = R sin(φ)

(44)

and thus,

σ − H = R sin(φ) − (p + H ) = −(p + H )cos2 (φ)

(45)

Taking into account the relations (40), the finite part of the bipotential writes

(34)

b = Hδv − (tan(ψ ) − tan(φ))(p + H )cos2 (φ)|δu|

Proof. Similarly to the proof of the Proposition 1, it is sufficient to show that, for any (ε˙, γ˙) ∈ Kϵ and (σ , τ ) ∈ Kσ we have

Hε˙ + |τ|

(43)

where p is the radius of the Mohr circle tangent to the cone Kσ and the tan(ψ ) coefficient K = tan(φ) .

In addition, the flow rule dictates that τ and γ˙ have the same sign. Therefore, Eq. (31) is established. Finally, for the vertex, τ vanishes and Eq. (31) is again fulfilled.□

⎧ tan(φ) − tan(ψ ) Hε˙ + |τ| ε˙ ⎪ ⎪ ∼ tan(ψ )tan(φ) b ((ε˙, γ˙), (σ , τ )) = ⎨ ⎪ if (ε˙, γ˙) ∈ Kϵ and (σ , τ ) ∈ Kσ ⎪ ⎩+ ∞ otherwise

(42)

(46)

Moreover, |δu| tan(ψ ) = δv holds on the velocity surface discontinuity, then we obtain:

⎞ ⎛ tan(ψ ) ⎛ tan(ψ ) ⎞ b=⎜ c + tan(φ) ⎜1 − ⎟ (p + H )cos2 (φ) ⎟ |δu| tan(φ) ⎠ ⎝ ⎠ ⎝ tan(φ)

(35)

First, the condition (σ , τ ) ∈ Kσ implies 4

(47)

International Journal of Non-Linear Mechanics 90 (2017) 1–10

M. Hamlaoui et al.

Finally, by introducing the coefficient K = identity tan(φ)cos2 (φ) =

sin(2φ) 2

tan(ψ ) tan(φ)

tion/compression is unique. Regarding the strip footing problem, according to Davis and Booker [19], the range of correct but not unique plastic limit load is insignificant. The finite element computations [39,37,32] agree with this statement.

and by using the

we obtain the desired relation (43).□

Remarks. 1. The bipotential depends on the stress field through p, the radius of the circle of Mohr, tangent to the cone Kσ at failure. 2. For the associated model, K=1, the bipotential (47) writes b = D = c|δu|. We recover the relation (42). 3. For incompressible materials, K=0, the bipotential (47) reduces to sin(2φ) b (p , δu ) = 2 (p + H )|δu|

4.2. Static approach Let (Oxy ) be a Cartesian frame in the plane of symmetry, σ1 be the major principal stress and σ2 be the minor one. Transversal normal stress is assumed to be the intermediate principal stress. θ, measuring the angle giving the inclination of the principal direction corresponding to σ1, with respect to the (Ox ) axis, results from the usual Mohr circle formulae (Fig. 5):

4. Collapse loads of the strip footing

σ = −p + R cos(2θ ),

4.1. Problem description

(48)

τ = R sin(2θ ) σ +σ

where the position of the circle centre is given by p = − 1 2 2 , and the σ −σ circle radius is R = 1 2 2 . When plastic yielding occurs, the Mohr circle is tangent to the cone which results in:

Consider the problem of punching rigid-perfectly plastic half-space using a rigid stamp subjected to a vertical force. The material is assumed to be weightless and follows the non-associated Mohr– Coulomb model (19)–(21). The study is carried out within the framework of plane strains. Numerical computations worked out by Hamlaoui et al. [32] using the Drucker–Prager model matching the Mohr–Coulomb criterion at plane strain [15] showed that the collapse mechanism is similar to Prandtl's one (4) in the sense that there is a unique wedge below the strip footing, unlike Hill's one. The size and inclination angles of the mechanism, may, however be rather different (see Figs. 2 and 3). On these grounds, a symmetric Prandtl-like mechanism (see Fig. 4) is considered. This mechanism is composed of three zones: (i) a rigid triangular block (OO′A) with base angles α and moving downwards with a constant velocity VP; (ii) a radial shear sector delimited by the lines (OA, OB ); and by a logarithmic spiral arc (AB) and (iii) a rigid wedge (OBC) with bases angles β moving upwards (passive zone). The mass soil below the line (ABC) is at rest. Recall that for the associated φ φ π π Mohr–Coulomb material α = 4 + 2 and β = 4 − 2 . For the non-associated Drucker–Prager soil, Hamlaoui et al. [32] showed that α is independent of the dilatancy angle and can be φ π approximated by 4 + 2 . Therefore, in the present work, we take φ π α = 4 + 2 . Unfortunately, no accurate estimates about the value of the base angle β are available. In this study we propose the trial value ψ π β = 4 − 2 (see Fig. 4). Remark. For non-associated models, the Drucker stability postulate is violated. Thus, the ultimate plastic load may not be unique and may depend on the initial state and/or the path loadings. However, for some problems, analytic proof and numerical verification of the uniqueness of the ultimate load are observed. For example, it has been shown in Cheng et al. [17], using numerical and closed analytical solutions, that the ultimate plastic load of the hollow sphere problem made by a nonassociated Drucker–Prager rigid plastic material is unique, whatever the elastic properties and initial states. Moreover, Boushine et al. [7] proved analytically that the plastic limit load of a non-associated Drucker–Prager specimen in plane strain subjected to uniform trac-

(49)

R = (H + p)sin(φ)

Stresses in the triangular wedge OBC. Consider the uniform stress region OBC (Fig. 6). In the absence of a surcharge stress the major principal stress σ1 = 0 and consequently, Mohr's circle is tangent to the criterion and tangent to the axis τ = 0 at the origin (0, 0). By simple geometric consideration p1 = R1 = (p1 + H )sin(φ) and thus

p1 + H =

H 1 − sin(φ)

(50)

Stress jump in the logspiral shear zone. The logspiral shear zone is bounded by the lines (OA = r0 , OB = r1) and the logspiral curve AB is defined by r (θ ) = r0 e θ tan(ψ ) in polar coordinates (r , θ ). In particular π φ ψ (Fig. 6) OB = r1 = r0 e( 2 − 2 + 2 )tan(ψ ) . We suppose that the stress is continuous through the line OB. The well-known formula of the infinitesimal stress jump formula in the logspiral zone writes [15]:

dθ = −

dp dp 1 1 cos(φ) = − ctan(φ) 2 R 2 p+H

(51)

Therefore, if we integrate over the logspiral zone we find that

∫p

p (θ ) 0

dp = −2 tan(φ) p+H

∫0

θ

dθ

(52)

As a result, the stress distribution in the logspiral shear zone at an arbitrary angle θ is given by

p (θ ) + H = e−2θ tan(φ) p0 + H

(53)

In particular,

p1 + H = e−(π + ψ − φ)tan(φ) p0 + H

(54)

Fig. 2. Collapse mechanism depicted by the isovalues of the vertical displacement increment at the numerical limit state for the Drucker–Prager model matching the Mohr–Coulomb criterion at plane strain, the ruin mechanism for the associate Mohr–Coulomb model is drawn with continued black line [32].

5

International Journal of Non-Linear Mechanics 90 (2017) 1–10

M. Hamlaoui et al.

Fig. 3. Collapse mechanism depicted by the isovalues of the vertical displacement increment at the numerical limit state for the Drucker–Prager model matching the Mohr–Coulomb criterion at plane strain, the ruin mechanism for the associate Mohr–Coulomb model is drawn with continued black line [32].

Principal stress in the wedge OO′A. We suppose that stress is continuous through the line OA (Fig. 6). Using the Mohr circle, we get the principal stresses in the OO′A region. The lower bound of the bearing capacity denoted by Ps is given by: (58)

σ2 = −Ps with

Ps = p0 + R0 = p0 + (p0 + H )sin(φ)

(59)

and thus

Fig. 4. The Prandtl like collapse mechanism for the non-associated material.

Ps + H = (1 + sin(φ))(p0 + H )

(60)

Finally, a lower estimate of the bearing capacity is given by:

⎡ ⎤ ⎛π φ⎞ Ps = H ⎢e(π + ψ − φ)tan(φ) tan2 ⎜ + ⎟ − 1⎥ ⎠ ⎝ 4 2 ⎣ ⎦

(61)

4.3. Kinematic approach Now, we turn our attention to the computation of the upper quasibound. For the compatible Prandtl-like collapse mechanism, the rate of the energy dissipation is a result of the velocity discontinuity on lines (OA, AB, BC ) and the deformation in the logspiral shear zone. It is worth recalling that the velocity is inclined at the angle ψ with lines of velocity discontinuity. Rate of the energy dissipation along the line OA. The rigid wedge (OO′A) under the strip sinks into the soil with a vertical velocity Vp equal to that of the punch. Let V0 be the velocity of the line (OA) and V0p the relative velocity inclined at the angle ψ with (OA). The norm of velocities V0 and V0p can be expressed in terms of the footing velocity Vp from the hodograph, depicted in Fig. 7. The compatibility of velocities yields:

Fig. 5. Mohr diagram at failure.

V0p Vp V0 = = ⎞ ⎛π ⎛π ⎛π φ⎞ φ⎞ sin ⎜ + ψ − ⎟ sin ⎜ + ⎟ sin ⎜ + ψ ⎟ ⎠ ⎝2 ⎝2 ⎝4 2⎠ 2⎠ In addition, it is easy to establish the following relation: Fig. 6. Stresses field in the half-plane.

Thus, we have

p (θ ) + H p (θ ) + H p0 + H = = e(π + ψ − φ −2θ )tan(φ) p1 + H p0 + H p1 + H

(55)

which leads to:

p (θ ) + H = e(π + ψ − φ −2θ )tan(φ)

H 1 − sin(φ)

(56)

It is worthy of note that on the line OA we have

p0 + H = e(π + ψ − φ)tan(φ)

H 1 − sin(φ)

(57)

Fig. 7. Hodograph of velocities.

6

(62)

International Journal of Non-Linear Mechanics 90 (2017) 1–10

M. Hamlaoui et al.

Fig. 9. Velocities in the active wedge OBC.

[15]. Rate of the energy dissipation along the spiral curve AB. It can be easily shown that the rate dissipation of energy, denoted by BAB, on the velocity discontinuity curve (AB ), is equal to Bf:

Fig. 8. Logspiral shear zone.

⎛π α⎞ L cos ⎜ + ⎟ = ⎝4 2 ⎠ 2r0

(63)

where r0 is the length of the line OA. The rate of dissipation on the velocity discontinuity line (OA), denoted by BOA, is simply obtained by multiplying the bipotential by r0:

⎡ ⎤ sin(2φ) BOA = r0 Vop cos(ψ ) ⎢Kc + (1 − K ) (p0 + H ) ⎥ ⎣ ⎦ 2

BAB =

(64)

Therefore,

BOA

⎛π φ⎞ sin ⎜ + ⎟ ⎡ ⎝4 2⎠ = r0 V0 cos(ψ ) ⎢Kc + (1 − K ) ⎛π ⎣ φ⎞ sin ⎜ + ψ − ⎟ ⎝4 2⎠ ⎤ sin(2φ) H e(π + ψ − φ)tan(φ) ⎥ ⎦ 2 1 − sin(φ)

(65)

⎡ ⎤ sin(2φ) H BBC = r0 V0 cos(ψ ) e(π + ψ − φ)tan(ψ ) ⎢Kc + (1 − K ) ⎥ ⎣ 2 1 − sin(φ) ⎦

(66)

Pext = Pk L Vp

Pint = BOA + Bf + BAB + BBC = r0 V0 (1 − sin(ψ )) ⎡ ⎞ ⎛ ⎛π φ⎞ ⎢ ⎟ ⎜ sin(ψ )sin ⎜ + ⎟ ⎝4 2⎠ ψ⎞ ⎢ (π + ψ − φ)tan(ψ ) 2 ⎛ π 1 ⎟ ⎜ e tan + + − 1 ⎟ ⎜ ⎢ ⎟ ⎞ ⎛π ⎠ 1 − sin(ψ ) ⎜ ⎝4 2 φ ⎢ ⎟ ⎜ sin ⎜ + ψ − ⎟ ⎢⎣ ⎝4 2⎠ ⎠ ⎝

⎡ sin(2θ ) He(π + ψ − φ) tan(φ) 2θ (tan(ψ )−tan(φ)) ⎤ + r0 v0 ⎢ (1 − K ) e ⎥ dθ ⎣ ⎦ 2 1 − sin(φ) (67) The energy dissipated in the central shear zone is obtained by integration:

B1 = r0 v0

∫0

(72)

The rate of dissipation of the internal energy is given by

[Kce 2θ tan(ψ ) ] dθ

φ ψ (π − + ) 2 2 2

(71)

An upper estimate of the plastic limit load. The upper bound of the bearing state of the footing is obtained by equating the rate of internal dissipation of energy to the external rate work furnished by the punch load P. The later is given by:

It turns out that

dB1 (θ ) = r0 v0

(70)

Rate of the energy dissipation in the active zone OBC. The velocity is continuous across the line (OB ) so that that the wedge BGF moves as a rigid body with the velocity V1 = V0 e(π + ψ − φ). The vector V1 is orthogonal to BF and inclined at the angle ψ to the line (OB ). The dissipation of internal energy is owing to the discontinuity of velocity on line (BC ) separating the wedge in motion and the material at rest (Fig. 9). The rate of dissipation of energy along line (BC ) is obtained by multiplying the bipotential by the length of BC:

Rate of the energy dissipation in the logspiral shear zone OAB. As in classical limit analysis [15,16], the logspiral shear zone is modelled by a series of rigid elementary triangles, as shown in Fig. 8. Each triangle moves as a rigid body with the hortoradial velocity v (θ ) = V0 e θ tan(ψ ) . The rate of dissipation of energy of an infinitesimal triangle reads:

⎡ ⎤ sin(2θ ) dB1 (θ ) = r (θ ) ⎢Kc + (1 − K ) (p (θ ) + H ) ⎥ v (θ ) dθ ⎣ ⎦ 2

⎤ r0 v0 ⎡ Kc (e(π + ψ − φ)tan(ψ ) − 1) ⎥ ⎢ ⎦ 2 ⎣ tan(ψ ) ⎤ sin(2θ ) r0 v0 ⎡ H + ⎢ (1 − K ) ⎥ 2 ⎣ 2 1 − sin(φ) ⎦ ⎤ ⎡ e(π + ψ − φ)tan(φ) (e(π + ψ − φ)(tan(ψ )−tan(φ)) − 1) ⎥ ⎢ ⎦ ⎣ tan(ψ ) − tan(φ)

Kce 2θ tan(ψ ) dθ

⎡ sin(2θ ) ⎢ (1 − K ) ⎣ 2

+ r0 v0

He(π + ψ − φ)tan(φ) 1 − sin(φ)

∫0

φ ψ (π − + ) 2 2 2

⎤ ⎥ ⎥ ⎥ ⎥ ⎥⎦

(1 − K )

⎤ e 2θ (tan(ψ )−tan(φ)) ⎥ dθ ⎦

(68)

e(π + ψ − φ)tan(φ)

Finally, we obtain:

⎤ r v ⎡ Kc (e(π + ψ − φ)tan(ψ ) − 1) ⎥ Bf = 0 0 ⎢ ⎦ 2 ⎣ tan(ψ ) ⎤ sin(2θ ) r0 v0 ⎡ H + ⎢ (1 − K ) ⎥ 2 ⎣ 2 1 − sin(φ) ⎦ ⎤ r0 v0 ⎡ e(π + ψ − φ)tan(φ) (e(π + ψ − φ)(tan(ψ )−tan(φ)) − 1) ⎥ ⎢ ⎦ 2 ⎣ tan(ψ ) − tan(φ)

⎞ ⎛ ⎛π φ⎞ ⎟ ⎜ (tan(φ) − tan(ψ )) sin ⎜ + ⎟ cos(ψ ) ⎝4 2⎠ ⎟ ⎜ + r0 V0 (1 + sin(φ)) ⎜ + 1⎟ ⎛π φ⎞ ⎟ ⎜ sin ⎜ + ψ − ⎟ ⎝4 2⎠ ⎠ ⎝ + r0 V0 (1 + sin(φ))((tan(φ) − tan(φ))cos(ψ ) − 1)

e(π + ψ − φ)tan(ψ )

(73)

Finally, an upper estimate of the ultimate load writes:

⎛π φ⎞ sin ⎜ + ψ − ⎟ ⎝4 2⎠ Pk = H (P + P2 ) ⎞ 1 ⎛π ⎞ ⎛π cos ⎜ ⎟ sin ⎜ − ψ ⎟ ⎠ ⎝4⎠ ⎝2 (69)

where

Let us notice that the computation of the rate of dissipation in the logspiral zone may be achieved by modelling the logspiral zone as an assembly of homothetic thin and rigid layers rotating about the point O 7

(74)

International Journal of Non-Linear Mechanics 90 (2017) 1–10

M. Hamlaoui et al.

⎡ ⎛π ψ ⎞⎤ P1 = (1 − sin(ψ )) ⎢e(π + ψ − φ)tan(ψ ) tan2 ⎜ + ⎟ ⎥ ⎝4 2 ⎠⎦ ⎣ ⎡ ⎞⎤ ⎛ ⎛π φ⎞ ⎢ ⎟⎥ ⎜ sin(ψ )sin ⎜ + ⎟ ⎝4 2⎠ ⎢ 1 ⎟⎥ ⎜ + (1 − sin(ψ )) ⎢ − 1 ⎟⎥ ⎜ ⎞ ⎛π ψ 1 − sin( ) φ ⎢ ⎟⎥ ⎜ sin ⎜ + ψ − ⎟ ⎢⎣ ⎠ ⎝ 4 2 ⎠ ⎥⎦ ⎝

Table 1 Plastic limit loads for the strip footing; φ = 35 , ψ ∈ {15, 25, 30, 35}.

(75)

and

⎤ ⎡ ⎛π φ⎞ ⎥ ⎢ (tan(φ) − tan(ψ ))sin ⎜ + ⎟ cos(ψ ) ⎝4 2⎠ ⎥ ⎢ P2 = (1 + sin(φ)) ⎢ + 1⎥ ⎛π φ⎞ ⎥ ⎢ sin ⎜ + ψ − ⎟ ⎝4 2⎠ ⎦ ⎣ e(π + ψ − φ)tan(φ)

This section is devoted to the comparison between the analytic estimates Ps and Pk against the Drescher and Detournay formula [25] and the numerical results worked out by the authors in [32] and many others [3,5,37]. First, let us recall that, for finite element computations, the plastic collapse load is reached when the reaction force beneath the punch ceases to grow even if the imposed displacement continues to increase. Also, recall that the exact collapse load of a strip footing for the Mohr–Coulomb criterion with a normal flow rule is given by the wellknown Prandtl-Hill formula scaled by L and c [15]:

⎛ ⎞ ⎛π φ⎞ Pu = cot(φ) ⎜e π tan(φ) tan2 ⎜ + ⎟ − 1⎟ ⎝4 cL 2⎠ ⎝ ⎠

c* =

cos(φ) cos (ψ ) tan(φ) 1 − sin(φ)sin(ψ )

cos(φ)cos(ψ ) c 1 − sin(φ)sin(ψ )

(77)

(78)

(79)

The Drescher–Detournay model stems from the work of Davis [18] who found that, for non-associated models, velocity characteristics do not coincide with the stress characteristics (slip lines), and on velocity characteristics, the shear and normal stresses satisfy a Mohr–Coulomb like criterion with parameters c* and φ*:

τ + σn tan(φ*) = c*

Numerical result [32]

Drescher and Detournay [25]

35/35 35/30 35/25 35/20 35/15 35/10

46.1236 43.305 40.6534 38.159 35.8125 33.605

46.1236 44.9536 42.6257 39.9977 37.3771 34.8739

46.41 46.15 43.74 39.65 36.01 32.8071

46.1236 45.2817 43.0335 39.8199 36.0531 32.0707

φ /ψ (°)

Static estimate

Upper estimate

Numerical result [32]

Drescher and Detournay [25]

30/30 30/25 30/20 30/15 30/10

30.1396 28.5736 27.0845 25.6686 24.3223

30.1396 29.5307 28.2718 26.8048 25.3065

30.81 30.35 29.31 28.06 25.32

30.1396 29.7098 28.5365 26.8105 24.7204

bearing capacity  u (77) for the associated model, whereas the FEM result is greater than the exact solution. Let us recall here that numerical computations delivered in [32] have been carried out with the Durcker–Prager criterion matching the Mohr–Coulomb criterion at plane strain [15]. Moreover, for the friction angle φ = 20°, Berga's associate solution is very close to (77), whereas Bousshine's solution is strangely inaccurate for associated and non-associated cases. It should be stressed that a dispersion of results can be noted for very low values of the dilatancy angle. We remark here that inaccuracy in numerical computation occurs and is a result of the oscillations exhibited in load-settlement curves, making it difficult to define the collapse load. In [32], the limit load was defined as the average of the fluctuating values in the plastic regime. It is worth noting that, in several works in the literature [39,27,40,52] using the Mohr–Coulomb model, oscillations are observed in the load-footing settlement, even for high friction angles φ, and the increase of φ and mesh refinement leads to an increase in the intensity of fluctuations. With the high degree of non-associativity (φ = 20°, ψ = 0°), the Krabbenhoft et al. approximations are far from other estimates. From Figs. 12 and 13, it can be concluded, that for moderate and low friction (φ ≤ 20 ), the finite element solution overestimates the plastic limit load, the D–D model and the kinematic quasi-bound Pk are closed and follow the same tendency, while the static bound Ps seems to be a bad approximation for almost all values of the dilatancy angle ψ. Thus the improvement of the trial stress field is necessary.

where L is the width of the foundation and c is cohesion of the material. As mentioned before, to our knowledge no analytical solution giving the bearing capacity of a strip footing with non-associated models is available. Drescher and Detournay [25] proposed a simple method to calculate limit loads for the non-associated Mohr–Coulomb criterion. They advocated replacing the parameters c and φ by c* and φ* respectively in the Prandtl–Hill solution (77), given by:

tan(φ*) =

Upper estimate

(76)

5. Discussion

u =

Static estimate

Table 2 Plastic limit loads for the strip footing; φ = 30 , ψ ∈ {10, 15, 25, 30}.

+ (1 + sin(φ))((tan(φ) − tan(φ))cos(ψ ) − 1)

e(π + ψ − φ)tan(ψ )

φ /ψ (°)

6. Conclusion In this paper the bipotential theory and the slip-line method have been used to derive estimates of the ultimate loads of rigid strip resting on a non-standard Mohr–Coulomb half-plane. The punch is subjected to a normal pressure and the weight of substrate is neglected. The major contributions of this paper are as follows:

(80)

The Drescher–Detournay solution will be denoted as the ‘D–D’ solution in the following text. The comparison between the values of (Pk,Ps) with the D–D solution and the FEM results for the non-associated strip footing are reported in Tables 1–4. For more convenience, and in order to facilitate interpretation of results, the curves of Pk, Ps, D–D and numerical ultimate loads normalized by the Prandtl–Hill solution (77) are plotted in Figs. 10–13. For the friction angle φ = 20° and different dilatancy angles ψ, finite element results for smooth and rigid footing derived by Berga [3] and Bousshine [5] using the finite element method within the bipotential framework are compared to our analytical estimates. First, notice that (Pk,Ps) and the D–D formula lead to the exact

• • • 8

The bipotential, representing the plastic dissipation of the nonassociated Mohr–Coulomb criterion is delivered, and its expression along a velocity discontinuity surface is provided. Analytical static and kinematic quasi-bounds are given based upon a statically and plastically trial stress field and a Prandtl-like collapse mechanism. It should be pointed out that these analytic estimates can be considered as the first analytic ones available in the literature. The obtained results are compared to existing numerical computations [32] and the formula proposed by Drescher and Detournay

International Journal of Non-Linear Mechanics 90 (2017) 1–10

M. Hamlaoui et al.

Table 3 Plastic limit loads for a rough footing; φ = 20 , ψ ∈ {0, 5, 10, 15, 20}.

φ /ψ (°)

Static estimate

Upper estimate

Numerical result [32]

Drescher and Detournay [25]

Berga [3]

Bousshine [5]

Krabbenhoft et al. [37]

20/20 20/15 20/10 20/5 20/0

14.8347 14.285 13.7525 13.2367 12.737

14.8347 14.6404 14.2146 13.6931 13.1393

15.146 15.049 14.794 14.29 13.68

14.8347 14.7011 14.3232 13.7409 12.9977

14.938 – – – 13.5886

14.5528 – 13.6628 – 13.2029

15.2 – – 14.1 –

Table 4 Plastic limit loads for a rough footing; φ = 15 , ψ ∈ {0, 5, 10, 15}.

φ /ψ (°)

Static estimate

Upper estimate

Numerical result [32]

Drescher and Detournay [25]

15/15 15/10 15/5 15/0

10.9765 10.6366 10.3045 9.98008

10.9765 10.86 10.599 10.273

11.20 11.154 11.011 10.6654

10.9765 10.8971 10.6689 10.3104

Fig. 12. Dimensionless collapse loads; φ = 20 , ψ ∈ {0, 5, 10, 15, 20}.

Fig. 10. Dimensionless collapse loads; φ = 35 , ψ ∈ {10, 15, 25, 30, 35}.

Fig. 13. Dimensionless collapse loads; φ = 15 , ψ ∈ {0, 5, 10, 15}.

the collapse load. However, the D–D formula cannot be improved, whereas the bipotential approach proposed in this paper can be improved by considering more relevant mechanisms and stress field. This work constitutes a first promising attempt to apply bipotential theory to the study of the plastic limit load of the rigid punch indenting a non-associated half-plane, and opens the way to new analytical investigations. Further researches should be carried out to improve the results presented in this paper. There are several possibilities: (i) considering unknown bases angles α and β of wedges in the Prandtl-like mechanism and then performing optimization of the rate of dissipation and (ii) finding more refined trial (S.A.) stress fields in order to improve both the static and the kinematic quasi-bound. These are topics for future work.

Fig. 11. Dimensionless collapse loads; φ = 30 , ψ ∈ {10, 15, 25, 30}.

[25]. In most cases, a good agreement between the kinematic quasibound and D–D formula is observed, whereas the static bound seems to underestimate the collapse load. Moreover, for a high degree of non-associativity, numerical computations must be considered with caution. The main conclusion to be drawn here is that the D–D formula and the kinematic quasi-bound deliver reliable estimates of 9

International Journal of Non-Linear Mechanics 90 (2017) 1–10

M. Hamlaoui et al.

implicit standard material approach, Mech. Struct. Mach. 19 (1991) 301–325. [25] A. Drescher, E. Detournay, Limit load in translational mechanisms for associative and nonassociative materials, Géotechnique 43 (4) (1993) 443–456. [27] H.L. Erickson, A. Drescher, Bearing capacity of circular footings, J. Geotech. Geoenviron. Eng. 128 (2002) 38–43. [28] W. Fenchel, On conjugate convex functions, Can. J. Math. 1 (1949) 73–77. [29] S. Frydman, H.J. Burd, Numerical studies of bearing-capacity factor Nc, J. Geotech. Geoenviron. Eng. 123 (1997) 20–29. [30] D.V. Griffiths, Computation of bearing capacity factors using finite elements, Géotechnique 32 (1982) 195–202. [32] M. Hamlaoui, A. Oueslati, B. Lamri, G. de Saxcé, Finite element analysis of the plastic limit load of strip foundations with non-associated Drucker–Prager model, Eur. J. Environ. Civ. Eng. 19 (2015) 1179–1201. [33] M. Hjiaj, J. Fortin, G. Saxcé, A complete stress update algorithm for the nonassociated Drucker–Prager model including treatment of the apex, Int. J. Eng. Sci. 41 (2003) 1109–1143. [34] M. Hjiaj, A.V. Lyamin, S.W. Sloan, Numerical limit analysis solutions for the bearing capacity factor Ny, Int. J. Solids Struct. 42 (2005) 1681–1704. [35] R. Hill, Mathematical Theory of Plasticity, Oxford University Press, London, 1950. [36] K. Krabbenhft, A.V. Lyamin, S.W. Sloan, Formulation and solution of some plasticity problems as conic programs, Int. J. Solids Struct. 44 (2007) 1533–1549. [37] K. Krabbenhoft, M.R. Karim, A.V. Lyamin, S.W. Sloan, Associated computational plasticity schemes for nonassociated frictional materials, Int. J. Numer. Methods Eng. 90 (2012) 1089–1117. [38] J. Lemaitre, Formulation and identification of damage kinetic constitutive equation, in: D. Krajcinovic and j. Lemaitre (Eds.), Continuum Damage Mechanics, Theory and Application, CISM Courses and Lectures, vol. 295, Springer-Verlag, NY, 1987. [39] D. Loukidis, R. Salgado, Bearing capacity of strip and circular footings in sand using finite elements, Comput. Geotech. 36 (5) (2009) 871–879. [40] N. Manoharan, S.P. Dasgupta, Bearing capacity of surface footings by finite elements, Comput. Struct. 54 (1995) 563–586. [41] A.A. Markov, On variational principles on theory of plasticity, Prok. Mat. Mekh. 11 (1947) 339–350. [42] J.J. Moreau, Fonctionnelles convexes, Istituto Poligrafico e Zecca dello Stato, Rome, 2003. [43] J.J. Moreau, La notion de sur-potentiel et les liaisons unilatérales en élastostatique, C.R. Acad. Sci. Sér. A 267 (1968) 954–957. [44] D.H. Nguyen, principes variationnels à deux champs pour le matériau de SaintVenant–Von Mises, Bull. Acad. Sci. Pol. Sér. Sci. Tech. 14 (7/8) (1976) 325–331. [45] P.D. Panagiatopoulos, Inequality Problems in Mechanics and Applications, Birkhauser, Boston, 1975. [46] L. Prandtl, Über die Häete plastischer Körper, Nachr. Ges. Wiss. Gött., Math.-Phys. Kl. (1920) 74–85. [47] D.M. Potts, L. Zdravkovic, Finite Element Analysis in Geotechnical Engineering: Applications, Thomas Telford Ltd., London, 2001. [48] R.T. Rockafellar, Convex Analysis, Princeton University Press, Princeton, 1970. [49] S.W. Sloan, P.W. Kleeman, Upper bound limit analysis using discontinuous velocity fields, Comput. Methods Appl. Mech. Eng. 127 (1995) 293–314. [50] M. Save, C. Massonnet, G. de Saxcé, Plastic Limit Analysis of Plates, Shells and Disks, Elsevier, Amsterdam, 1997. [52] J. Yin, Y. Wang, A. Selvadurai, Influence of nonassociativity on the bearing capacity of a strip footing, J. Geotech. Geoenviron. Eng. 127 (2001) 985–989. [53] N. Zouian, I. Pontes Filho, L. Borges, L. Moutada Costa, Plastic collapse in nonassociated hardening materials with application to Cam-clay, Int. J. Solids Struct. 44 (2007) 4382–4398.

References [1] A. Berga, Mathematical and numerical modeling of the non-associated plasticity of soils – part 1: the boundary value problem, Int. J. Non-Linear Mech. 47 (2012) 26–35. [2] A. Berga, Mathematical and numerical modeling of the non-associated plasticity of soils – part 2: finite element analysis, Int. J. Non-Linear Mech. 47 (2012) 36–45. [3] A. Berga, calcul élastoplastique des sols à lois non-associées par éléments finis basé sur l′approche des matériaux standards implicites (Ph.D. thesis), Université de Technologie de Compiègne, 1993. [4] C. Bouby, G. de Saxcé, J.-B. Tritsch, A comparison between analytical calculations of the shakedown load by the bipotential approach and step-by-step computations for elastoplastic materials with nonlinear hardening, Int. J. Solids Struct. 43 (2006) 2670–2692. [5] L. Bousshine, Simulation numérique des processus de poinçonnement, de coupe et de formage des métaux et des matériaux granulaires (Ph.D. thesis), Faculté Polytechnique de Mons, 1994. [6] L. Bousshine, A. Chaaba, G. De Saxcé, Plastic limit load of plane frames with frictional contact supports, Int. J. Mech. Sci. 44 (2002) 2189–2216. [7] L. Bousshine, A. Chaaba, G. de Saxcé, Softening in stress–strain curve for Drucker– Prager non-associated plasticity, Int. J. Plast. 17 (1) (2001) 21–46. [8] M. Buliga, G. de Saxcé, C. Vallée, A variational formulation for constitutive laws described by bipotentials, Math. Mech. Solids 18 (2013) 78–90. [11] R. Casciaro, A. Di Carlo, G. Valente, Plane stress limit analysis by finite elements, in: Proceedings of the International Symposium on Discrete Methods in Engineering, 1974, pp. 171–181. [12] R. Casciaro, A. Di Carlo, Mixed finite element models in limit analysis, in: J.T. Oden (Ed.), Computational Methods in Nonlinear Mechanics, Austin, Texas, 1974, pp. 171–181. [13] A. Chaaba, L. Bousshine, G. de Saxcé, Kinematic limit analysis of non associated perfectly plastic material by the bipotential approach and finite element method, J. Appl. Mech. – Trans. ASME 77 (2010) 031016. [14] E. Christiansen, Computation of limit loads, Int. J. Numer. Methods Eng. 17 (1981) 1547–1570. [15] W.F. Chen, Limit Analysis and Soil Plasticity, Elsevier, New York, 1975. [16] W.F. Chen, Limit Analysis in Soil Plasticity, Elsevier, Amsterdam, 1990. [17] L. Cheng, J. Yun, A. Oueslati, G. de Saxcé, D. Kondo, Plastic limit state of the hollow sphere model with non-associated Drucker–Prager material under isotropic loading, Comput. Mater. Sci. 62 (2012) 210–215. [18] E.H. Davis, in: I.K. Lee, (Ed.), Theories of Plasticity and Failure of Soil Masses. Soil Mechanics: Selected Topics, Butterworth's, London, 1968, pp. 341–380. [19] E.H. Davis, J.R. Booker, Some adaptations of classical plasticity theory for soil stability problems, in: A.C. Palmer (Ed.), Proceedings of the Symposium on the Role of Plasticity in Soil Mechanics, Cambridge University Press, Cambridge, UK, 1973, pp. 24–41. [20] P. de Buhan, Plasticité et calcul à la rupture. Presses de l′Ecole Nationale des Ponts et Chaussées, 2007. [21] G. de Saxcé, Une généralisation de l′inégalité de Frenchel et ses applications aux lois constitutives. C.R. Acad. Sci. t.314 (série 2) (1992) 125–129. [22] G. de Saxcé, L. Bousshine, The limit analysis theorems for the implicit standard materials: application to the unilateral contact with dry friction and the non associated flow rules in soils and rocks, Int. J. Mech. Sci. 40 (4) (1998) 387–398. [23] G. de Saxcé, L. Bousshine, Implicit Standard Materials, CISM Courses and Lectures, vol. 295, Springer-Verlag, NY; 2002. [24] G. de Saxcé, Z.Q. Feng, New inequality and functional for contact with friction: the

10