## Abstract

The non-coaxiality of the directions of principal stresses and principal plastic strain rates in granular soils under stress rotations has long been observed and recognized in soil tests using both simple shear and hollow cylinder apparatuses. A few constitutive theories have also been proposed in the literature to account for the effect of stress rotations and the subsequent non-coaxial soil behaviour, particularly in the context of shear band analysis. However, the lack of corresponding general numerical methods makes it difficult to investigate the influence of non-coaxial stress–strain behaviour on the results of geotechnical boundary value problems. This paper presents a numerical evaluation of a class of non-coaxial, elastic–plastic models that are developed by combining the conventional plastic potential theory and the double shearing theory. The general non-coaxial constitutive theories are first formulated and then a finite element implementation of the theories is carried out. To evaluate the non-coaxial theories, the problem of simple shear of soils is chosen to investigate the predicted behaviour of soils under simple shear loading conditions where the axes of principal stresses rotate. In particular, the influence of initial stress states and the degree of non-coaxiality are examined. It is found that the numerical results predicted using the non-coaxial model are in general agreement with the experimental observations reported in the literature.

## 1. Introduction

Laboratory element testing is essential for developing useful soil models and measuring soil parameters for geotechnical analysis and design. Among many existing testing devices, two have been developed to impose a simple shear condition on soil specimens. The first simple shear device was initially developed by Kjellman (1951), and subsequently modified by Bjerrum & Landva (1966) at the Norwegian Geotechnical Institute. It uses a cylindrical specimen enclosed in a wire-reinforced membrane. The second simple shear device, described by Roscoe (1953) at Cambridge, uses a parallelepiped, square specimen enclosed by rigid platens. Since then other simple shear devices have been developed and used in geotechnical research and practice (e.g. Matsuoka 1974; Oda & Konishi 1974; Allersma 1982; Ochiai *et al*. 1983; Airey *et al*. 1985; Joer *et al*. 1998).

An important feature of simple shear testing is that the direction of the principal stress rotates as the shear stress is applied on the specimen. Roscoe *et al*. (1967) showed that the principal axes of strain rate and of stress are generally not coincident during simple shear tests of sand. Figure 1 shows the experimental results reported by Roscoe (1970). The rotations of the principal stress and the principal plastic strain rate are non-coaxial particularly at the early stage of loading. However, the axes tend to become coincident at large shear strains. Oda & Konishi (1974) carried out two-dimensional simple shear tests on an assembly of a cylinder (made of photoelastic material) packed at random to simulate the shear deformation of sand, and derived a similar conclusion to that of Roscoe *et al*. (1967). Drescher & de Josselin de Jong (1972) carried out an experimental micro-mechanical study using a photoelastic disc assembly as a two-dimensional analogue of granular media and presented further experimental evidence for non-coaxiality. In addition, other researchers have also observed that the direction of the principal stress deviates from that of the principal strain increment under rotation of the principal stress axes (Arthur *et al*. 1980; Gutierrez *et al*. 1991; Joer *et al*. 1998). These all point to the fact that the plastic behaviour of granular material is generally non-coaxial when subjected to stress rotations. Apart from the experimental research, non-coaxial plastic flow rules have also been derived theoretically from a micro-mechanical study of granular materials (Christoffersen *et al*. 1981).

Non-coaxiality has been studied theoretically within the context of fully developed flows (i.e. steady shear flow) in perfectly plastic granular materials. An early kinematic model for granular flow was developed by de Josselin de Jong (1958, 1971). His so-called ‘double sliding, free rotating model’ for planar flow was based on the assumption of shear flow occurring along two surfaces where the available shear resistance has been exhausted. Using the same concept of double shearing, Spencer (1964) established a set of kinematic equations termed as the ‘double shearing model’ but with a different rotation term from that of de Josselin de Jong (1971). Although these two models have a common physical basis, they are slightly different in nature (Mandl & Luque 1970; Spencer 1982; Harris 1995; Savage & Lockner 1997). In particular, the ‘double sliding, free rotating model’ of de Josselin de Jong is indeterminate in nature and additional assumptions are required before it can be used to solve boundary value problems. For example Teunissen & Vermeer (1988) proposed a version of de Josselin de Jong's model by assuming that a minimum amount of work is stored as elastic energy. In this way, the model becomes determinate, but it has not yet been shown to be able to reproduce the non-coaxial behaviour as observed in simple shear tests (Roscoe 1970).

The classical potential theory may also predict a non-coaxial behaviour if the plastic potential is assumed to be an anisotropic function of the stress tensor. However, this approach is unable to model the dependence of plastic strain rates on the stress rate tangential to the yield surface. This important feature has been observed in many soil experiments for many years (e.g. Tsutsumi & Hashiguchi 2005).

The validity of any proposed non-coaxial, elastic–plastic models has to be assessed by comparing their predicted behaviour with experimental data (such as that obtained from simple shear tests). Existing numerical analyses of the simple shear test are mainly based on the coaxial theories, in which the principal axes of stress and plastic strain rate are assumed to coincide (Potts *et al*. 1987; Budhu & Britto 1987; Vermeer 1990). This assumption clearly contradicts the experimental data shown in figure 1. In view of the above discussions, the objectives of the present paper are as follows:

To present the theory of a class of non-coaxial, elastic–plastic soil models with its links established with both the plastic potential theory and the double shearing theory;

To carry out a finite element implementation of the proposed non-coaxial models;

To use a simple shear problem as an example to assess numerically the predictive capability of the proposed soil models in light of experimental observations.

In contrast to the approach used in the standard double shearing model, which is based on kinematic assumptions, this paper sets out to establish a general non-coaxial flow rule by extending the conventional plastic potential theory (Hill 1950) to account for non-coaxiality. Moreover, the scope of the modelling in the paper includes the range of pre-failure deformation and, therefore, represents an extension to the original double shearing models that were proposed for describing fully developed plastic flow only (Spencer 1964, 1982; Mandl & Luque 1970; Mehrabadi & Cowin 1978). In agreement with Anand (1983) and Savage & Lockner (1997), the results of this present study suggest that it is necessary to relax the original kinematic hypothesis of the coincidence of stress and velocity characteristics in order to allow the double shearing concept to be used more successfully in the range of pre-failure deformation.

## 2. A class of non-coaxial plasticity theories

For simplicity, granular material is often modelled as an elastic perfectly plastic material. Given that the double shearing theory of Spencer (1964, 1982) has been developed and discussed largely in the context of perfect plasticity and plane strain conditions, this paper will also focus on a perfectly plastic formulation under plane strain conditions. However, it should be pointed out that extension to isotropic or kinematic hardening plasticity and non-plane strain cases is possible, which will be addressed elsewhere. All stresses are assumed to be effective stresses and tensile stresses are taken as positive.

An incremental stress–strain relationship for an elastic–plastic material in plane strain can be written as follows:(2.1)where(2.2)Under plane strain conditions we have *ϵ*_{xx}=0. The out-of-plane stress *σ*_{xx} is assumed to be the intermediate principal stress. is the plastic strain rate. Although the stress rate in equation (2.1) is strictly an objective stress rate, the rotation term may be neglected for small strain problems (as will be considered in this paper). This has been the standard assumption used in almost all the major existing geotechnical analyses.

The elasticity matrix can be shown to be as follows:(2.3)where *E* is Young modulus and *ν* is Poisson's ratio.

To synthesize the double shearing theory and the plastic potential models, we closely follow Spencer (1982) and Harris (1993, 1995) in proposing the formulation of a class of non-coaxial models that take the following general form in terms of plastic strain rate :(2.4)where *f* and *g* are the yield function and plastic potential, respectively, and *α* are scalar functions and *α* is dimensionless. The advantage of using the flow rule equation (2.4) is that it reduces to the conventional plastic potential theory when *α* is set to zero. As will be explained later, the flow rule equation (2.4) can also be made as identical to the double shearing theory if *α* is given by equation (2.7).

The vector * t* is defined as(2.5)where

*θ*

_{σ}is the angle between the minor principal compressive stress direction and the

*x*-axis, defined by(2.6)Equation (2.4) is a generalization of the double-shearing model and classical plastic potential theory, to which it reduces when

*α*≡0. In the case when

*α*is not equal to zero, the strain rate and stress

*σ*are not coaxial. The plastic strain rate given by the first term of the equation (2.4) is the same as that from the conventional coaxial plastic theory, which is normal to a plastic potential surface. The second term is the non-coaxial, stress-rate dependent, plastic strain component that is tangential to the yield locus in the deviatoric plane.

In the original double shearing theory of Spencer (1964, 1982) and Mehrabadi & Cowin (1978), the stress and velocity characteristic directions are known to coincide. It can be shown (Harris 1993) that this property of the original double-shearing models leads to the following expression for *α*:(2.7)where *ϕ* is the friction angle and *ψ* is the dilation angle of the material. In other words, the models of Spencer (1964) and Mehrabad & Cowin (1978) may be recovered from equation (2.4) by using equation (2.7) to determine *α*. In addition, it can also be shown that for plane strain conditions the present flow rule equation (2.4) becomes identical in form to the yield-vertex flow rule of Rudnicki & Rice (1975) if we assume a pressure-dependent expression for *α* such that *α*=*R/h*_{1} where and *h*_{1} is a plastic modulus used by Rudnicki & Rice (1975) to account for non-coaxial plastic strain rates.

For a non-coaxial behaviour, the directions of principal stress and principal plastic strain rate are different. The definition of , which is the angle between the minor principal compressive plastic strain rate and the *x*-axis, is expressed by equation (2.8) and shown in figure 2.(2.8)

## 3. Incremental non-coaxial stress–strain relationship

Equation (2.4) represents a general plane strain, non-coaxial theory that combines the double shearing theory with the conventional plastic potential theory. It is obvious from equation (2.4) that the non-coaxial plastic strain components are controlled by the value of *α* and stress rate (as can be expressed in terms of stress rate ). In order to implement the non-coaxial theory equation (2.4) into a finite element program for use to solve general boundary value problems in geotechnical engineering, it is necessary to derive a general incremental stress–strain relationship.

By combining equations (2.5) and (2.6), we obtain(3.1)in which the matrix **H**_{s} is defined by(3.2)where(3.3)(3.4)(3.5)By using the consistency condition, the constitutive equation of a non-coaxial model defined by equation (2.4) for an elastic perfectly plastic material can be shown to be(3.6)where(3.7)and the modified elastic stiffness matrix is given by(3.8)(3.9)in which [* I*] is the identity matrix.

The Mohr–Coulomb yield criterion has been widely used to model granular soils. In plane strain, it takes the following form(3.10)where *c* is known as cohesion. As shown by Spencer (1982), the Mohr–Coulomb yield condition (3.10) may be graphically represented by a cone in the space, as illustrated in figure 3*a*. The section defined by (*σ*_{xx}+*σ*_{yy})/2=constant in this space is shown in figure 3*b*. The principal stress angle *θ*_{σ} and the plastic deviatoric strain components defined in equation (2.4) are also plotted in figure 3*b*.

It is noted that the general non-coaxial plastic flow rule, defined by equation (2.4) and shown graphically in figure 3, is capable of reproducing much of the experimental findings on the non-coaxiality of the directions of principal stresses and plastic strain rates. Figure 4 presents the experimental results obtained by Gutierrez & Ishihara (2000) using hollow cylinder tests by varying the directions of principal stresses. The figure shows the plastic strain increment directions obtained during tests involving pure principal stress rotation at different constant level of mobilized friction angles. As clearly shown in figure 4, the unit plastic strain increments (or rates) deviate substantially from the radial directions which correspond to the orientations of the principal stresses. For example, for the test at a constant shear stress corresponding to a mobilized friction angle of 20°, the average angle of non-coaxiality is about 30°. It is also evident that the degree of non-coaxiality decreases with increasing shear stress level and this is consistent with the simple shear data reported by Roscoe (1970).

## 4. Finite element implementation

To evaluate the performance of the proposed non-coaxial stress–strain theories, it is essential to implement them into a finite element program to allow for their detailed evaluation and wide applications. In this study, the finite element program, known as ABAQUS (2001), has been chosen. The proposed stress–strain law was incorporated using its user subroutine UMAT.

To overcome the sharp vertices at the lode angle *θ*=±30° in the Mohr–Coulomb surface, the method of Owen & Hinton (1980) was adopted. Zienkiewicz & Pande (1977) proposed a hyperbolic approximation to remove the singularity at the tip of the surface. Following Zienkiewicz & Pande (1977), Abbo (1997) obtained a close approximation to the straight line which defines the Mohr–Coulomb yield surface. The rounded hyperbolic Mohr–Coulomb yield function (Abbo 1997) is adopted in our own subroutine for the analysis. This yield function is defined as(4.1)where(4.2)(4.3)(4.4)where the lode angle *θ* is defined by(4.5)with(4.6)The parameter *η* in equation (4.1) is used to adjust the closeness between the hyperbolic surface and the original Mohr–Coulomb surface. The original Mohr–Coulomb yield function can be recovered by setting *η* to zero. In practice, however, setting *η*=0.05*c* cot*ϕ* is found to give results which are almost identical to those from the original Mohr–Coulomb model (Abbo 1997).

The corresponding plastic potential, which is similar to the yield function and takes into account of the effect of dilation angle *ψ*, is given below:(4.7)where(4.8)(4.9)

## 5. Numerical analysis of simple shear problems

One of the purposes of this study is to check the validity of a class of new non-coaxial models using simple shear problems. The ideal simple shear problem has been solved by Hansen (1961) analytically with coaxial plasticity for zero dilation angle. In this paper, the shear stress–strain curves were obtained numerically by imposing the simple shear boundary conditions to an isoparametric, eight-noded, plane strain finite element. This assumption has also been used in the finite element studies of Potts *et al*. (1987). A reduced integration rule is used because it is known to perform better for problems involving incompressible behaviour (Yu & Netherton 2000). Newton's method was used for solving the nonlinear finite element equations. The Riks method was adopted in the case of unstable problems that exhibit softening. Modified Euler scheme with substepping (Sloan 1987; Abbo 1997) was implemented in the subroutine to integrate the elastic–plastic stress–strain equations to update the stresses at each integration point.

For clarity of presentation, a purely frictional material will be considered in the numerical analyses. A constant uniform stress *σ*_{yy}=200 kPa was used in all analyses, and the initial horizontal stress given by *σ*_{xx}=*K*_{0}*σ*_{yy} is varied.

## 6. Results and discussion

### (a) Coaxial plasticity (*α*=0)

The elastic–plastic model defined by equation (2.4) becomes coaxial when *α* is equal to zero. In order to validate the finite element implementation of elastic–plastic models, the finite element calculations for the coaxial plasticity are discussed first.

To check the accuracy of the finite element formulation and solution procedures, numerical analyses of simple shear problems are carried out with the same set of soil parameters as those used by Hansen (1961) in deriving his analytical solutions. Figure 5 shows a direct comparison between the numerical and analytical solutions for the stress ratio—shear strain relationship. *K* in the horizontal axis in figure 5 is defined by Hansen (1961) as *K*=*E*/(1−*ν*^{2}). It is evident that the numerical results are practically identical with the analytical results.

For a simple shear problem, the horizontal plane may be regarded as a slip line (i.e. a velocity characteristic or velocity discontinuity) at large strains (e.g. ultimate failure states) when deformation is plastic. For a purely frictional soil, Davis (1968) first showed that on a slip line the Mohr–Coulomb failure criterion implies the following ratio of shear stress to normal stress:(6.1)For the special case of a non-dilatant soil (*ψ*=0), the above ultimate stress ratio simply reduces to sin *ϕ*.

Figure 6 shows the shear stress–strain curves of several cases during simple shear. The friction angle for all cases is 35°, and elastic modulus is 26 000 kPa, Poisson's ratio is 0.3. When initial stress ratio *K*_{0} is equal to 0.43, the ultimate stress ratio is 0.579 for the dilation angle of 0°, 0.608 for the dilation angle of 5° and 0.71 for the dilation angle of 35°. These numerical values agree well with those from Davis's analytical solution (6.1).

The initial *K*_{0} value and the angle of dilation affect the location of the peak value and shape of the stress–strain curve between peak condition and ultimate conditions. For non-associated flows when *ψ*<*ϕ*, it is possible that the coaxial, perfectly plastic material model exhibits strain softening in simple shear, purely as a consequence of a high *K*_{0} and the applied boundary conditions. In other words, the combination of a high *K*_{0} and a non-associated flow rule can produce an apparent strain-softening response, as also noted by Vermeer (1990). Figure 6 presents this strain-softening phenomenon. Hansen (1961) note that the value of the stress ratio at the peak point is(6.2)The peak values for the dilation angle of 0° and 5° are the same and are equal to 0.7 as predicted from equation (6.2).

On the other hand, for an associated flow rule with *ϕ*=*ψ*, slip lines coincide with stress characteristics (when the Mohr–Coulomb condition is satisfied). Therefore, on the horizontal plane, the ultimate ratio of the shear stress to normal stress would also be tan *ϕ* as for the peak ratio. The values *σ*_{xy}/σ_{yy} for different initial *K*_{0} values are 0.7.

In summary, for coaxial plasticity when *α* is equal to zero, the numerical results agree very well with the analytical solutions. This testifies the correctness of the finite element implementation of the constitutive models.

### (b) Non-coaxial plasticity (*α*≠0)

#### (i) Cases with *α*<0—the original double shearing theory

To obtain the original double shearing theory (Spencer 1964, 1982; Mehrabadi & Cowin 1978) from the general non-coaxial theory equation (2.4), the non-coaxial coefficient *α* may be calculated using equation (2.7). In reality, the dilation angle *ψ* for a soil is smaller than the angle of friction *ϕ*. Thus equation (2.7) suggests that *α* will be a negative value.

Presented in figure 7 are the numerical results for the ratio of shear stress to normal stress against shear strain. For a non-associated flow rule with *ψ*<ϕ, the numerical results suggest that an increasing plastic shear strain is accompanied by a decreasing shear stress. It is noted that the analytical study of a simple shear problem reported by Spencer (1986, 2003) using his incompressible double shearing theory also indicates that shear flow takes place with a decreasing shear stress.

#### (ii) Cases with *α*>0—the generalized double shearing theory

As previously shown by Anand (1983), it is necessary to abandon the original kinematic hypothesis that slip lines coincide with stress characteristics to allow the double shearing concept to be used more successfully in the range of pre-failure deformation. Such a generalization has been made by Anand (1983) to apply the double shearing concept to shear band analysis. Further experimental evidence in supporting this modification has been provided recently by Savage & Lockner (1997) whose experiments on simple shear deformation suggest that the slip lines in the velocity field do not generally coincide with the Coulomb shears in the stress field.

Using the analysis of Harris (1993), it then follows that the non-coaxial coefficient *α* can take a positive value if the stress and velocity characteristic directions are different. We will refer to this resulting model as the generalized double shearing theory. To check the influence of the value of *α* on the predicted behaviour, the values of 0.02 and 0.05 have been used in numerical calculations. The material properties used for seven different cases are detailed in table 1.

Presented in figures 8–10 are the stress–strain curves for the case with *K*_{0}=0.43 (i.e. a normally consolidated sample). It is clear from the figures that an increasing value of the non-coaxial coefficient leads to a softer stress–strain response and a smaller dilation. However, the ultimate stress ratios for different *α* values tend to approach the value of the coaxial plasticity when *α*=0 at large shear strains.

In the simple shear test, the principal stress axes rotate during the increase of shear stress and the stress ratio *σ*_{xy}/*σ*_{yy} on the horizontal plane is related to the inclination angle *θ*_{σ} of the principal stress axis to the horizontal direction (Oda & Konishi 1974). At the ultimate failure when the deformation is plastic, however, any horizontal plane is a velocity characteristic and is therefore inclined at (45°+*ψ*/2) to the direction of the principal stress axis (Davis 1968).

With coaxial plasticity, the directions of principal stress and principal plastic strain are coincident. The changes of directions for principal stresses and principal plastic strain increments for Cases 1, 2 and 3 are shown in figures 11–13. When *α* is equal to zero (coaxial plasticity), the ultimate value is 45° for Case 1, 47.5° for Case 2 and 62.4° for Case 3, which are consistent with the theoretical predictions. For non-coaxial plasticity, however, figures 11–13 show clearly that the difference between orientations of principal stresses and principal plastic strain increment enlarges as the non-coaxial coefficients increase, and decreases with increasing shear strains. In general, the pattern of this predicted non-coaxial behaviour is consistent with simple shear experimental data, shown in figure 1, as reported by Roscoe (1970).

Hansen (1961) and Potts *et al*. (1987) noted that for coaxial plasticity the horizontal stress *σ*_{xy}=*K*_{r}*σ*_{yy} at the ultimate stress state is defined by(6.3)For the special case *K*_{0}=*K*_{r}, the minor principal compressive stress is already orientated at (45°+*ψ*/2) to the horizontal direction at the peak point, and the stress–strain relationship for the element becomes bilinear. This special situation corresponds to Case 4 where *K*_{r} is equal to 1.98. This case has been studied by Potts *et al*. (1987) using a coaxial plasticity model. Here further calculations are carried out for non-coaxial plasticity with non-zero *α* values. Figure 14 presents the stress ratio–shear strain and volumetric shear strain relations for this case. The numerical results suggest that the directions of principal stress and principal plastic strain rate coincide once first yield occurs.

For Cases 5 and 6 with a high *K*_{0} value and a non-associated flow rule, a strain softening response is observed (figures 15 and 16). As discussed by Vermeer (1990) for a coaxial model, this softening is not due to the material strength reduction but comes entirely from the decrease of the horizontal compressive stress in the soil when a non-associated flow rule is used. This strain softening phenomenon is absent as long as an associated flow rule is employed, as for Case 7 and shown in figure 17. However, in all the cases the ultimate stress ratio values appear to approach the coaxial value at large shear strains. In addition, the inclusion of non-coaxial behaviour (i.e. positive *α* value) tends to delay the onset of this non-associated strain softening and make the rate of softening slower. From these results, it can be concluded that the non-associated strain softening effect depends on initial stress states and the flow rule as well as the degree of non-coaxiality of the soil.

The changes of principal stress and plastic strain rate directions are shown in figures 18–20 for Cases 5, 6 and 7 respectively. For coaxial plasticity when *α* is zero, it is noted that the ultimate values at large shear strains are identical to the analytical solutions, namely 45° for Case 5, 47.5° for Case 6 and 62.5° for Case 7, respectively.

The numerical results for the principal stress and strain rate orientations suggest that the difference between these two axes is the largest at small shear strains and tends to decrease with increasing shear strain. This difference also increases with the increasing value of the non-coaxial coefficient *α*. Another important numerical result from figures 11–13, 18–20 is that the direction of principal plastic strain rate is always ‘in advance’ to that of principal stress. This phenomenon has been observed in experiments and reported by many researchers (e.g. Gutierrez *et al*. 1991; Joer *et al*. 1998).

In order to apply the theory of non-coaxial plasticity in practical geotechnical problems, *α* value needs to be determined for a given material. The current research suggests that the simple shear test could be used for this purpose by comparing numerical simulations with experimental results. Of course, other stress rotation soil tests, such as the hollow cylinder test, may also be used to determine the non-coaxial coefficient *α*.

## 7. Concluding remarks

Non-coaxiality refers to the deviation of the directions of principal stresses and principal plastic strain rates, and generally occurs during loading involving rotation of the principal stress directions. The simple shear test allows this non-coaxial behaviour to be observed and studied as the principal stresses rotate during shear loading. This paper presents and evaluates a general non-coaxial elastic-plastic model that incorporates both the conventional coaxial plastic potential theory and the double shearing theory as special cases. The theory has, for the first time, been numerically implemented into a finite element program. The results from numerical analyses of the simple shear test are found to be in general agreement with the experimental data reported in the literature.

For simple shear testing of a normally consolidated sample (i.e. *K*_{0}<1), the non-coaxial theory predicts lower shear to normal stress ratios than a coaxial theory, and the difference decreases with an increasing shear strain. However, for heavily overconsolidated soils with a high *K*_{0} value and a non-associated flow rule, strain-softening behaviour is observed. In these cases, the non-coaxial theory tends to delay the onset of the softening response. In agreement with experiments, the numerical analyses reported in this paper suggest that the difference between the directions of principal stresses and plastic strain rates is largest at the initial stage of loading and becomes smaller when shear strains increase. In addition, the numerical results are also consistent with the experimental observation that the direction of the principal plastic strain rate is always ‘in advance’ to that of the principal stress.

Finally, it is noted that the non-coaxial theory presented in the paper has been formulated in an isotropic framework. This is justified as many granular materials are isotropic on the macroscopic scale. As pointed out by Spencer (1982), any theory for anisotropic materials will necessarily be more complicated than the corresponding theory for isotropic materials, and so it is natural to focus on the isotropic case in the first instance. The non-coaxial theory described in the paper can, however, be readily extended to include the effect of both inherent and strain-induced anisotropy by adopting an anisotropic yield function and a kinematic hardening law, as has been carried out by Papamichos & Vardoulakis (1995) for the yield vertex theory.

## Acknowledgments

The work reported in this paper was supported by an EPSRC grant and the authors are grateful for this support. The authors are also very grateful to Professor A. J. M. Spencer FRS for his encouragement and stimulating discussions on the double shearing theory. Discussions with Dr D. Harris and Professor I. F. Collins at the initial stage of the project are also much appreciated. Finally the authors wish to thank three referees for their constructive and helpful comments.

## Footnotes

- Received December 17, 2004.
- Accepted October 11, 2005.

- © 2005 The Royal Society