Simulations of Fine-Meshed Biaxial Tests with Barodesy

Recent experimental studies showed that shear band development starts at the beginning of triaxial tests. In experimental testing, it is impossible to obtain a soil sample with a homogeneous void ratio. Therefore, a homogeneous deformation, i.e., an element test, is questionable well before the peak. In this article we carry out finite element simulations of fine-meshed biaxial tests with the constitutive model barodesy, where the stress rate is formulated as a function of stress, stretching and void ratio. The initial void ratio in the simulations is normally distributed over all elements in a narrow range. In this article, we evaluate the pre-peak shear band development. We further compare stress paths and stress-strain curves of the biaxial test of relevant elements (e.g., inand outside the shear band) with the results of the average response of all elements. We show how the response in an element test differs from the average response of the fine-meshed test. We present the resulting potential for understanding (early) shear band development and stress-strain behaviour in a biaxial test: The inhomogeneous void ratio distribution in a sample favours early shear band development. This effect is modelled with barodesy. The obtained stress paths and stress-strain curves show that the maximum deviatoric stress is higher in the element test than it is in the average response of the fine-meshed test.


Introduction
Soil is a material which consists of different grain sizes and grain shapes.The density or void ratio in a soil sample is naturally scattered and is thus not constant.Even for an ideal material consisting of spheres of the same size, inhomogeneities are to be expected due to inhomogeneous sample mounting.Experiments, e.g., triaxial tests, are carried out to investigate soil behaviour.In the evaluation of stress and strain, it is common to assume the sample to be a single homogeneous element with rectilinear deformation.In experimental testing, however, it is impossible to obtain a completely homogeneous sample with constant void ratio.In addition, the deformation of the sample in the experiment is not rectilinear: lateral barreling of the sample and/or development of shear band(s) are observed.
Finite element simulations on homogeneous samples lead to homogeneous deformations and not to development of shear bands.However, the pronounced localization in a biaxial test results in different stress-strain curves in-and outside the shear band.Thus, a single element test obviously cannot reproduce a realistic stress-strain response, e.g., [1].To obtain shear bands, at least one weak element must be included in the simulation, e.g., [2][3][4].Tomographic investigations on experiments by Desrues et al. [5] show that shear bands develop from the very beginning of triaxial tests.Therefore, a homogeneous deformation, i.e., an element test, is questionable well before the peak.
We perform finite element simulations of fine-meshed biaxial tests with barodesy [6].Barodesy is a material model which includes the void ratio as a state variable.In this article, the initial void ratio is normally distributed over all elements in a narrow range, comparable to [7].Thus, the natural scattering of density is modelled.The randomly distributed void ratio leads to a scatter of dilatancy and therefore in a scatter of peak strength in the elements.Hence, patterns of shear bands develop from the very beginning of the biaxial test.We evaluate stress-strain curves, stress paths and the pattern of the shear bands with ongoing deviatoric strain [8].
The innovative part of this research is to scatter the void ratio in combination with an evaluation of the average response of stress and strain of the fine-meshed test.The aim of this work is to find out whether and in which cases the approximation of a soil sample by a single element is justified and valid.In this context, it is possible to distinguish what has to be modelled by the constitutive relation and what by the numerical approach.We therefore compare the different response of a single element test with the average response of a fine-meshed test.

Constitutive Model-Barodesy
Barodesy is a material model that differs from conventional elasto-plastic models: it does not distinguish between elastic and plastic deformations and does not require expressions such as yield function, plastic potential or flow rule.The stress rate is formulated as a tensorial function of the current stress, stretching and the void ratio, i.e., σ = h(σ, D, e).Thus, barodesy has certain similarities with hypoplasticity [9,10].Appendix A summarizes all equations of barodesy.

Critical State Soil Mechanics
The following concepts from Critical State Soil Mechanics are included in the mathematical formulation of barodesy [6]:

•
A critical stress state locus [11,12], which is a line, the Critical State Line (CSL), in the p -q plot: with M = 6 sin ϕ c /(3 − sin ϕ c ) for axisymmetric triaxial compression, ϕ c is the critical friction angle.In three dimensional stress space, the critical stress locus of barodesy practically coincides with the Matsuoka-Nakai failure criterion [11], see Figure 1.
• A stress dependent critical void ratio e c (CSL from [13]) in the p -e plot, in order to distinguish between normally to slightly overconsolidated soil (e > e c ) and highly overconsolidated soil (e < e c ):

•
The isotropic normal compression line (NCL) from [14], in order to define the loosest possible packing: σ * = 1 kPa is a reference stress, N and λ * are soil parameters from Critical State Soil Mechanics.

Barotropy and Pyknotropy
By including the actual void ratio e as a state variable in barodesy, effects such as pyknotropy (i.e., the dependence of stiffness and strength on density) and barotropy (i.e., the dependence of stiffness and strength on stress level) can be modelled.In Figure 2 consolidated drained (CD) triaxial tests with Weald clay (parameters in Table 1) with the initial void ratio e ini = const are simulated with barodesy.Highly overconsolidated samples (grayed area) dilate, normally and slightly overconsolidated samples contract, see Figure 2a.The stress states and the void ratios then reach the critical state.Highly overconsolidated samples achieve higher peak friction angles than the critical friction angle.The stress states reach the CSL after softening.
1 (b) mean stress -deviatoric stress plot: highly overconsolidated samples reach peak states which are higher than critical strength.(c) axial strain -stress ratio plot: highly overconsolidated samples reach higher stress ratios/mobilized friction angles than the critical friction angle.Similar results are obtained if the mean stress p is kept constant and the void ratio is varied: In Figure 3 consolidated drained (CD) triaxial tests with Weald clay with p ini = const are simulated with barodesy.The highly overconsolidated samples dilate, the normally and slightly overconsolidated samples contract, see Figure 3b.The lower the initial void ratio is, the higher is the maximum deviatoric stress q, see Figure 3b.Again the highly overconsolidated samples reach peak friction angles that are larger than the critical friction angle of 24 • , see Figure 3c.q (kPa) 1 stress plot: the lower the initial void ratio is, the higher is peak strength.(c) axial strain -stress ratio plot: highly overconsolidated samples reach higher stress ratios/mobilized friction angles than the critical friction angle.

Biaxial Tests
We simulate drained plane strain triaxial compression (biaxial) tests, i.e., ε 2 = 0 and σ 3 = 100 kPa = const.Biaxial tests on Weald clay with different overconsolidation ratios are investigated.The parameters of Weald clay are in Table 1.It is common to define the overconsolidation ratio (OCR = p e /p ) by means of the so-called Hvorslev's equivalent consolidation pressure p e = exp N−ln(1+e) λ * , divided by the actual mean stress p .p e is the value of mean stress on the isotropic normal consolidation line which refers to the current specific volume (1 + e).

Initial Conditions
We simulate biaxial tests, with a width of 37 mm and a height of 74 mm, see Figure 4a.Multiple elements and the stress boundary condition (σ 3 = const) allow barreling of the sample in direction 3. The shear band thickness in finite element simulations is mesh-dependent, as the mesh serves as an internal length parameter, e.g., [20].As the stress-strain curves in-and outside the shear band differ, the thickness certainly has an effect on the average response, see [21,22] and Figure 4b,c.
inside shear band A outside shear band B single element test ε q (%) Final shear band at global deviatoric strain of ε q = 8.8%, displayed as deviatoric strain, green: ε q = 22% to red: ε q = 53% (c) Stress-strain curves in-A and outside B the shear band, single element test.
Gylland et al. [23] report about shear band thicknesses in the range of 2 to 8 mm for clay.In our simulations, the number elements (5000) have been chosen in a way that the final shear band thickness is in the reported range.In Figure 4b the final shear band is about 3 mm.By determining the preferred shear band thickness, we expect realistic predictions in the average response of strain and strain.This approach is only applicable for small dimensions due to the large number of elements.To overcome the mesh-dependency, a regularization technique is required [2,20,21].
The simulations in this work are carried out with ABAQUS 6.14 [24].Four-node plane strain elements are used, see Figure 4b.The UMAT of barodesy [6] is available at SoilModels [25].Finite element applications with barodesy can be found in [4,26,27].
In our simulations the initial void ratio is randomly and normally distributed over all 5000 elements in a narrow range (The standard deviation in our simulations is 0.005.In the simulations by Schneider-Muntau et al. [3], one weak element with an increased void ratio of 0.02 was included to obtain a shear band.) for different overconsolidation ratios according to Table 2.For an OCR of 6, the void ratio distribution is shown in Figure 5a.The mean value of ēini = 0.526 is marked with a black line.The red lines show the standard deviation, i.e., ēini ± 0.005.In order to better estimate the size of this standard deviation, we mark ēini and ēini ± 0.005 in the mean stress-void ratio plot with black and red lines in Figure 5b.

Evaluation of Stress Paths in the Deviatoric Plane
The deviatoric direction of a stress path is indicated by the Lode angle α σ , Equation (4): α σ = 30 • refers to triaxial compression, α σ = −30 • refers to triaxial extension.Nakai [28] reports that stress paths of plane strain conditions lie within 0 • < α σ < 15 • .The experimental findings by Nakai have been carried out on normally consolidated Fujinomori clay.However, the values also provide an indication for overconsolidated samples.Investigations under plane strain conditions were carried out with barodesy in [12,26].

Evaluation of the Inclination of Final Shear Band
The higher the overconsolidation ratio is, the steeper is the inclination of the final shear band.In barodesy, for a higher overconsolidation ratio there is a higher dilatancy and thus a higher mobilized friction angle, see Figure 3.The angle of dilatancy ψ under plane strain conditions (e.g., [29]) is The orientation of the shear band is reported (e.g., [23,30]) to lie between Roscoe's solution [31] with and the Mohr-Coulomb solution with Arthur et al. [32] propose based on their experimental results.The angle ψ in Equations ( 6) and ( 7) is the angle of dilatancy at the peak of q. ϕ max m is the maximum mobilized friction angle, which is for the simulations in this article at the peak of q.

Results
In this Section we evaluate the simulations of fine-meshed biaxial tests.The results of the overconsolidation ratios 1.5, 3, 6 and 54 are discussed.We focus on early shear band development, stress-strain behaviour, stress paths in the deviatoric plane and inclination of the final shear band.All simulations are carried out with the material model barodesy [6].Further details concerning the model are in Section 2.1.Information about the initial conditions and evaluation of the biaxial tests are in Section 2.2.

Early Shear Band Development
In Figures 6-9 we show the development of shear bands at different global deviatoric strain ε q .Figure 6 corresponds to an OCR of 1.5, Figure 7 to an OCR of 3, Figure 8 to an OCR of 6 and Figure 9 to an OCR of 54.For all overconsolidation ratios, it is visible that shear bands develop well before the peak.The peak state is marked in the figures or figures captions.For example, Figure 6a-c show pre-peak pattern of shear bands, (d) is just at peak state and (e,f) show the final shear band.As expected, the smaller the OCR, the more pronounced is the barreling of the sample.For all OCRs a final shear band develops.The development of the final shear band has a strong influence on the shape of the average response of the stress-strain curve, especially in the post-peak range.
In Figures 10-13 the result of an element test (blue), as well as the average response of the fine-meshed test (red) are shown for the OCRs of 1.5, 3, 6 and 54.The average densities of the fine-meshed biaxial test and the single element test coincide for each OCR.The global deviatoric strain corresponding to Figures 6-9 is marked in Figures 10-13.In Figure 10, the average response of the fine-meshed test (red) and the element test (blue) give similar results until ε q ≈ 8%.An overconsolidation of 1.5 corresponds to a slightly overconsolidated specimen.Therefore, the element test does not show any softening.It is interesting to note that the average response of the fine-meshed test has a peak with subsequent softening.Biaxial tests on loose Hostun sand samples [33] show a similar stress-strain behaviour.For all OCRs, the peak and subsequent softening is more pronounced in the average response of the fine-meshed test than it is in the element test.
In the element test, the maximum deviatoric stress is higher than in the fine-meshed test.The lower the initial overconsolidation ratio is, the larger is the difference, see Table 3.For the OCR of 1.5, the difference is 9.64% and for the OCR of 54, 0.37%.In Table 3 the maximum mobilized friction angles ϕ max m in the element tests and in the average response of all elements are added for the different overconsolidation ratios.
--element test --average response of all elements (a) The global deviatoric strain ε q marked with (a-f) corresponds to Figure 7a-f.
--element test --average response of all elements ε q (%) q (kPa) After the peak, the average response of the fine-meshed test shows a pronounced decrease in deviatoric stress compared to the single element test.The pronounced decrease results from the varying stress-strain curves in the individual elements, during softening, see Figures 14a-17a.The stress-strain curves that are close to the element test are inside the shear band.The stress-strain curves that are similar to unloading curves are outside the final shear band, see also Figure 4c and [22].The average response shows a typical stress-strain curve of a biaxial test (e.g., [33]).Figures 14b-17b show the volumetric behaviour.Again the average response differs from the single element test.
--element test --all elements --average response of all elements --locus of critical stress states • peak in q

Stress Paths in the Deviatoric Plane
In Figures 14c-17c the Lode angles according to Equation ( 4) α σ = 0 • , 15 • and 30 • are marked.At failure, stress paths of the average response of all elements (red line) as well as the stress paths of the element test in the simulations do not leave the range of 0 • < α σ < 15 • for the OCR = 1.5, 3 and 6.The experimental findings by Nakai [28], that the Lode angles under plane strain compression are within 0 • < α σ < 15 • have been carried out on normally consolidated samples.The element test in Figure 14c reaches the locus of critical stress states of barodesy, whereas the average response reached a lower maximum mobilized friction angle which is inside the locus of critical stress states, see Figure 14d.
--element test --all elements --average response of all elements --locus of critical stress states • peak in q For the simulation with OCR = 54, the increase of σ 2 is larger than in the simulations with lower overconsolidation ratios.The higher the OCR is, the higher is the dilatancy.By keeping the lateral expansion to zero (ε 2 = 0), higher stresses σ 2 develop for the OCR of 54 in Figure 17c,d than for OCR of 1.5, 3 and 6, see Figures 14-16c,d.
As expected, the highly overconsolidated samples in Figures 16 and 17c,d include peak states which lie outside the locus of critical stress states; see also Figure 2b.

Inclination of the Final Shear Band
The shear band inclination is evaluated for the tests with ψ > 0 • , that are the highly overconsolidated samples with OCR 6 and 54.The orientations of the shear bands are estimated according to Equations ( 6)- (8), see Table 4.In the finite element simulations, the inclination of the shear band is evaluated just after the peak as soon as the inclination of the shear band is clearly visible, see Figures 8e and 9d.
For an OCR of 6, the measured angle of θ = 49.5 • lies between Roscoe's solution (θ R = 49.4 --element test --all elements --average response of all elements --locus of critical stress states • peak in q ε q (%) --element test --all elements --average response of all elements • peak in q --locus of critical stress states

Discussion and Conclusions
Soil is an inhomogeneous material as it consists of different grain sizes and grain shapes.A constant void ratio distribution over the entire sample is therefore not possible.Inhomogeneities in the samples density are also generated due to inhomogeneous sample mounting.

•
Inhomogeneities in the samples density may be one cause for early shear bands development.
Finite element simulations of fine-meshed biaxial test with the material model barodesy show this effect: In barodesy, the density is considered with the void ratio as state variable.A randomly distributed void ratio results in a scatter of dilatancy.Thus, pre-peak patterns of shear bands develop from the very beginning of the biaxial tests.

•
The stress-strain behaviour of a biaxial test-in particular the post-peak behaviour-cannot be modelled with a single element test.In our simulations, the number of elements has been chosen in a way that the final shear band thickness is in the experimental confirmed range.A realistic estimation of the shear band thickness results in a realistic average response of the stress-strain behaviour, due to the pronounced localization.

•
The stress paths in the deviatoric plane show reasonable results: The Lode angles of the stress paths for the OCRs 1.5, 3 and 6 are in the range of 0 • < α σ < 15 • , which is experimentally confirmed for plane strain failure.

•
The maximum deviatoric stress in the element test is higher than the maximum deviatoric stress of the average response of the fine-meshed test.The lower the initial overconsolidation ratio is, the larger is the deviation in peak strength.

•
The inclination of the final shear band is dependent on the overconsolidation ratio: The higher the overconsolidation ratio is, the larger is the inclination of the final shear band.

Figure 2 .
Figure 2. Barotropy in barodesy, figure modified from [19]: simulations with barodesy of CD triaxial tests with Weald clay with e ini = const.(a) mean stress -void ratio plot: highly overconsolidated samples (grayed area) dilate, normally and slightly overconsolidated samples contract under shearing.(b)mean stress -deviatoric stress plot: highly overconsolidated samples reach peak states which are higher than critical strength.(c) axial strain -stress ratio plot: highly overconsolidated samples reach higher stress ratios/mobilized friction angles than the critical friction angle.

Figure 3 .
Figure 3. Pyknotropy in barodesy: simulations with barodesy of CD triaxial tests with Weald clay with p ini = const.(a)mean stress -void ratio plot: highly overconsolidated samples (grayed area) dilate, normally and slightly overconsolidated samples contract under shearing.(b) axial strain -deviatoric stress plot: the lower the initial void ratio is, the higher is peak strength.(c) axial strain -stress ratio plot: highly overconsolidated samples reach higher stress ratios/mobilized friction angles than the critical friction angle.

1 Figure 5 .
Figure 5.For OCR = 6 follows the mean value ēini = 0.526, black line.The standard deviation is 0.005.The red lines mark the standard deviation, i.e., ēini ± 0.005.(a) The initial void ratio is normally distributed over all 5000 elements.(b) shows the range in the e-p plot with ēini and ēini ± 0.005.

Table 1 .
Critical state soil mechanics parameters used for the calibration of barodesy.

Table 2 .
Initial conditions for the FE simulations.

Table 3 .
Maximum deviatoric stress q max in the element test and in the average response of all elements and maximum mobilized friction angles ϕ max m in the element test and in the average response of all elements.

q max Difference in q ϕ max m Difference in ϕ m
--element test --average response of all elements