Visco Drucker-Prager Model: Triaxial Driver versus Semi-Analytical Solution

Problem description

This example uses the Triaxial Driver to simulate a triaxial compression test of a Visco Drucker-Prager solid. Constant lateral confining stress together with loading/unloading axial strain periods are imposed. Imposed axial strain range are high enough for allowing visco-plastic yield in both loading and unloading period. This complicated scenario is used for verifying the numerical convergence and accuracy of the Visco Drucker-Prager constitutive model implemented in GEOS.

Semi analytical result for axial stress variation \Delta\sigma_{V} and lateral strain variation \Delta\varepsilon_{V} can be established for the imposed triaxial boundary conditions following the theoretical basis of the Perzyna time dependent approach presented by (Runesson et al. 1999) as:

\Delta\sigma_{V} = (\Delta\varepsilon_{V}-\Delta\lambda\frac{b^\prime-3}{3}) E

\Delta\varepsilon_{H} = \Delta\varepsilon_{V} - \frac{\Delta\sigma_{V}}{2\mu} + \frac{3}{2}\Delta\lambda

where E and \mu are the elastic Young and shear moduli. The visco-plastic multiplier \Delta\lambda can be approximated by:

\Delta\lambda = \frac{\Delta t}{t_*} \frac{F}{3\mu + Kbb^\prime + h}

in which \Delta t is the time increment, t_* is the relaxation time, F is the stress function defining the visco-plastic yield surface, K is the elastic bulk modulus, b is the frictional parameter defining the visco-plastic yield surface, b^\prime is the dilation parameter defining the plastic potential and h is the hardening rate. These solutions are applied only when plastic yield condition is satisfied. The cohesion parameter defining the plastic yield surface is updated with stress change as

\Delta a = h \Delta\lambda

These solutions were established for a positive shear stress q = -(\sigma_{V} - \sigma_{H}) (negative sign convention for compression stress). For the case when the plastic yield occurs at a negative shear stress, we have

\Delta\sigma_{V} = (\Delta\varepsilon_{V}-\Delta\lambda\frac{b^\prime+3}{3}) E

\Delta\varepsilon_{H} = \Delta\varepsilon_{V} - \frac{\Delta\sigma_{V}}{2\mu} - \frac{3}{2}\Delta\lambda

These solutions are implemented in a Python script associated to this example for verifying GEOS results.

Input files

This benchmark example uses two GEOS xml files that are located at:




It also uses a set of table files located at:


A Python script for the semi-analytical solutions presented above as well as for post-processing the GEOS results is provided at:


For this example, we focus on the Task and the Constitutive tags.


The imposed axial strain loading/unloading periods, the constant lateral confining stress as well as the initial stress are defined in the Task block as

      output="ViscoDruckerPragerResults.txt" />

Constitutive laws

The elasto-visco-plastic parameters are defined as


All constitutive parameters such as density, viscosity, and bulk and shear moduli are specified in the International System of Units.

A comparison between GEOS results and semi-analytical results

The simulation results are saved in a text file, named ViscoDruckerPragerResults.txt. A comparison between the results given by the TriaxialDriver solver in GEOS and the approximated semi-analytical results presented above is shown below. Interestingly we observed that the Duvaut-Lions approach implemented in GEOS can fit perfectly with the Perzyna approach that was considered for deriving the analytical results. This consistency between these time dependence approaches is because of the linear hardening law of the considered constitutive model as already discussed by (Runesson et al. 1999) .

(Source code)


To go further

Feedback on this example

For any feedback on this example, please submit a GitHub issue on the project’s GitHub page.