Skip to main content

The well-known Hardening Soil Model is used worldwide to model geotechnical engineering problems, such as deep excavations and embankment construction, because it captures several features of real soil behavior for both sandy soils and clays and silts.

This webinar briefly introduces the Hardening Soil Model theory and formulation. The creation of a SIGMA/W analysis using the Hardening Soil Model to simulate the triaxial test is included. You will be taken through a step-by-step procedure to calibrate the model and identify the model constants based on drained triaxial test results. The simulated results are compared to the laboratory data.



Curtis Kelln


34 min

See more on demand videos


Find out more about Seequent's civil solution

Learn more

Video Transcript:

<v Curtis>Welcome to this GeoStudio</v>

Material Model webinar series

on the Hardening Soil Model.

My name is Curtis Kelln,

I’m the Director of Research and Development

with the GeoSlope Engineering team at Seequent.

Today’s webinar will be approximately 30 minutes long.

Attendees can ask questions using the chat feature

and I will respond to these questions

by email as quickly as possible.

A recording of the webinar will be available,

so that you can review the demonstration, at a later time.

GeoStudio was developed for geotechnical engineers

and earth scientists.

As such, it comprises a range of products,

that can be used independently

or together to solve a wide array of problems.

Today’s webinar will focus on SIGMA/W.

Those looking to learn more about the products,

including background theory, available features,

and typical modeling scenarios,

can find an extensive library of resources

on the Geoslope website.

Here, you can find a tutorial videos,

examples with detailed explanations

and reference books on heat and mass transfer modeling,

statics stress rain modeling, slope stability modeling,

and dynamics stress strain modeling with GeoStudio.

To begin this webinar, let us consider a soil specimen

subject to an experimental procedure,

and a numerical simulation

of that experimental procedure.

And in experimental approach,

there’s a real soil sample,

from which specimens can be extracted

and subjected to various tests,

such as triaxial and odometer loading.

The data collected from these tests is then interpreted

and the interpretation has to be done in a manner

that is specific to a particular constitutive of model.

The parametrization procedure

produces the models constants.

These constants are then used as inputs for the model.

The simulated results are obtained

by solving the physics subject to boundary conditions

and a constitutive law.

And finally, the simulated results are compared

to the corresponding laboratory results

to verify the parameterization.

In today’s webinar, we will review the key elements

of the hardening soil constitutive model.

Then, a step-by-step procedure will be presented

for parameterizing the model.

The procedure will then be applied to data

from numerous triaxial tests on sandy soil.

A demonstration of the steps from modeling a triaxial test

and SIGMA will be provided,

and the simulated and measured results will be compared.

Some concluding remarks will also be provided

at the end of the webinar.

The Hardening Soil Model was proposed by Schanz,

Vermeer and Bonnier in 1999.

It was formulated on the framework

of the Classical Theory of Plasticity.

The formulation comprises three stress-dependent

stiffness quantities and two yield surfaces.

One that resembles a hexagonal cone

and a second that resembles the tip of an ellipsoid.

The shear surface predominantly governs the response

to changes in deviatoric stress

while the cap predominantly governs

the response to increments in mean effective stress.

Failure is controlled by a Mohr-Coulomb failure law

hence the hexagonal shape of the yield surface.

The model implemented in Sigma

actually comprises a smoother cap

than what is depicted in this figure.

The reference book should be consulted

to understand the details of the formulation.

In the triaxial stress state,

the total stress path is a straight line

that reaches the Mohr-Coulomb envelope

out a specific failure point.

The deviatoric stress,

axial strain response of the soil

under this loading condition

is estimated by hyperbolic curve,

as will be discussed shortly.

The characteristics of this curve

are determined by three

stress-dependent stiffness quantities,

including the initial tangent stiffness,

E subscript I,

the half failure secant stiffness

E subscript 50

and an unloading reloading stiffness E subscript UR.

In order to have a closed yield surface,

a cap function is defined in the hardening soil model.

This cap has an elliptical shape

with an aspect ratio of M

or the Greek letter mu.

In the PQ stress base,

the stress state of any normally consolidated sample

is located on the cap.

That’s at the tip of this blue stress path

that intersects the ellipse.

The stress strain behavior of the sample

in an odometer test is defined

by a tension stiffness called the odometer stiffness

or E subscript OED.

Like the other stiffness quantities,

the odometer stiffness is not constant,

but it’s value depends on the stress of the specimen.

To review the key ingredients of the hardening soil model,

Let’s consider a theoretical triaxial response.

The axial strain of the sample in this condition

is to defined by a hyperbolic function.

The hyperbolic curve continues

until the Mohr-Coulomb failure stress is reached,

Q subscript F.

The failure occurs before the curve reaches

a asymptotic tonic value.

The ratio of the failure to asymptotic deviatoric stresses

is R subscript F.

The behavior of the sample at the beginning of loading

is modeled by the initial stiffness

E subscript I.

Since this quantity is difficult to measure,

its value is presented in terms of

the Half failure Secant stiffness,

E subscript 50.

The half failure secant stiffness

is a stress dependent stiffness

that depends on the minimum principle stress,

Sigma subscript three.

A similar stress dependency is also considered

in the Hardening soil model

for unloading reloading stiffness.

The model’s constants are now highlighted in gray.

In an odometer test on normally consolidated soil,

the slope of the stress path

is a function of the coefficient of lateral earth pressure,

K nought superscript NC.

Moreover, the cap yield surface and PQ space is an ellipse

with a horizontal major axes.

In a non cohesive soil,

the elliptical yield surface is centered around the origin.

The odometer stiffness is also a stress dependence stiffness

presented in terms of the vertical stress,

Sigma one in a one dimensional odometer tests.

Notice that the dependence on Sigma one

is different from the other moduli

that were dependent on Sigma three.

The model’s constants are once again

highlighted in light gray.

To summarize,

the hardening soil model requires

three stiffness Constants,

a reference stress,

and an exponent M.

The strength constants and failure ratio

are subscript F,

the coefficient of lateral earth pressure.

Some other parameters will also be needed

to describe the initial condition of a specific soil sample.

These include the initial void ratio,

the unit weight of the soil,

and an isotropic over consolidation ratio,

which we’ll touch on briefly in the subsequent slides.

We can now discuss a step-by-step parametrization procedure

for determining the Hardening soil model Constants

from a number of drain triaxial tests.

Step one involves determination

of the effective friction angle

and the effect of cohesion C.

The failure stress from numerous triaxial tests

at different confining pressures

are theoretically located

on the Mohr-Coulomb failure line in PQ space.

As such the slope of the line,

which is 10 beta,

as well as the Y intercept alpha,

is calculated using the method of least squares

where N is the number of data points.

Alternatively, the trendline functionality

and programs like Excel

can be used instead of the least squares formulas.

The effective friction angle and effective cohesion

of the soil can be found by changing the coordinate system

to a shear normal stress base.

The equations corresponding to that transformation

are shown here.

The second step of the procedure

is to consider a reference stress

and to estimate reference stiffnesses.

One of the available triaxial tests

is considered as a reference test.

The confining pressure of this test

is then called the reference stress.

The failure of deviatoric stress for this particular test

is determined based on Mohr-Coulomb strength constants

determined in the previous step.

The half failure referenced stiffness

is the slope of a straight line that passes through

the point carrying 50% of the failure deviatoric stress.

The least squared methods is applied again

on the unloading reloading branch

of the stress strain curve

to calculate E subscript UR.

In the third step of the procedure,

the stress dependency exponent M is estimated.

In stress strain curves of drain triaxial tests

with different confining pressures,

half failure stresses are first calculated

using the method described in step two.

A new natural log space can then be defined

that represents the half failure stiffnesses of each test

in the form of a single point.

The point corresponding to the reference test

is located on the origin.

The stress dependency exponent M

is the slope of a trend line

that passes through these points.

The similar steps can be taken for unloading reloading

in order to calculate

the stress dependency exponent twice.

In the hardening soil model however,

it is assumed that the constant M

is the same for both half failure

and unloading reloading stiffnesses.

The fourth step of the parameterization procedure

is to calculate the failure ratio, RF.

This constant has a value between 0.5 and one,

and is usually assumed to be 0.9.

A more accurate value for this constant can be calculated

based on the results of the triaxial tests.

First, the stress strain hyperbolic function

is rewritten in terms of the failure ratio,

R subscript F.

This expression introduces a new coordinate system

in which each hyperbolic curve

is projected to a straight line.

The failure ratio is the slope of these parallel lines.

In the second last step of the parametrization procedure,

the odometer reference stiffness

and the coefficient of lateral earth pressure are estimated.

The odometer test result is required for the step.

However, these values can be approximated

in the case where the appropriate laboratory result

is not available.

As mentioned earlier, the odometer reference stiffness

is defined in a one-dimensional compression test.

By definition, it is the tangent stiffness

at the point where the vertical stress

is equal to the reference stress.

In addition, the coefficient of lateral earth pressure

for a normally consolidated soil, K nought NC,

is the ratio of lateral stress, Sigma three,

to vertical stress, Sigma one, in the odometer test.

These two constants must be approximated

in the absence of laboratory data.

For instance, the value of the odometer reference stiffness

is usually greater than 50%

of E 50 rough and less than 80% of E 50 rough.

And also the coefficient of lateral earth pressure

can be estimated using the well-known Jackie’s formula,

which relies on the friction angle of the soil.

One additional step in the privatization procedure

is required if the soil is over consolidated,

and that is the determination of the isotropic over

consolidation ratio,

which is a ratio equivalent isotropic stresses.

The denominator of this expression

is simply the confining pressure of the triaxial test.

In contrast, the numerator, P subscript C,

defines the size of the yield locus

passing through the yields stresses QY and PY

As such, P subscript C

must be calculated from the equation

of the yield function.

The yields stresses are generally determined

from plots of deviatoric stress versus axial strain

and void ratio versus natural log mean effect of stress.

And in both of these plots, the yield point,

that is that transition from

elastic to elastic plastic behavior

it deduce by an inflection point

in the measured laboratory data.

That brings us to the constant mu,

which is in the equation for the cap yield function.

It is calculated by the material model in the solver

from this expression here.

We elaborate on this expression in the reference book,

but for the purpose of this webinar,

it suffices to say that theta is

Q/P stress ratio

that can be calculated from the user input,

K nought NC,

which we discussed on the previous slide

and lambda/Kappa is calculated from

failure’s ratio


this ratio of reference stiffnesses.

Now the entire procedure for determining OCR

can be simplified if the stress path is instead

that of isotropic compression, because in that case,

the stress path tracks along

the mean effect of stress axes.

Meaning that the deviatoric stress is zero

and PC is equal to PY.

To demonstrate the parametrization procedure

we are going to consider 12 drain triaxial compression tests

on Ottawa sand.

These tests were performed at different confining pressures

ranging from 50 to 600 kilopascals.

The effective friction angle and the cohesion of the soil

are calculated from the first step

of the parametrization procedure.

The deviatoric compression stress in 12 triaxial tests

continued until failure was achieved.

The failure state of each task can be represented

by a point in PQ space.

The failure line that passes through these failure points

has a slope of 1.18

and a Y intercept of zero.

Therefore the effect of friction angle

as calculated as 29.6 degrees

while the cohesion is calculated as zero,

since alpha is zero.

At the second step of the procedure,

the test with the confining stress of 100 kPA

is considered as the reference test.

The half failure deviatoric stress for this reference test

is calculated to be 97.6 kilopascals.

According to the measured stress strain curve,

the axial strain corresponded into this half failure

stress is equal to 0.0055.

This means that the half failure referenced stiffness

for the soil is 17,745 kilopascals.

Again, that is simply the slope of this line.

Finally, the unloading reloading branch

of the measured stress drain curve

gives us the value of the unloading reloading

reference stiffness.

And the slope of that line

is estimated to be 45,000 kilopascals.

The third step of the parametrization procedure

estimates the stress dependency exponent M.

As shown here this constant can be estimated

from either the unloading reloading

or the half failure stiffnesses in natural log,

natural log space.

The slope of these two lines is almost the same

and their average is 0.680.

And this is what we consider to be the constant M.

As can be seen here,

the proposed expressions can adequately estimate

the unloading reloading and half failure stiffnesses

in terms of their corresponding confining stress

Sigma three.

The proposed curve, however,

is more consistent with the measured results

for the unloading reloading stiffnesses

than for the half failure stiffnesses.

This is commonly the case.

The failure ratio, R subscript F

is calculated at the fourth step.

The linear projections of the hyperbolic curves

are shown in the suggested space in the left graph.

The near vertical branches of these curves

represent unloading reloading phases,

where Q approaches zero

causing the ratio of QF/Q

to approach infinity.

These near vertical branches are not included

in the calculation of RF,

rather it is the sloping portions

of the curves that are described

by the equation of a line

where the slope of each line is our subscript F.

The graph on the right shows the slopes

for all 12 triaxial tests.

The average of 0.941 is then estimated

for the failure ratio of the soil.

For the odometer related constant in step five,

we apply the approximate method

of calculating the odometer reference stiffness.

According to the value

of the half failure reference stiffness,

which was calculated at the second step,

we estimate the odometer reference stiffness to be 115,

excuse me, 11,500 kilopascals.

The Jackie’s formula is also used to calculate the

coefficient of lateral earth pressure

in terms of the effect of friction angle

and from that relationship,

we obtained 0.506.

The results the parametrization procedure

for the Ottawa sand can be summarized as follows.

Our two strength constants gave us a friction angle

of 29.6 and a cohesion of zero.

The reference dress was 100 kilopascals,

and the reference stiffnesses were calculated accordingly.

The stress dependency exponent and the failure ratio

were estimated as 0.68

and 0.941 respectively,

and the coefficient of lateral earth pressure

was obtained from Jackie’s formula

to be 0.506.

Now that we have reviewed the step-by-step parametrization

procedure for the Hardening soil model,

I will describe how to model a Triaxial test in Sigma,

and compare the numerical results

with the previously mentioned lab data.

The numerical model configuration

is based on the geometry and boundary conditions

of the triaxial test conducted in the laboratory.

A typical triaxial sample has a cylindrical shape

with a diameter of five centimeters

and a height of 10 centimeters.

Due to axis symmetric shape of this geometry

with respect to the central vertical axis,

we use a 2D axis symmetric geometry and Sigma.

The geometry is also symmetric

with respect to the horizontal axis.

So we can model the top one quarter of the specimen.

An eight nodes single element model

is considered for this resulting rectangle.

The material definition is of course the key subject

of this webinar.

This step is performed in the software

using the defined materials dialog box.

The material model is set to the Hardening soil model,

the sample parameters

including the initial void ratio

and over consolidation ratio are entered.

It should be noted that the effect of the soil self weight

can be removed by setting the unit weight to zero.

And lastly, the constitutive model constants

are entered based on what was previously obtained

during the parametrization procedure

for the Ottawa sand.

The triaxial tests are simulated

by firstly considering the isotropic confining pressure.

And secondly, the deviatoric pressure.

To do this, the project file is set up

with an analysis tree structure for each test.

In each branch,

we see the confining pressure phase

acting as a parent,

which provides initial data

to the subsequent displacement control deviatoric phase.

In test A one as an example,

a constant normal stress of 50 kPA

is applied on both the right and top sides of the model

using a normal tangential stress boundary condition.

And in the deviatoric phase of this test,

a displacement of two centimeters

that is equivalent to 40% strain,

is applied on the top side of the model.

The displacement spline data function

contains a small plip

that just before

a time and a lapse time

of 1.2,

and this is the unloading reloading phase.

Let’s just take a few minutes to explore

some of the finer details of these ostensibly simple

numerical analysis.

First, we’ll start at the top of the project definition.

And I will point out that you can right click on here

and add multiple geometries to the same project file.

In this analysis, of course we have only one single geometry

and it is of the type 2D Axis symmetric.

The next noteworthy thing under our geometry sits

our confining stress analysis

for all 12 triaxial tests

and then unlike what is shown

in the PowerPoint presentation,

we actually have two cases

for each deviatoric loading condition.

One that is monotonic,

meaning strained without unload reload,

and the other with the LUR or Load Unload Reload branch.

All three of the analysis in each branch

are low deformation.

Now the first analysis

is the one in which we apply the confining pressure

of 50 all the way through to 600 kilopascals.

And it is an analysis of the type load deformation.

We’re not doing an in-situ or gravity activation analysis



dimensions of the specimen

are so negligible

that the body load or the gravity load

doesn’t contribute to the stresses in any meaningful way.

The other thing that’s worth noting is that

our deviatoric phases

are also low deformation,

and the initial stresses and poor pressures

come from the parent analysis.

We also have the ability to set

the final pour water pressure conditions.

And you can see in this case,

it’s the same as the initial conditions,

meaning that there’s no increment or decrement

in the pour water pressure.

Had these been different due to, for example,

specifying the pour pressures

from another geostudio analysis like Sigma CAP

or the use of a water table or

special function to find the pressure heads,

then not change in pour water pressure

would be reflected

in the deformation response,

the simulated deformation response of the specimen.

I should also note while I’m here,

that there are these two options to reset

displacements and strains

and reset state variables.

These are not on automatically

when a low deformation analysis is added

to another low deformation analysis.

And so they were not set in this case.

Now that means that any of the deformations accumulated

during the isotropic compression

will carry through

into the deviatoric loading phase.

And so in actuality,

it would have been better modeling practice

to toggle the reset displacements and strains,

and these are accumulated values to toggle that on.

If you’re using a material model that has state

like the hardening soil model in the confining phase,

then you can also reset the state variables.

Now, in this case,

what we’ve done is we’ve used

an isotropic elastic material model

in the compression phase,

and we’ve given it a very high stiffness quantity

of 15 million kilopascals,

meaning that the simulated deformations

from this phase are negligible.

And therefore we need not reset

the accumulated displacements and strains.

The other point I’ll make is that

when you deal with large project files

containing multiple analysis,

there are right click options in the solve manager

to uncheck all,

to check children, which picks a single branch,

that’s the zero one branch.

And then there are a couple other options.

Similarly, if you go up into the solve manager,

it’s often more convenient

if these branches are quite elaborate,

in this case there are only two analysis

that are parenting to that,

the confining stress analysis.

But we can imagine a scenario

where they’re eight, 12, and so on analysis underneath

the parent analysis.

In which case right clicking up here

and queuing an analysis branch for solving

is more convenient.

You can also

end queuing,

you can clear, and you can also

queue one particular analysis.

Now, in this case,

I’m going to right click and


all analysis

and then solve.

And once the solution is obtained for any one of these,

we can switch to the results view.

And it is here where we can

go to draw graph



the response through various graphs,

such as deviatoric stress versus axial strain,

deviatoric stress versus mean effective stress

and void ratio versus P.

You can see that there’s a unload branch in there.

There’s the monotonic loading case.

You can also multi-select these graphs,

which of course gives a rather strange view

since the axes are not all the same,

but I just wanted to point out

that you can copy this graph data

and then go to a spreadsheet program like Excel,

where you can then


and paste the data into the spreadsheet.

This will then make for easier comparison

with the measured laboratory data

If you turn out to the PowerPoint presentation,

the simulation results of six triaxial trust tests

with confining pressures ranging from 50 to 300 kilopascals

are compared here with laboratory results.

Colored curves are laboratory results.

And the black lines are the results

from the Sigma simulation

and analytical solution is also shown on these plots

as gray dash lines.

The full compatibility of the analytical

and simulated solutions verifies the implementation

of the hardening soil model in SIGMA/W.

In addition,

the acceptable consistency

between the lab data and the simulator response

confirms that the hardening soil model

with carefully parameterized material constants

can capture the stress strain response of the Ottawa sand.

The same trends can be observed for the next six tests

that had confining pressures ranging from

350 kilopascals to 600 kilopascals.

To conclude, the parametrization procedure

of the hardening soil model

was provided in five straightforward steps

and a sixth optional step for over consolidated soil.

It was shown that the results from

multiple drain triaxial tests

and at least one odometer tests

are required to parameterize

the hardening soil model.

In the absence of odometer data,

approximations can be used to obtain

the reference odometer stiffness

and the K nought NC

or normally compressed earth pressure coefficient inputs.

The material constants were then used

in numerical simulations of the lab tests

and the numerical results were shown

to compare favorably with both the laboratory data

and analytical solutions.

During the presentation,

I referenced these top two


Shanz, Vermeer and Bannier,

that’s the original model formulation.

I also want to remind you that you can consult

our stress strain modeling

with GeoStudio 2021 reference book,

to understand the details of the formulation.

And then this is the paper that comprises

the laboratory tests on the Ottawa sand.

I’ll also remind you that

the website comprises a number of resources,

including example files.

One noteworthy example is that of

the tieback wall in Berlin,

it is completed using both the Mohr-Coulomb

and hardening soil models

and those results are compared.

With that, we’ve reached the end of the webinar.

Our recording will be available to view online.

Please take time to complete the short survey

that appears on your screen

so that we know the types of webinars

that are of interest to you.

Thank you very much for joining us.

The video transcript gets copy and pasted here