_{ cap }

Capillary force due to water between two grains of soil

DDiameter of soil's particles

u_{ a }

Pressure in the air phase

s_{ e }

Suction of air entry

u_{ w }

Pressure in the water phase

σ_{c}

Capillary stress

A_{ REV }

Cross-sectional area of the REV in the plane normal to the direction observed

γWater-air surface tension

eVoid ratio

RA non-physical parameter to fix yield function

sSuction

pTotal mean stress

p′Effective mean stress

qDeviatoric stress

TTemperature

M(T)Slope of critical state line influenced by temperature

p′_{c}(σ

_{c}, T)

Preconsolidation stress influenced by capillary stress and temperature

a_{1}, a

_{2}

Parameters defining thermal softening curve

αThermal expansion parameter

λ(σ_{c})

Stiffness parameter for changes in effective mean stress for virgin states of the soil influenced by capillary stress

κElastic stiffness parameter for changes in effective mean stress

kParameter controlling the rate of increase of preconsolidation stress with suction

e_{ 0 }

Initial void ratio

GShear modulus

dɛ_{v}

^{e}

Increment of elastic volumetric strain induced by changes in effective mean stress or temperature

dɛ_{v}

^{p}

Increment of plastic volumetric strain induced by changes in effective mean stress or temperature

dɛ_{s}

^{e}

Increment of elastic deviatoric strain induced by changes in deviatoric stress

dɛ_{s}

^{p}

Increment of plastic deviatoric strain induced by changes in deviatoric stress

## 1 Introduction

Thermohydromechanical (THM) behaviour of unsaturated soils is a subject of interest for civil engineers, for example, in the use of geothermal energy or for nuclear waste disposal. Over the past three decades, the study of THM behaviour and non-isothermal mass flow in porous media has been a key area of research in environmental geomechanics.

Many studies have been carried out in the last forty years concerning the impact of temperature changes on saturated soil behaviour during isotropic and oedometric loadings. Studies had for object the influence of temperature on preconsolidation stress (Eriksson, 1989; Tidfors and Sàllfors, 1989), on stiffness parameters (c_{s} or κ and c_{c} or λ) (Burghignoli et al., 2000; Campanella and Mitchell, 1968; Cekerevac and Laloui, 2004; Despax, 1976; Eriksson, 1989; Fleureau, 1972; Graham et al., 2001; Tidfors and Sàllfors, 1989) and the evolution of thermal conductivity on a thermo-hydro-mechanical path (Djéran-Maigre, 1993; Djéran-Maigre et al., 1997).

Concerning the failure criterion, defined by the slope M of the critical state line in the (p′, q) stress plane, the experimental results obtained by several researchers show some divergence: Belanteur et al. (1997), Burghignoli et al. (2000), De Bryn (1999), Despax (1976), Hicher (1974), Hueckel and Pellegrini (1989), Hueckel and Baldi (1990) noted a reduction of the slope M with temperature. Other results were obtained by Cekerevac and Laloui (2004), Graham et al. (2001, 2004), Hueckel and Pellegrini (1989) and highlighted the absence of influence of temperature on M. All these results led us to be particularly attentive concerning the kind of material tested, the initial state and the way of preparation of the specimens. Hueckel and Baldi (1990) and Cekerevac and Laloui (2004) showed that the reduction of the strength of soils with temperature is amplified for low values of the OCR.

In the case of unsaturated soils, there are few results in the literature (Recordon, 1993; Saix, 1991; Saix et al., 2000; Romero, 1999; Weibe et al., 1998). On triaxial paths, Romero (1999) on Boom clay and Weibe et al. (1998) on samples of compacted sand-expansive clay mixtures, showed a contraction of the yield surface with temperature at any suction. On the other hand, the works of Jamin et al. (2004) showed a thermo-dilatancy of the yield surface for samples of sand-loam mixture. These observations highlight the major role of the clay fraction of the soil on the yield surface change (El-Youssoufi, 2006).

### 1.1 Modelling

Two points of view have been proposed in literature to model the THM behaviour of unsaturated soils, one based on independent stress analysis (Alonso et al., 1990; Mathyas and Radhakrishna, 1968), and the other, on the effective stress concept.

In the latter case, it is admitted that, as for saturated soils, the behaviour of unsaturated soils depends only on one parameter, which is the combination of stress and suction. This concept, in association with a hardening model of soil with suction, presents many advantages, in particular because it assures a continuous passage from saturated to unsaturated medium. The validity of the model highly depends on the correct definition of the capillary stress which is the contribution of suction to the effective stress (Aboubekr, 1995; Biarez et al., 1993; Modaressi and Aboubekr, 1994; Loret and Khalili, 2002; Taibi, 1994). Under non-isothermal conditions, the expression of the constitutive laws for unsaturated soils is formulated in the same way as under isothermal conditions. Lewis et al. (1986, 1998) have developed a finite elements model for studying the hydromechanical behaviour and heat transfer in saturated and unsaturated porous media.

Recently, two new models, one based on a generalised effective stress concept (Loret and Khalili, 2002), and the other one based on the independent variable approach (Wenhua et al., 2004), have been developped. In their model, Loret and Khalili (2002) propose a thermohydromechanical model for unsaturated soils, based on the effective stress definition of Bishop, with the χ parameter being a function of suction. On the other hand, for Wenhua et al. (2004), the THM model for unsaturated soils is an extension of the thermomechanical models for saturated soils proposed by Hueckel and Borsetto (1990) and Cui et al. (2000).

## 2 Thermohydromechanical model for unsaturated soil

### 2.1 Effective stress in unsaturated soils

A simple way to define effective stresses is to start from the expression of the intergranular forces between two particles in an idealized medium at the microscopic scale (Aboubekr, 1995; Biarez et al., 1993; Taibi, 1994). The method consists in calculating the capillary force F_{cap} due to water between two grains of soil modelled by rigid spheres with the same diameter D (Fig. 1).

In a representative elementary volume (REV), the capillary stress σ_{c} in any direction is defined by:

$${\sigma}_{c}=\frac{\sum {F}_{cap}}{{A}_{REV}}$$ | (1) |

_{cap}is the sum of the forces acting in the direction; A

_{REV}is the cross-sectional area of the REV in the plane normal to this direction; σ

_{c}is a function of s the suction, of the void ratio e of the material, of the diameter D of the spheres (Biarez et al., 1993; Taibi 1994) and expressed as:

$${\sigma}_{c}=\frac{2\pi \gamma}{K\left(e\right){D}^{2}}\left[2D+\frac{3\gamma -\sqrt{4\gamma Ds+9{\gamma}^{2}}}{s}\right]$$ | (2) |

$$K\left(e\right)=0.32{e}^{2}+4.06e+0.11$$ | (3) |

Fig. 2 illustrates the transient between the two cases in the (s, σ_{c}) plane. For the highest suctions, capillary stress reaches a horizontal asymptote, where σ_{c} is equal to σ_{cmax} defined as:

$${\sigma}_{c\text{\hspace{0.17em}}\mathrm{max}}=\frac{4\pi \gamma}{K\left(e\right)D}$$ | (4) |

The capillary stress appears as an isotropic pressure, that replaces the matrix suction when the water phase is discontinuous. The expression of the effective stress tensor becomes:

$$\left[{\sigma}^{\prime}\right]=\left[{\sigma}_{\text{net}}\right]+{\sigma}_{c}\left[I\right]$$ | (5) |

When considering a real soil with different diameters of particles, the determination of the capillary stress is difficult, because of the need to choose a “characteristic diameter” D equivalent to the real grain-size distribution. In practice, this diameter is a parameter of the model and must be derived from experimental data (Biarez et al., 1993, 1994).

### 2.2 The elastoplastic general framework

The purpose of the model presented in this article is to associate the capillary stress defined for unsaturated soils in the previous part with an elastic-plastic behaviour law like Modified Cam-Clay. The yield function used in this model is based on the Khalili and Loret (2001) yield function in which suction is replaced by the capillary stress σ_{c}. Moreover, the increase of preconsolidation stress with suction is proposed and defined as a function p′_{c}(σ_{c}) = f((p_{c0}, k, σ_{c}) where k is a parameter depending of soils. The thermal effect is coupled in the model with the thermal softening curve as explained by Wenhua et al. (2004).

Concerning the influence of temperature on M, the slope of the critical state line in the plane (p′,q), experimental results (De Bryn, 1999; Despax, 1976) show a decrease of M and q_{max} with temperature. These results lead us to introduce a function defining the evolution of M with temperature as follows:

$$M\left(T\right)=M\left(T-{T}_{ref}\right)+aT$$ | (6) |

_{ref}is a reference temperature and a, a negative empirical parameter.

The yield function is described by two distinct expressions depending on the $\frac{{{p}^{\prime}}_{c}}{R}$ ratio, as follows:

$$\frac{{q}^{2}}{M{\left(T\right)}^{2}{p}^{\prime}}+{p}^{\prime}-\frac{2}{R}{{p}^{\prime}}_{c}\left({\sigma}_{c},T\right)\text{\hspace{0.17em}}\text{for}\text{\hspace{0.17em}}{p}^{\prime}<\frac{{{p}^{\prime}}_{c}\left({\sigma}_{c},T\right)}{R}$$ | (7) |

$${\left(R-1\right)}^{2}\frac{{q}^{2}}{M{\left(T\right)}^{2}{p}^{\prime}}+{p}^{\prime}-\frac{2}{R}{{p}^{\prime}}_{c}\left({\sigma}_{c},T\right)+\left(\frac{2}{R}-1\right)\frac{{{p}^{\prime}}_{c}^{2}\left({\sigma}_{c},T\right)}{{p}^{\prime}}\text{\hspace{0.17em}}\text{for}\text{\hspace{0.17em}}{p}^{\prime}>\frac{{{p}^{\prime}}_{c}\left({\sigma}_{c},T\right)}{R}$$ | (8) |

_{c}(σ

_{c},T) is the preconsolidation pressure influenced by capillary stress and temperature, R is a non-physical parameter used to fix the yield function (R = 2 for Cam-Clay, R = 1.6 for Loret and Khalili) and to derive the consolidation line and critical state line from one another using the experimental data. The continuity of yield function is preserved at $\frac{{{p}^{\prime}}_{c}\left({\sigma}_{c},T\right)}{R}$ value. The yield function is kept elliptical as in the Cam-Clay model and normal to the p′ axis for p′ = p′

_{c}(σ

_{c},T).

### 2.3 The thermal and the capillary stress softening curves

Many authors have carried out thermomechanical tests consisting in oedometric and triaxial loadings, under constant temperature, on overconsolidated samples in a saturated state. They have observed a nonlinear decrease in the preconsolidation stress when the temperature increases (Eriksson, 1989; Ghembaza, 2004; Tidfors and Sàllfors, 1989). This form was adopted to define a yield surface in the (p′, T) plane by Hueckel and Borsetto (1990) and also by Wenhua et al. (2004) who used it in their model as thermal softening curve (TS). This expression has the following form:

$${{p}^{\prime}}_{c}\left(T\right)={{p}^{\prime}}_{c}+A\left(\Delta T\right)$$ | (9) |

$$A\left(\Delta T\right)={a}_{1}\Delta T+{a}_{2}\Delta T\left|\Delta T\right|$$ | (10) |

_{1}and a

_{2}are coefficients depending on the thermal sensitivity of the soil, ΔT is the increase in temperature from a reference temperature, p′

_{c}is the preconsolidation pressure for saturated state and a temperature of reference. The term A(ΔT) is independent of the initial value of the preconsolidation pressure. Suction and capillary stress have a great influence on unsaturated soils behaviour and preconsolidation pressure p′

_{c}. To take into account this increase in p′

_{c}with suction, a parameter k is introduced in association with σ

_{c}. This equation gives the increase in p′

_{c}as a function of capillary stress, starting from an initial preconsolidation pressure under saturated conditions p′

_{c0}. This equation is written as follows:

$${{p}^{\prime}}_{c}\left({\sigma}_{c}\right)={{p}^{\prime}}_{c0}+k\left({\sigma}_{c}-{s}_{e}\right)$$ | (11) |

_{e}is the suction of air entry and k is an empirical parameter derived from appropriate experiments to quantify the increasing of preconsolidation stress with capillary stress. This phenomenon takes place only for suctions higher than the suction of air entry. The plastic slope is influenced by capillary stress and takes the form λ(σ

_{c}) = f(λ

_{0}, σ

_{c}), where λ

_{0}is the plastic slope for saturated soils.

To couple the effect of temperature with that of suction, the relationship is derived from that of Wenhua et al. (2004) in which suction has been replaced by capillary stress. The relationship between the preconsolidation pressure p′_{c}(σ_{c}, T), the capillary stress and temperature T in unsaturated soils becomes:

$${{p}^{\prime}}_{c}\left({\sigma}_{c},T\right)={{p}^{\prime}}_{c0}\left(T\right)+k\left({\sigma}_{c}-{s}_{e}\right)$$ | (12) |

For this coupling between temperature and capillary stress for unsaturated soils, p′_{c0}(T) must first be calculated from p′_{c0} using the TS relation then p′_{c}(σ_{c}, T).

The response of the volumetric elastic strain to a change in p′ and temperature T becomes:

$$\left(\text{d}{\epsilon}_{v}^{e}\right)=\frac{k}{1+{e}_{0}}\text{\hspace{0.17em}}\frac{d{p}^{\prime}}{{p}^{\prime}}+\alpha \Delta T$$ | (13) |

$$\text{d}{\epsilon}_{s}^{e}=\frac{1}{3G}\text{d}q$$ | (14) |

### 2.4 Plastic behaviour

The set of constitutive equations for the plastic part consists in the following incremental relations. For an isotropic test in which a soil sample, at a given suction, is submitted to mechanical load increments along the normally consolidated (NC) line, the incremental volumetric plastic strain is given by:

$$\text{d}{\epsilon}_{v}^{p}=\frac{\lambda \left({\sigma}_{c}\right)-\kappa}{1+{e}_{0}}\text{\hspace{0.17em}}\frac{\text{d}{p}^{\prime}}{{p}^{\prime}}$$ | (15) |

For suction values higher than s_{e}, p′_{c}(σ_{c}) increases faster than p′ and this leads to hardening of the soil. Consequently, only elastic strains occur.

Concerning NC or slightly overconsolidated (OC) materials at ambient temperature, an increase in temperature induces volumetric plastic strains due to activation of Thermal Softening mechanism (TS). If p′_{c}(σ_{c}, T) > p′, the preconsolidation stress p′_{c}(σ_{c}, T) decreases when temperature increases and p′_{c}(σ_{c}, T) is equal to the effective mean stress p′ in some cases. For these particular cases, the preconsolidation stress decrease is stopped because p′_{c}(σ_{c}, T) cannot be lower than p′. Plastic strains appear and are calculated as indicated below:

$$\text{d}{\epsilon}_{v}^{p}=-\frac{\left(\lambda \left({\sigma}_{c}\right)-\kappa \right)\left({a}_{1}+2{a}_{2}\Delta T\right)}{\left(1+{e}_{0}\right){{p}^{\prime}}_{c}\left({\sigma}_{c}\right)}\Delta T$$ | (16) |

_{1}and a

_{2}are the coefficients depending on thermal sensitivity of the soil used in relation (10). Considering an associated flow-rule with respect to effective stresses and the yield function written in relations (7) and (8), the incremental plastic deviatoric strain is obtained by:

$$\frac{\text{d}{\epsilon}_{s}^{p}}{\text{d}{\epsilon}_{v}^{p}}=\frac{\frac{2q}{M{\left(T\right)}^{2}{p}^{\prime}}}{1-\frac{{q}^{2}}{M{\left(T\right)}^{2}{{p}^{\prime}}^{2}}}\text{\hspace{0.17em}}\text{for}\text{\hspace{0.17em}}{p}^{\prime}<\frac{{{\text{p}}^{\prime}}_{\text{c}}\left({\sigma}_{c},T\right)}{R}\text{\hspace{0.17em}}$$ | (17) |

$$\frac{\text{d}{\epsilon}_{s}^{p}}{\text{d}{\epsilon}_{v}^{p}}=\frac{\frac{2q}{M{\left(T\right)}^{2}{p}^{\prime}}{\left(R-1\right)}^{2}}{1-\frac{{\left(R-1\right)}^{2}{q}^{2}}{M{\left(T\right)}^{2}{{p}^{\prime}}^{2}}-\left(\frac{2}{R}-1\right)\frac{{{p}^{\prime}}_{c}\left({\sigma}_{c},T\right)}{{{p}^{\prime}}^{2}}}\text{\hspace{0.17em}}\text{for}\text{\hspace{0.17em}}{p}^{\prime}>\frac{{{p}^{\prime}}_{c}\left({\sigma}_{c},T\right)}{R}$$ | (18) |

## 3 Simulations of typical phenomena of unsaturated soils

Table 1 gathers the parameters used for simulations below.

### 3.1 Collapse path

In unsaturated soil, wetting under constant net mean stress may result in collapse (volumetric contraction). When an unsaturated soil is loaded then wetted, the preconsolidation stress p′_{c}(σ_{c}, T) decreases faster than the effective mean stress p′. As long as p′_{c}(σ_{c}, T) is higher than p′, elastic swelling is observed and may be calculated using relation (13). When p′_{c}(σ_{c}, T) is equal to p′, the yield surface evolves to produce only irreversible contraction calculated with relation (16).

Fig. 3 illustrates collapse. The corresponding stress path is shown in Fig. 3a. The soil is first dried to a suction of 5 MPa (larger than the air entry suction s_{e} = 1 MPa). Afterwards, mechanical loading is applied up to p = 1000 kPa with p′ = p′_{c}(σ_{c}), followed by a wetting path (Fig. 3b). During this wetting step, capillary stress decreases like p′_{c}(σ_{c}) and irreversible contraction appears (Fig. 3c) highlighting the collapse phenomenon. This phenomenon takes place as long as s > s_{e}. When s < s_{e}, elastic swelling is observed (Fig. 3c).

### 3.2 Coupled drying-heating loading

Fig. 4 shows a coupled drying-heating path with a constant total stress (p = 200kPa) on normally consolidated samples in the (T, e), (ln s, e), (ln p′, e); (s, p′_{c}) and (T, p′_{c}) planes. During the drying step, the same observations as those described in section 3.1 for the drying path can be done. After the drying step, a heating path is simulated by a temperature increase. At the beginning of the heating, the behaviour is elastic as evidenced by the yield surface position (Fig. 4b), and the preconsolidation stress p′_{c}(σ_{c}) which is larger than the effective mean stress p′. The behaviour is first elastic with thermal expansion and the preconsolidation stress decreases according to the relation 20.

As long as temperature continues to increase, p′_{c}(σ_{c}, T) is equal to p′, temperature continues to increase and results in irreversible deformations and irreversible contraction as shown in Fig. 4b. In this case, p′_{c}(σ_{c}, T) keeps a constant value, as shown in Fig. 4e.

This simulation on a normally consolidated sample shows the predominant effect of the drying path in the saturated zone, which induces the larger part of total and plastic strains. Concerning reversible strains, the part linked to the drying path in the unsaturated zone is larger than the part linked to heating but stays of the same order of magnitude (Fig. 4a and b).

### 3.3 Deviatoric path under constant suction and imposed temperature

Fig. 5 shows the effect of temperature on the coefficient M, the slope of the Critical State Line in the (T, M) plane (Despax, 1976). When the sample is heated from 22 °C to 80 °C, M decreases from 0.8 to 0.4.

Fig. 6 shows the results of two triaxial tests at reference temperature (T = 22 °C) and high temperature (T = 80 °C) under constant net mean stress (p = 800 kPa) and suction (s = 10 MPa), higher than air entry suction value (s_{e} = 1000 kPa) in the (ɛ_{1}, q), (p′, q), (ɛ_{1}, e), (ln p′, e) planes. Both samples have been first isotropically consolidated in saturated state at ambient temperature up to p = 800 kPa. Then, suction (s = 10 MPa) has been applied and one of it has been heated to 80 °C before starting triaxial tests.

For the sample at ambient temperature, the path is first elastic then elastoplastic. This phenomenon is due to suction and capillary stress that increase the preconsolidation stress p′_{c}(σ_{c}) and consequently the yield surface. Elastic strains are calculated by relation (13). When the yield surface is activated, elastic and plastic strains occur that are estimated by the increase in the preconsolidation stress and flow rules described in section 2. In the (q, ɛ_{1}) plane, the elastic part of the strains for the test at 80 °C is smaller than those calculated for the test at reference temperature. For the sample heated before the triaxial test, the preconsolidation stress is reduced with the yield surface contraction. The behaviour is only elastoplastic (Fig. 6a), because with thermal expansion, the preconsolidation stress decreases according to the relation (12) and p′_{c}(σ_{c}, T) is equal to p′ when triaxial test begins. Both samples reach a plateau for different values because the slope of the Critical State Line decreases and the failure criterion is reached for a lower effective mean stress p′.

For sample at ambient temperature, the void ratio first follows a straight line (slope κ) when the behaviour is elastic, and it reaches the critical state line when the yield surface is activated and the stresses are close to the failure criterion.

## 4 Comparison of model predictions with experimental results

### 4.1 Drying and wetting cycle

Fig. 7 shows a drying-wetting cycle on a mixture of kaolin and sand prepared at a water content of 1.5 w_{L} (where w_{L} is the liquid limit): specimens were first dried up to 10 MPa of suction, then they were wetted (Ghembaza, 2004). The air entry suction s_{e} is almost equal to 0.3 MPa. As shown in Fig. 7, the model has the capacity to simulate this kind of path. For s > s_{e}, the increase in suction has a limited influence on void ratio and volumetric strain. Moreover, the irreversible strains on the drying and wetting paths are well simulated.

#### 4.1.1 Deviatoric path under constant suction

Fig. 8 shows the results of triaxial tests on the same kaolin-sand mixture under constant suction (Ghembaza, 2004). The samples are first consolidated to an effective mean stress p′ = 1250 kPa under saturated condition, then suction is applied. After stabilization under the required suctions (3 MPa and 8.5 MPa), deviatoric loading begins.

## 5 Concluding remarks

The modified Cam-Clay model has been chosen as a plastic driver to keep the model for unsaturated soils and thermal effects as simple as possible. The replacement of suction by capillary stress in Terzaghi's effective stress definition is used to extend relations defined for saturated soils to unsaturated soils. Moreover, this capillary stress is based on a microstructural model and calculated from attraction forces due to water menisci.

The introduction of the thermal softening curve to describe the thermal effect on soils in the elastoplastic model allows to reproduce some coupled paths like a drying path under mechanical stress with transition between saturated and unsaturated states, a heating path under constant suction or a deviatoric path with imposed suction and temperature.

The typical compression phenomena observed in saturated normally or slightly overconsolidated soils, or the swelling of overconsolidated soils (OCR > R) submitted to coupled thermomechanical loading have been highlighted in this study.

In the (p′, q) plane, the failure criterion is intrinsic and its slope M remains constant for a given saturated or unsaturated soil at ambient temperature. This criterion is affected by temperature, and the introduction of the decreasing function M = f(T) in the model allows to reproduce the decrease in shear strength during a triaxial test at high temperature.