An intergranular strain concept for material models formulated as rate equations.

Summary The intergranular strain concept was originally developed to capture the small‐strain behaviour of the soil with hypoplastic models. A change of the deformation direction leads to an increase of the material stiffness. To obtain elastic behaviour for smallstrains, only the elastic part of the material stiffness matrix is used. Two different approaches for an application of this concept to nonhypoplastic models are presented in this article. These approaches differ in the determination of the elastic stress response, which is used for reversible deformations. The first approach determines an elastic response from the original material model, and the second one uses an additional elastic model. Both approaches are applied on barodesy. The simulations are compared with experimental results and with simulations using hypoplastic models with the original intergranular strain concept.


INTRODUCTION
The development of material models is a major field of research in soil mechanics. One specific issue is modelling the behaviour of soil in cyclic loading, ie, the deformation of soil, which is loaded, unloaded, and reloaded several times. The reloading behaviour differs significantly from the behaviour at first loading. Additionally, any reloading does not simply follow the reverse unloading path. We outline this peculiar behaviour qualitatively in Figure 1, starting with the reloading behaviour in triaxial tests in Figure 1A. The difficulties of the majority of constitutive relations to model cyclic loading are already revealed at the first reloading: with basic elastoplastic models, reloading is purely elastic, so all subsequent cycles of loading and unloading do not cause any additional irreversible deformation ( Figure 1B). With basic hypoplastic models 1,2 and barodesy, 3,4 however, it proves difficult to distinguish the initial loading from subsequent reloadings, and hence too large irreversible deformations are predicted. An accumulation of deformation, the so-called ratcheting effect, is obtained in small stress cycles ( Figure 1C). The problem related to reloading is inherently linked with the so-called small-strain stiffness. FIGURE 1 Loading and reloading in triaxial tests. It is assumed that at points 1 and 3 the same stress T, the same stretching D, and the same void ratio e occur. A, Realistic soil behaviour schematically. After the load reversal at point 3, the stiffness for reloading is increased compared with monotonic loading (1)(2). B, Basic elastoplastic models do not distinguish between unloading and reloading stiffness. C, Basic hypoplastic and barodetic models of the typeT = f(T, D, e) do not distinguish between monotonic loading and reloading and therefore underestimate the reloading stiffness [Colour figure can be viewed at wileyonlinelibrary.com]

FIGURE 2
The role of memory in soil behaviour (schematically): Isotropic compression of soil. For (cyclic) reloading, the stiffness at p 1 is increased compared with that for monotonic loading (dashed line) passing through the void ratio e 1 . Thus, the void ratio e 1 and the stress p 1 are not sufficient to describe both, monotonic loading and the reloading behaviour. Figure slightly modified from Herle 5 [Colour figure can be viewed at wileyonlinelibrary.com] not possible. This shortage is negligible in many cases, but it proves deleterious for reloading and for so-called small-strain problems. 5,6 In Figure 2, the isotropic compression behaviour of soil is shown schematically. For (cyclic) reloading, the stiffness at p 1 is strongly increased compared with that for monotonic loading passing through the same void ratio e 1 (dashed line). Thus, the void ratio e 1 and the stress p 1 are not sufficient to describe both monotonic loading and reloading behaviours. The same conclusion can be drawn if we consider triaxial compression of soil; see Figure 1C. After the load reversal at point 3, the stiffness for reloading is increased compared with that for initial loading (1)(2). It is assumed that at points 1 and 3, the same stress T, the same stretching D, and the same void ratio e occur. Again, the void ratio e and the stress T are not sufficient to describe both, initial loading and the reloading behaviour.
To improve the memory of the soil model, one needs to consider what the independent variables in a constitutive relation of the rate type are. If a model does not include any additional information about the deformation history, it has to be enhanced with special extensions to capture the small-strain behaviour. For elastoplastic models, the small strain overlay model 6 is a possible approach. It uses a strain history surface as memory. To model small-strain effects and thus to overcome the ratcheting effect in hypoplastic models, Niemunis and Herle 7 introduced the intergranular strain concept. In this model, the so-called intergranular strain is introduced as an additional state variable that acts as memory of former deformations.
This concept was extended by several authors, eg, Wegener and Herle, 8 Mašín 9 or Fuentes and Triantafyllidis, 10 to improve its performance and to capture more soil features. Wegener and Herle 8 added an additional material parameter to control the transition between the small strain and the large strain model. This parameter is used to reduce the predicted strain accumulation for cyclic loading. Mašín 9 added an anisotropic model for the linear part of the material model and allowed for a nonlinear dependence of the very smallstrain shear modulus on the mean effective stress. Fuentes and Triantafyllidis 10 introduced the intergranular strain anisotropy (ISA) model, which incorporates a yield surface within FIGURE 3 Simulated oedometric compression tests using hypoplasticity 16 show some shortcomings of the intergranular strain concept for different reloading cycles. Parameters of Hochstetten sand 7 starting at isotropic stress state with p 0 = 100kPa the intergranular strain space. This allowed for improved predictions of the memory effects upon reloading paths, which also improved predictions of the strain accumulation in cyclic tests and reduced the overshooting effects. In its original form, the intergranular strain concept is designed for hypoplastic models. We propose two possible extensions of the intergranular strain concept to make it applicable to other material models of the rate type. For simplicity, we only use the basic intergranular strain concept as introduced by Niemunis and Herle 7 without any modifications; however, the above-mentioned improved versions can also be used with our presented approaches. In the last section, we give a brief outlook about a possible application of the improved versions. We only consider models that are positively homogeneous of degree one in the stretching D, ie, rate-independent models, although the original model could also be used for models considering rate effects as the visco-hypoplastic model. 11 The good performance and acceptance of the intergranular strain model is shown by a wide range of applications. [12][13][14] However, the original intergranular strain concept also has some shortcomings; eg, for different load cycles, using the same set of parameters, effects of ratcheting, and overshooting can appear. These effects are shown in Figure 3 for simulations of different load cycles in oedometric compression tests using hypoplasticity with the original intergranular strain concept. As the concept uses a hypoelastic model where a complete reversibility of stress upon a closed strain cycle is not ensured, the concept is not suitable for high-cyclic problems. Wichtmann et al 15 recently compared the behaviour of three sophisticated models for cyclic loading, including the intergranular strain concept and the ISA model. For all models, they showed advantages and shortcomings under different conditions; eg, hypoplastic models using the intergranular strain concept perform better than the ISA model for the accumulation rate under stress-and strain-controlled undrained cycles, whereas the opposite is true for cyclic mobility effects.
The approaches presented here do not represent an improvement of the original intergranular strain concept. They will serve as an application of the intergranular strain concept to other models, which opens the possibility to improve their small-strain behaviour. Moreover, making the intergranular strain concept accessible for a wider range of constitutive models can be a motivation to create new developments and to finally improve the concept itself.

Symbols and notation
In this article, stress is denoted by T and stretching by D. Compression is defined as negative. Second-order tensors are written in bold capital letters (eg, X), ||X|| = √ tr X 2 is the Euclidean norm of a symmetric tensor X, and tr X is the sum of the diagonal components of X. The superscript 0 marks a normalized tensor, ie, X 0 = X∕||X||. Any symmetric second-order tensor can be written as vector with the principal values X v = [X 1 , X 2 , X 3 ]. We use this to display tensors in figures; however, we do not use the notation X v , as it is implicitly clear that X is shown as a vector in these figures. Bold calligraphic letters denote tensors of the fourth order (eg, ). We use different kinds of tensor operations using the Einstein summation convention. In particular, the indices follow the lexicographic order X ⊗ Y = X i Y kl , X ∶ Y = X i Y i , and  ∶ D =  i kl D kl . We employ the unit tensor of second order I with I i = i and the fourth-order tensor  with  i kl = ik l using the Kronecker delta i .
Stresses are considered as effective ones; the normally used dash is omitted. The objective stress rate is denoted byT. The stretching tensor D is the symmetric part of the velocity gradient. The void ratio e is the ratio of the volume of the voids V p to the volume of the solids V s . The additional state variables of the material models are collected in the vector Q. The effective mean stress is denoted by p = − 1 3 tr T. For a conventional triaxial compression or an oedometric compression test, the axial stress is denoted with T 1 and the radial stress is denoted with T 2 = T 3 . The associated strains are 1 and 2 = 3 . Deviatoric stresses and strains are denoted by q = −(T 1 − T 3 ) and q = − 2 3 ( 1 − 3 ).

INTERGRANULAR STRAIN CONCEPT
Before we present our approaches for the extension of the intergranular strain concept, we shortly outline the original concept and investigate some features that are useful for understanding our proposed extensions. For a detailed derivation of the basic concept, we refer to the original publication. 7 We start with a general formulation of hypoplasticitẙ In this formulation, the terms of the hypoplastic equation are collected in terms linear in D, ie, ( ∶ D), and terms that are not linear in D, ie, N||D||. In hypoplasticity, the tensor N does not depend on D, and the use of ||D|| leads to a nonlinear behaviour. 17 To model a stiffness increase subsequent to a change in deformation direction, an additional tensorial state variable, the so-called intergranular strain tensor , is introduced as memory of former deformations. We need the normalized tensor which can be interpreted as direction in the principal stretching space and the magnitude of the intergranular strain normalized with the maximum intergranular strain magnitude R. The intergranular strain tensor follows an evolution equation which must be integrated together with the evolution Equation (1) of the stress rate. Here, R and r are material parameters. The stress rate is given byT with the material stiffness matrix The scalar variables m T , m R , and are material parameters, so, in total, five additional material parameters are needed. For a monotonic deformation, the maximum intergranular strain is reached, hence = 1. A sharp change of the deformation direction following a continuous monotonic deformation yields the maximum increase of stiffness. In case of a full reversal with D 0 = − 0 , the material behaviour is elastic, as only the linear part of the stiffness matrix is used, which is increased by a scalar factor m R . For an orthogonal change of direction with D 0 ⟂ 0 (or 0 ∶ D 0 = 0) and = 1, Equation (5) with Equation (6) reduces toT = m T  ∶ D. For monotonic continuation of deformation with D 0 = 0 (or 0 ∶ D 0 = 1) and = 1, we get the original material model as given in Equation (1). The scalar factor is used to perform an interpolation with change of the magnitude of . A useful graphical representation to visualize the directional response of a constitutive model are response envelopes. 18 The stress response ΔT for the stretching D is visualized by multiplying the stress rateT with a single time increment Δt as a scaling factor. The distance between the initial state T 0 and the point on the response envelope is a measure of the incremental stiffness; see Figure 4A. In Figure 4B, the intergranular strain concept is used. The directions of D 0 A and 0 are the same, ie, the stiffness is the same compared with that of the original material model. The direction of D 0 B differs from 0 , and thus the stiffness is increased. Figure 5 shows the results for response envelopes simulated with hypoplasticity in the version of von Wolffersdorff 16 using the material parameters for Hostun sand given by Niemunis and Herle. 7 The arrows show the former deformation direction, which is stored in the intergranular strain tensor 0 . The response envelopes for hypoplastic models are ellipses, 18 the shape of which is given by the linear part . The centre of the ellipse is shifted from the initial stress point by the nonlinear part N. A response envelope with its centre at the initial stress point maps elastic behaviour, as the stiffnesses in reversal directions are equal. As small-strain behaviour is assumed to be elastic, for reversal directions ( 0 ∶ D 0 ≤ 0), only the part linear in D is used and scaled by m R and m T , respectively. For directions with 0 ∶ D 0 > 0, the nonlinear part is faded in, to obtain the original material behaviour.
As it is needed for further discussion, we combine Equations (5) and (6), with the following abbreviations: Equation (7) depends on three different tensors A, B, and C. The expression 0 ∶ D 0 = cos serves as interpolation function to consider different angles between the direction of intergranular strain 0 and the direction of actual deformation D 0 ; see Figure 4B. With = 0 ∶ D 0 = cos , and D = ||D||D 0 , we can write We get = 1 for monotonic deformation with 0 = D 0 and = 0 for perpendicular deformation with 0 ⟂ D 0 . For reverse of deformation, = −1, and the second line in Equation (7) is used.
The concept of intergranular strain is designed for hypoplastic models. Its original formulation is strongly based on the structure of hypoplastic models, which consist of terms linear and nonlinear in D as shown in Equation (1). The part linear in D is used to calculate the elastic stiffness within the small strain range. For some constitutive models of the rate typeT = f(T, D, Q), (10) this linear part cannot be explicitly extracted from the model formulation.
To be able to apply the intergranular strain concept to such material models without an explicit linear part, we need to offer a consistent elastic stress response for small strains. This is the motivation for creating a more general approach for an intergranular strain extension that can easily be used with other constitutive models of the rate type.

BARODESY
We use barodesy to show the application of our intergranular strain approaches. Barodesy, introduced by Kolymbas, 19 is a material model that has similarities to hypoplasticity. As in hypoplasticity, the stress rate is formulated as a function of stress, stretching, and void ratio:T = f(T, D, e). Barodesy 4 also includes concepts from critical state soil mechanics as a stress-dependent critical void ratio. In barodesy, critical stress states almost coincide with predictions by Matsuoka-Nakai. 20 Compared with elastoplastic models, barodesy and hypoplasticity do not distinguish between elastic and plastic strain. In contrast to hypoplasticity, barodesy cannot be directly split into parts that are linear and nonlinear in D. The full set of equations for barodesy for clay can be found in the Appendix .

INTERNAL ELASTIC MODEL
The intergranular response is assumed to be elastic. 7 In hypoplasticity, this is modelled with the linear part of the material model. In a more general approach, we need an appropriate stress rate that delivers equal stress rates for reversal stretchingsT ie,T el is an odd function with respect to D. We denote byT el the incrementally elastic stress rate. The stress rate T can be decomposed in an incrementally elastic partT el and a nonelastic partT ne with The incrementally elastic stress rateT el is the odd part ofT with respect to D. In the case of hypoplasticity,T el coincides with  ∶ D, which is the part of Equation (1) that is linear in D. Note, however, that for more complex material models, T el is not necessarily linear in D. In our first approach, the given material model is directly used to compute the elastic stress rate. The resulting model will be called internal elastic model (IEM).

Decomposition of the material model
In Equation (7), the three tensors A, B, and C were identified, where evaluations of the elastic and the nonelastic parts of the constitutive model are necessary. The incrementally elastic stress rate in Equation (11) represents the odd part of the constitutive Equation (10) The odd part is then obtained asT The rest of the stress rate in direction of D is the nonelastic part, which is also the even part of the constitutive relation with respect to D:T The incrementally elastic part of the stress rateT el is the first expression A in Equation (7). In expression C, the nonelastic part of the material modelT ne is used.
Expression B in Equation (7) considers different directions of the actual deformation and the intergranular strain, to perform a smooth transition between the tangential and the reversal stiffness increase. For hypoplasticity, B is the evaluation of the linear part of the constitutive model in the direction of the intergranular strain. For our formulation, the odd part of the constitutive equation in direction of intergranular strain is used. Setting D = 0 in Equation (10) and denotingT = f(T, 0 , Q), we can use the same procedure as before: and furtherT el = This finally provides a possible extension of the intergranular strain concept within the IEM approach. Employing the abbreviations in Equations (8) and (9), we get Note that we are using two different interpolation functions 1 and 2 for more flexibility. Both functions depend on the variable = arccos( 0 ∶ D 0 ). At = 90 • , ie, 0 ∶ D 0 = 0, they must vanish to ensure a smooth transition between the two cases in Equation (18). These functions are 1 = 2 = for hypoplasticity; see Equation (9).
Applying these ideas to the hypoplasticity given in Equation (1), we get the expressions For hypoplasticity, the new concept thus reduces to the original intergranular strain concept. The IEM approach models an increase of stiffness and an elastic response in a small-strain range for a given constitutive model using the incrementally elastic part of the model. For constant stretching, Equation (18) leads to the response of the original constitutive model. We need four evaluations of the material model at each time step and each integration point to calculate the different parts of the model.

Application of the IEM approach to barodesy
We apply the IEM approach on barodesy for clay as an example of its use. 4 As the IEM approach uses the behaviour of the given material model to calculate an elastic response and to increase the material stiffness, it also reproduces all its peculiarities, which must be considered to reproduce a realistic behaviour.
The intergranular strain concept using the IEM approach can be regarded as a blow-up of the original response envelope. Thus, the original shape strongly influences the resulting one. The elliptical shape of response envelopes in hypoplasticity is only an assumption, as the exact shape is not known from experiments. 22,23 Response envelopes of barodesy are roughly elliptical but not smooth for isotropic unloading; see Figure 6. Through the interpolation process, this irregularity can also appear in other parts of the intergranular strain response envelope, where the original response envelope is smooth. This is a result of the directional interpolation of the partsT el andT ne with 1 for 0 ∶ D > 0 in Equation (18). The intergranular strain response with 1 = 2 = cos for > 90 • as shown in Figure 6 is not realistic, because a small change of the direction D 0 to the direction 0 leads to a decrease of stiffness. Henceforth, we call this interpolation scheme ipA.
A simple way to adjust the directional stress response is a change in the interpolation function . In our case, we want to get horizontal tangents in the interpolation functions 1 and 2 for = 0 • and = 90 • to reach a smooth transition between the two cases in Equation (18). This can be done by using and 2 = − 1 . The horizontal tangent for = 90 • , as shown in in Figure 7, now leads to a smooth transition between the parts of the model. This interpolation is denoted by ipB.
Using the ipB interpolation scheme, we get a more realistic shape of the intergranular strain response envelope for barodesy, which is shown by the black solid lines in Figure 6.
For other material models, other interpolation functions 1 and 2 may be more suitable. This can be investigated by simulation of response envelopes for different and different initial stress states.

EXTERNAL ELASTIC MODEL
Another approach to get an elastic response for intergranular strain is the following: We define an external elastic model (EEM) and perform an interpolation between the original constitutive model and the external model using the framework of the original intergranular strain concept. Here, as an example, we use an elastic model with the stiffness matrix with the elastic shear modulus G and the Poisson ratio .

FIGURE 6
Response envelopes using barodesy with the parameters in Table 1 and the internal elastic model (IEM) approach for different initial stress states T 0 and different former deformations with = 1. ipA: The arrows show the intergranular strain direction 0 [Colour figure can be viewed at wileyonlinelibrary.com] FIGURE 7 Different interpolation functions 1 and 2 used for Equation (18) and the simulations in Figure 6. ipA: 1 = 2 = cos ; ipB: 1 = − 2 = cos 2 We start with the original intergranular strain concept given in Equation (7). The model response should be elastic within the small-strain range. For hypoplasticity, the elastic response is achieved with the part linear in D. We replace the linear part  with the elastic stiffness matrix  e . Hence, for a reversal of the deformation, only the elastic model is used and scaled up with the factor m R . With ongoing constant stretching, the enlarged stress rate has to reduce to the stress rate of the original model, denoted byT m . We can thus replace the nonlinear part N||D|| with the stress rate of the constitutive modelT m minus the stress rate resulting from the elastic model. A smooth transition between the different models is guaranteed by using the interpolation framework of the original intergranular strain concept. With these assumptions, we get a formulation for the objective stress rateT employing the abbreviations from Equation (8) and the interpolation function = 1 = 2 = cos :

Determination of the elastic model
The main task in this approach is determining the elastic model parameters -in this case, G and from Equation (21) which should be consistent with the original model, and finding a proper calibration. Only if the elastic model is consistent with the original model, a proper interpolation and calibration of the intergranular strain parameters is possible. We can determine an incremental bulk modulus K and a shear modulus G from the original material model by calculating the model response for isotropic and isochoric deformation. For this hypoelastic model, we use the assumption that the intergranular strain behaviour is independent from the stress ratio, so we can perform these calculations at an isotropic stress state.
Calculating the stress rate for isotropic loading, the incremental bulk modulus can be obtained as with the time rate of the effective mean stress . p and the volumetric stretching tr D. The same can be done for isotropic unloading, resulting in K iu . We take then K = 1 2 (K il +K iu ), which is an average value of incremental loading and unloading stiffness from the original model.
The shear modulus G links the deviatoric strain with the deviatoric stress and can be estimated as stiffness for deformations with constant volume. We need to calculate the stiffness for an isochoric deformation, ie, tr D = 0. For an isotropic material, the shear stiffness is the same in all (isochoric) deviatoric directions. Thus, one calculation is enough to assess the shear modulus, as it can be adjusted by calibration of the stiffness parameters m R and m T to fit the initial small-strain shear stiffness G 0 . For a calculation under axially symmetric conditions and in principal stresses, we get The elastic stiffness matrix can be written as Note that Equations (21) and (25) are the same with K = 2G(1+ ) 3(1−2 ) . In the general case, we need three evaluations of the material model to get the elastic parameters and one evaluation to get the material response considering intergranular strain at the actual state. So again, we need four evaluations of the material model at each time step and each integration point, if the parameters cannot be extracted analytically from the model formulation.

Calibration of the EEM model
For a full reversal of the deformation direction, the stiffness should be increased by the scalar factor m R . For a calibration of the stiffness increase, the pressure-dependent small-strain shear stiffness G 0 is used, which is known from laboratory tests, eg, Mašín 24 or estimated by empirical relations, eg, Hardin and Black 25 and Biarez and Hicher. 26 Using the formulation of the elastic stiffness matrix from Equation (21) with = 3K−2G 6K+2G , the calibration procedure presented by Mašín 21 for hypoplasticity is applicable. We get  e = 2G (26) with After a full reversal of deformation direction with = 1, Equation (22) yields With Equation (26), we get For the calibration of m R , we just need to calculate the incremental elastic shear modulus G for the same pressure, for which the small-strain shear modulus G 0 is known.

Application of the EEM approach to barodesy
We apply the EEM approach to barodesy for clay as an example. By predefining the respective boundary conditions, we can directly calculate the incremental elastic parameters at an isotropic state for barodesy. This reduces computational effort as not the entire determination of the elastic parameters has to be done at each time step.
As for hypoplasticity, the response envelope of barodesy at the isotropic state does not change its shape by changing the overconsolidation ratio (OCR). A different OCR only leads to a shift of the centre of the response envelope away from the initial stress point T 0 . In barodesy, the incremental shear modulus and the average bulk modulus stay the same for a different OCR. At normally consolidated isotropic states in barodesy, the bulk modulus for isotropic loading is directly given with the parameter * as with K * il = 1 * . For isotropic unloading, we use Equation (A1) inserting T 0 = R 0 = − 1 √ 3 I and tr D 0 = √ 3. The stress rate is only dependent on the actual mean pressure p and the volumetric strain rate for isotropic unloading with 3||D||. The parameters c 3 and c 4 are constants that can be calculated from the critical state material parameters; see Appendix . We thus get the bulk modulus for isotropic unloading withT = − . p · I.
Using the average value K * = 1 2 (K * il + K * iu ) for the linear elastic model (25), we get The incremental average bulk modulus K = pK * only depends on the actual mean pressure. The shear modulus G can be calculated with an isochoric deformation D 11 = −2D 22 , D 22 = D 33 (here, tr D 0 = 0) and Equation (24):

FIGURE 8
Response envelopes using barodesy with the parameters in Table 1  3p as and The shear and bulk moduli for barodesy for clay at the isotropic state only depend on the actual mean pressure. The stiffness coefficients K * and G * only depend on material parameters and can be used to calculate the linear elastic material stiffness matrix. Thus, Poisson ratio here is also independent of the pressure and can be calculated as Using this calibration, we obtain a hypoelastic model with a pressure-dependent shear modulus and a constant Poisson ratio. The stiffness matrix according to Equation (26) reads with the stiffness factor e = p. Figure 8 shows response envelopes simulated with barodesy using the EEM approach. The dashed lines show the original response envelopes and the grey solid lines the responses from the linear elastic model centred at the initial stress state. The intergranular strain response is a combination of both model responses using the interpolation in Equation (22).

COMPARISON AND SIMULATIONS
In this section, we show and discuss some simulation results using both proposed approaches in combination with hypoplasticity in the versions of von Wolffersdorff (v.W.) and Mašín 21 (clay hypoplasticity) as well as with barodesy. Those results are compared with experimental data and with simulations using the original intergranular strain concept.

Simulation results
In a first step, we apply both approaches to hypoplasticity. The response envelopes using v.W. hypoplasticity for different intergranular strain directions and different initial stress states are shown in Figure 9. The IEM approach in combination with hypoplasticity yields exactly the same results as the original intergranular strain concept. The linear part in the v.W. hypoplastic model is dependent on the actual stress state, while the elastic model used for the EEM only depends on the mean stress and is independent of the deviatoric stress. This leads to different orientations of the response envelopes for nonisotropic stress states. In clay hypoplasticity, the linear part is also not effected by the deviatoric stresses. 21 Combining this model with the IEM or the EEM approach leads to the same results as with the original intergranular strain concept.
We use triaxial tests on reconstituted London clay from Mašín, 24 namely, PhM17, PhM19, and PhM21, for a comparison of the simulations with experimental data. A discussion of these tests is done by Mašín. 27 We use barodesy for clay and clay hypoplasticity for the simulations. The employed parameters for London clay are given in Table 2. The boundary conditions for the simulations are shown in Table 3. The simulations start with a reverse intergranular strain direction.
For the simulations with hypoplasticity, we use m R = 4.5. For the EEM approach, G * 0 = G 0 p can be assessed using linear regression as shown by Mašín. 27 We get m R,EEM = 4.0 for barodesy using Equations (37) and (30) to fit the hypoplastic calibration of test PhM19. 27 As the IEM approach uses the response of the actual material model, the parameters m R and m T may differ from those of the EEM approach. We use m R,IEM = 4.3 for barodesy to get an equal stiffness for test PhM19. For all simulations, m T = 0.5m R . The other intergranular strain parameters are the default values from Mašín, 27 which are r = 0.2, R = 10 −4 and = 6. Figure 10 shows the degradation of the shear stiffness for the tests PhM17, PhM19, and PhM21. As the shear stiffness is calibrated for test PhM19, all models show the same results within the small-strain range. Using clay hypoplasticity, the proposed approaches show the same results as the original intergranular strain concept ( Figure 10A.) With a change of the pressure and the stress states, the results of barodesy differ slightly. For test PhM21, the initial pressure of 450kPa was quite different to the pressure of 110kPa used for the calibration with test PhM19. According to Mašín,27 this may explain the difference between the simulations and the experimental results. Overall, however, with both approaches, barodesy is able to model those pressure effects and represent the experimental results in the same quality as hypoplasticity.
Standard laboratory tests are simulated using the proposed extensions of the intergranular strain concept shown in Figures 11 and 12. Both approaches, IEM and EEM, show almost the same results and present a clear improvement compared with the original material model, as the reloading stiffness is increased. The shortcomings of the original intergranular strain formulation also emerge by applying our approaches. The second and the third reloading loops for the isotropic compression test in Figure 12B show ratcheting and overshooting effects, respectively.

Application to the modified Cam-Clay model
The proposed concepts can also be used for elastoplastic models. Basic elastoplastic models with isotropic hardening, such as the modified Cam-Clay model 28 (MCCM), behave elastically for unloading and reloading. Such models do not show ratcheting effects as hypoplastic or barodetic models, but they fail to model the hysteretic behaviour that occurs due to reloading, as shown in Figure 2. The here proposed small-strain extensions yield an increase of the stiffness for small strains and the hysteretic loops as shown in Figure 1.
The MCCM, enhanced with our intergranular strain approaches, is used for the simulations shown in Figure 13. The material parameters for MCCM are M = 0.88, v = 3.49, = 0.24, = 0.044, and = 0.3. The parameters for the elastic stiffness matrix for the EEM approach can be calculated as explained in Section 5.1 or directly be taken from the model p . For the normally consolidated simulations, it holds that p 0 = p ini . The parameters for the intergranular strain extension are set to the same values as for hypoplasticity. The IEM approach and the EEM approach show here virtually the same results.

General remarks
Small-strain stiffness can be modelled with the IEM approach by increasing the elastic stiffness that is obtained using the original model. As an advantage of the IEM approach compared with the EEM approach, no additional model, including a calibration of additional model parameters, is necessary. The reproduced model behaviour is directly given by the original material model. As shown, however, using the original model response can lead to unrealistic responses because of the directional interpolation. This can be corrected by adjusting the interpolation function for a specific material model. For hypoplastic models using the formulation of Equation (1), the IEM approach leads exactly to the original intergranular strain concept; see Equation (19).
With the EEM approach, an additional external model is added, which is used to model the elastic behaviour in small-strain range. This can be a simple elastic model, as presented in this article, but it can also be any other hypoelastic model. 7 Using the EEM approach, it is possible to separate the small-strain behaviour from the large-strain model. For example, concepts of anisotropy as proposed by Mašín and Rott 29 could easily be incorporated for the small-strain range by changing the external model.
Both proposed approaches are able to use the basic intergranular strain concept as proposed by Niemunis and Herle 7 and also any improved version. For an application of the approach of Wegener and Herle, 8 to reduce the strain accumulation by changing the fade-in of the nonlinear model, only the exponent in the third term of Equations (18) and (22) has to be changed. As the ISA model 10 is also applicable to hypoplastic models using the -N formulation, 30 the possibility of applying the ISA model to other material models is provided by the approaches presented here. However, it has to be mentioned that the shortcomings of the original intergranular strain concept are also reproduced by the here presented approaches; see Figure 12.
Both approaches need in general four calls of the material model in each time step, so they need essentially the same computing time. The elastic parameters for the EEM approach can for some models directly be calculated, as shown for barodesy and the MCCM, so the same computational effort as for the original intergranular strain concept is needed.

SUMMARY
Two different approaches for extending the intergranular strain concept to nonhypoplastic models are presented in this article. The first approach -IEM -uses the properties of the original constitutive model and calculates an elastic response, which is used for the stiffness increase in the small-strain range. For the second approach -EEM -an external elastic model is used for the calculation of the small-strain stiffness, and an interpolation between the original model and the external model is performed. For both approaches, an application to a specific material model -barodesy -is shown. For the simulated standard laboratory tests, both approaches show almost the same results and are an improvement compared with the original model behaviour. Both approaches also work for elastoplastic material models, which is demonstrated for the MCCM. The original intergranular strain concept is a well-established concept to improve small-strain behaviour in hypoplastic models, but it still also has its shortcomings. However, the possibility to apply the intergranular strain concept to more constitutive models will make it available for a wider range of constitutive modellers and may also help to develop improvements for the concept.