# Microporomechanics study of anisotropy of ASR under loading

Laurent Charpin, Alain Ehrlacher

► To cite this version:

Laurent Charpin, Alain Ehrlacher. Microporomechanics study of anisotropy of ASR under loading. Cement and Concrete Research, 2014, 63, pp.143-157. 10.1016/j.cemconres.2014.05.009. hal-01157339

HAL Id: hal-01157339

https://enpc.hal.science/hal-01157339

Submitted on 16 Jun 2021

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
Microporomechanics study of anisotropy of ASR under loading

Laurent Charpin<sup>a,*</sup>, Alain Ehrlacher<sup>a</sup>

<sup>a</sup>Laboratoire Navier, ENPC, 6-8 av Blaise Pascal, Cité Descartes, Champs-sur-Marne, 77455 Marne-la-Vallée Cedex 2, France

**Abstract**

In this article, we introduce a new micromechanical model for alkali-silica reaction. Our idea was to build a model with the following characteristics. First, the model has to be simple enough to be used to compute damage under loading and chemical attack at the level of each element in a structure code. Second, its parameters must be easy to identify on available alkali-silica reaction lab experiments. We have chosen to model the behavior of a concrete containing aggregates such that most of the damage occurs in the cement paste. Using micromechanics and an energy criterion, the model remains analytical except for the minimization of the energy. The parameters were identified on Multon's triaxial experiments and good results were obtained for compressive loadings up to 10 MPa.

**Keywords:** Alkali-aggregate reaction, Anisotropy, Micromechanics, Energy criterion

**1. Introduction**

Concrete is vastly used for construction and the life span of buildings is an important parameter of their rentability. Therefore durability issues such as freeze-thaw, carbonation or endogenous reactions need to be addressed. The alkali-silica reaction is an endogenous reaction discovered in the US in the 40's by Stanton [66]. It can develop several years after building. The reaction is externally visible because of a characteristic surface cracking which tends to align with the main compression directions and also through reaction products which sometimes flow from these cracks. The reaction also induces irreversible expansions which can lead to complications when using the structure. Various kinds of buildings can be affected. They all share some characteristics. First, their aggregates must contain some reactive silica, in the sense that it can dissolve in the concrete interstitial solution. Second, the concrete must contain alkali, usually provided by the cement. Third, there must be a water supply to the structure, since water enhances the transport of chemicals inside the concrete and is absorbed by the alkali-silica gels, leading to their swelling.

How to stop this reaction is an open problem. Today, a lot of effort has been made to make sure that no new building will be affected. For affected structures, mechanical techniques exist to release the stresses due to the reaction or try to stop the expansion mechanically. However, these solutions are expensive. It is therefore necessary to have simulation tools available to estimate the remaining life-span of affected structures and the efficiency of repair

solutions. In this paper we give some elements for the numerical simulation of ASR.

ASR modeling can roughly be divided in three branches. First, chemistry and transport models studying at the aggregate or structure scale, the transport of chemicals and water, so as to predict the amounts of gel produced at reactive sites, that is close to or in the aggregates. These models directly simulate the chemical reaction with various simplifying assumptions [59, 58]. Second, mechanical models at the structure scale which are usually implemented in more general concrete durability codes accounting for various important phenomena such as creep, shrinkage and macroscopic damage [45, 2, 63, 33, 17]. Finally, microscopic scale models which can be analytical or numerical and try, by a fine description of the microstructure, to determine the mechanical consequences of the presence of a swelling gel in the concrete porosity close to the reactive sites [44, 16, 22, 64, 49, 4]. The aim of these models is at the same time to help studying the affected structures and help understanding the physics of the degradation. It seems to us that analytical models are more likely to help for structure computations since their simplicity makes more obvious the influence of their parameters. Also, we can hope that such a model could serve as a way to compute degradation in each element in structure-size codes. That is why we have focused our research and this article on such models.

Therefore, for our model to be satisfactory, it must be such that its parameters are easily identifiable in the lab on concrete samples, some by direct testing, some by fitting of expansion curves. The model must be simple enough to keep short computation times if we want to couple it with a structure-size code. It must also reproduce expansion

<footnote>*Corresponding author,
now at EDF R&D, Les Renardières, 77818 Morêt-sur-Loing Cedex,
France</footnote>

<footnote>Email address: laurent.charpin@edf.fr (Laurent Charpin)</footnote>
and reduction of mechanical properties of attacked con-
crete samples under mechanical loading.

We will first motivate our study by recalling some works
discussing the anisotropy of alkali-silica reaction. Then,
we will present our thermodynamical description of the
attacked concrete. We will continue by describing the
microporomechanics equations we use to describe the mi-
crostructure of the concrete and extract needed informa-
tion. Then, we will explain the behavior of the model on
simple examples first, and by identifying its parameters on
lab experiments of alkali-silica reaction.

## 2. Evidence and modeling of the anisotropy of ASR

The question of anisotropy of ASR, as is well explained
in [37, 50, 39], was first raised by the orientation of the
macroscopic cracks on affected structured. These cracks
tend to align with the principal compression direction in
the structure (which means that the normal to the crack is
orthogonal to the principal compression direction), which
can be due to reinforcing bars or weight. At a sample-size
level, in the lab, anisotropy is also seen, sometimes even in
the case of free-swelling of the sample. The directions of
casting, as underlined in [39], play a major role, since they
influence the orientation of the aggregates and the zones of
higher porosity around the aggregates. These phenomena
alone can induce anisotropy factors up to 4 [39].

However in this article, we are interested in the anisotropy
of swelling which is due to mechanical loads. As pointed
by Hobbs, [37], applying stresses on concrete during attack
can modify expansion a lot. The mechanical restraint also
can modify expansion, as explained and measured by Bi-
nal, Kawamura and Iwahori, and Berra et al. [7, 38, 6],
using an experimental set-up introduced by Ferraris et al.
[27].

The effect of stress and restraint was originally under-
stood as a redistribution of expansion from the compressed
to the less compressed directions. Macroscopic expansion
models by Larive and Multon [39, 50] are based on this
idea.

Some models were also developed at a macroscopic
scale to be able to compute deformations of attacked struc-
tures, for example in the works of Sellier, Capra and Bour-
nazel, or Léger et al. [64, 9, 42], generally without be-
ing too precise about what is happening at the micro-
scopic level, except in the microscopic model by Sellier
from which are extracted some information used in the
macroscopic model.

More refined mechanical models can be built, if the
description of what is happening microscopically is im-
proved. An attempt based on micromechanics was made
by Lemarchand et al. [43], where the initial load influences
the filling of initially existing cracks and hence, the expan-
sion. The model can explain the tendency of Larive's test
results under uniaxial load [39], in terms of redistribution

of the expansion at intermediate load, and reduction of the
total (volumetric) expansion at high loads.

The most advanced effort in understanding the effect of
an external load on the mechanical consequences of ASR
seems to be the numerical model of Dunant and Scrivener,
presented in [22]. In his second paper, the focus is placed
on anisotropy [23]. The author gives interesting experi-
mental results of expansion under uniaxial loading. The
cracking/damage pattern in the aggregate (where the gel is
located, in pockets) and in the cement paste is influenced
by the loading. However, the model fails at high loads
due to an artificial coalescence of cracks in 2d. This model
was recently improved by Giorla et al. [31, 30] who in-
troduced creep in the cement paste. The effect of particle
shape on the anisotropy of free expansion was successfully
reproduced as well as many features of ASR. The effects
of creep and damage (using different criteria) are investi-
gated. This model is very advanced but few comparisons
with experiments are available at the time.

Finally it seems the most complete experimental work
concerning the expansion of concrete samples under load-
ing is Multon's [50, 51], since a variety of compressive loads
are used (0 MPa, 10 MPa, 20 MPa) for the loading direc-
tion, in conjunction with a restriction of the radial ex-
pansion by steel rings of two different thicknesses (3 mm
rings, 5 mm rings). This set-up induces a triaxial stress
state which influences cracking in a more complicated way
than the usual uniaxial compression [39, 23]. Giorla also
developed a new device for ASR under triaxial load which
is very promising [31].

Our goal is, as stated in our introduction, to improve
the micromechanical modeling of ASR. We will present a
model based on micromechanics, poromechanics, and an
energy fracture criterion. Therefore, we will start our pre-
sentation of the model by explaining the thermodynamic
framework we used for our model.

## 3. Thermodynamic evolution of a concrete under-
going ASR

In this section we will present the basis of our model.
We will start by some geometrical simplifications of the
description of the concrete. Our model has been built in
the idea of studying concrete made of aggregates which
are sometimes classified in the literature as *fast reacting*
aggregates [32, 56]. For us the important point is that we
focus on aggregates in which the expansive products can
easily reach the interface with the cement paste. Therefore
the most important cracking phenomenon are considered
to be the decohesion between the aggregate and the inter-
face first, and in a second stage the propagation of cracks
in the cement paste. Some authors focus on *slow reactive*
aggregates, in which the most important cracking mecha-
nism occurs inside the aggregates [21, 22, 62, 61, 31]. This
category of aggregates is of great interest, and the method-
ology developed in this article will hopefully be used in the
future for this kind of aggregates.
### 3.1. Microscopic description of a concrete undergoing ASR

In this section we introduce the geometrical description of our attacked concrete at a given time $t$. We will call $\Omega$ the actual volume of concrete that we study.

#### 3.1.1. Description of the concrete under attack before cracking

Our concrete, at the scale at which we consider it, is composed of two main components: grains and cement paste matrix, whose properties vary due to the attack. Let us detail their properties, beginning with the grains. They are an ensemble of aggregates and sand particles (which size are extremely variable, from microns to a few centimeters). Let us call their number $N_g$ in the volume $\Omega$ ($g$ stands for grain). All grains are considered spherical. The grain $i \in N_g$ is of radius $R^i$, and its current attack degree is called $\alpha^i(t)$. The attack degree represents the proportion of the radius of the grain which has undergone attack. Therefore the shell between radii $(1 - \alpha^i(t))R^i$ and $R^i$ will be called the *attacked zone*, while the sphere of radius $(1 - \alpha^i(t))R^i$ will be called the *sound zone*, at time $t$ (Fig. 1). The *sound zone* has the mechanical properties of

<figure><img src="image_2.png" /><figcaption>Figure 1: Schematic attack of an aggregate</figcaption></figure>

a sound aggregate. It is therefore taken as a linear elastic isotropic material, which stiffness tensor is called $C_a$. The *attacked zone* has been partly dissolved by the attack of ions coming from the cement paste. We do not describe this attack but only its mechanical consequences. We assume that it transforms the sound aggregate into a porous material of porosity $\rho^i$, whose poroelastic properties are a tensor of elasticity $C_p^i$, a Biot coefficient $b_p^i$ and a Biot compliance $M_p^i$, which are defined according to the habits in poromechanics as presented by Coussy [18], except for the Biot compliance, which is taken as the inverse of the usual Biot modulus. We get the following constitutive law

for the *attacked zone* of grain $i$:

$$ \left\{ \begin{array}{l} \underline{\underline{\sigma}} = C_p^i : \underline{\underline{\varepsilon}} - b_p^i p \underline{\underline{1}} \\ (\tilde{\varphi}^i - \rho^i) = b_p^i \underline{\underline{1}} : \underline{\underline{\varepsilon}} + M_p^i p \end{array} \right. \qquad (1) $$

Where $\underline{\underline{\sigma}}$ is the stress tensor, $\tilde{\varphi}^i$ is the deformed (due to strains $\underline{\underline{\varepsilon}}$ and pressure $p$) porosity of the *attacked zone*. We assume that the attack keeps the isotropy of the aggregate, therefore the Biot coefficient $b_p^i$ is a scalar. In our model the porosity $\rho^i$ can be different in different aggregates, but is constant in time. The attack progresses by increasing of the size of the attacked zone (described by $\alpha^i(t)$) only. The poroelastic coefficients of this zone are obtained using the Mori-Tanaka estimate [5], assuming the porosity is composed of spherical cavities. This simple homogenization step is not described here.

The grain is surrounded by an interfacial transition zone (ITZ). We think it is important to take it into account for two reasons. First, attempts to predict the stiffness of concretes have shown that if micromechanics models are used considering concrete as a two-phase (aggregates and cement paste) material, the mechanical properties were overestimated. Nielsen and Monteiro [54] argued that it is better to consider three phases (aggregates, ITZ, and cement paste) to predict the mechanical properties. It confirmed mechanically observations that were made long before by Farran for example [25] about the modification of the packing of cement grains close to the aggregates leading to different composition and porosities in this zone. There is discussions about the dependence of the properties of this zone on the total aggregate volume fraction [35] based on generalizations of the self-consistent scheme of Christensen and Lo [60, 15], or on the aggregate size [36], and whether it is possible to model it as a uniform zone or not [36, 10, 53]. The second reason why we need to take the ITZ into account is related to ASR. The ITZ plays the role of a reservoir for gel escaping from the aggregate, limiting the pressure increase. Not willing to put too much detail in the ITZ which is not the focus of our study, we choose to model it as a homogeneous porous medium of constant thickness $l_c$ and porosity $\rho_{itz}$ (both independent on the grain considered). The properties are written in the same manner as for the *attacked zone*: a tensor of elasticity $C_t$, a Biot coefficient $b_t$ and a Biot compliance $M_t$ in the following constitutive law:

$$ \left\{ \begin{array}{l} \underline{\underline{\sigma}} = C_t : \underline{\underline{\varepsilon}} - b_t p \underline{\underline{1}} \\ (\tilde{\varphi}_{itz} - \rho_{itz}) = b_t \underline{\underline{1}} : \underline{\underline{\varepsilon}} + M_t p \end{array} \right. \qquad (2) $$

Where $\tilde{\varphi}_{itz}$ is the deformed (due to strains $\underline{\underline{\varepsilon}}$ and pressure $p$) porosity of the ITZ. The poroelastic coefficients relative to the ITZ are also determined using the Mori-Tanaka estimate.

Now that we have described in detail the geometry of our concrete during attack, we will proceed with the description of cracking.
### 3.1.2. Description of cracking

The concrete described in the previous paragraph is going to change morphology due to the high pressures which will develop in the pores of the *attacked zone* and *ITZ*. Let us describe the changing microstructure by a number of damage parameters for each grain *i*. First, a decohesion parameter $d^i(t)$ which can take values of 0 and 1 only and is an increasing function of time. At the beginning of the attack and as long as the aggregate remains bounded to the cement paste, $d^i(t) = 0$. When the interface at radius $R^i$, that is between the *attacked zone* and *ITZ* breaks, $d^i(t) = 1$. Second, crack size parameters we note collectively $x^i(t)$ which describe the size of the cracks developing in the cement paste. The choice of crack shape and exact meaning of these parameters will be explained in § 4.2, since some approximations will be related to specificities of the micromechanical description. The crack sizes are also increasing functions of time.

The geometrical properties of our sample has been described and we have defined the damage parameters accounting for the evolution of the microstructure. Let us sum things up defining the problem we want to solve:

(P): Find the evolution of the damage parameters $d^i(t)$ and $x^i(t)$ due to the ASR, represented by the set of attack degrees $\alpha^i(t)$.

To simplify the resolution of this problem, let us introduce aggregate *families*.

### 3.1.3. Definition of the aggregate families

An aggregate *family* is a group of aggregates which are assumed to have identical evolutions. Usually it will correspond to a given aggregate size, but a given aggregate size can also be divided in *families*. Let us call their number $N_f$. We assume that the grains of *family* *i* are of size $R^i$. The volume fraction of grains of site *i* is called $f^i$, and the number of grains per unit volume in *family* *i* writes $N^i = f^i/\frac{4\pi}{3}(R^i)^3$. Aggregates which are in the same *family* will follow the same evolution in terms of attack and cracking. The vicinity of each aggregate contain various kinds of pores. First, pores in the aggregate which appear due to the chemical attack. As explained in § 3.1.1, they are located in the *attacked zone* between radius $(1 - \alpha^i) R^i$ and radius $R^i$, and represent a volume fraction $\rho^i$ of the region where they are located. That means that in the *family* *i*, there is a total porous volume in each aggregate equal to $\rho^i \frac{4\pi}{3}(R^i)^3(1 - (1 - \alpha^i)^3)$. The overall constitutive law has been written in Eq. 1. Second, pores in the *ITZ* which, due to its small thickness $l_c$ occupy a volume $\rho^{itz} 4\pi l_c (R^i)^2$ for each aggregate. This volume is available to the gel, which we think has threshold behavior and because of that, cannot penetrate farther in smaller pores of the cement paste [14]. The constitutive law of this zone was given in Eq. 2. Third, if decohesion has occurred, the gel is also able to occupy the space between the aggregate and the *ITZ* in each site, and if cracking has occurred, it can go in the cracks. These zones have no undeformed

volume, but gel under pressure can invade them.

In our model, we assume that all pores belonging to the same *family* are at the same pressure that we call $p^i$, which means that we consider that the alkali-silica reaction develops slowly enough so that the gel can flow to every part of the porosity without dissipating large amounts of energy.

To solve problem (P) in various cases (for example for different macroscopic loading conditions on the sample), let us introduce a virtual problem.

## 3.2. Virtual problem: loading with imposed strains and attack degree

As explained in the previous section, the chemical attack is described by a number of attacked degrees. Let us now take external loadings into account, since one of the goals of our model is to be able to estimate expansions of an attacked concrete under loading.

### 3.2.1. Total energy function

Let us now describe problem (P*) which will give an energetic framework for solving problem (P). Let us first be clearer about the loading parameters of our problem. Our problem is driven by a set of time-dependent parameters $\alpha^i(t)$ which represent the attack degree for each site present in the sample. The second loading parameter from now on will be the macroscopic strain $\underline{E}$. That means the displacements $\underline{E} \cdot \underline{x}$ are prescribed on the external boundary of $\Omega$. This is a particular external loading which is convenient for the presentation, but our model is suited for any macroscopic loading as will be shown in the examples section at the end of this article (section 5).

Now that our loading parameters are defined, let us choose a given loading state $(\underline{E}, \alpha^i)$ and write formally the total energy of our sample for a set of virtual damage parameters $(d^{i*}, x^{i*})$. This total energy is the sum of the elastic energy stored in all components of the system, plus the surface (or dissipated) energy due to crack creation:

$$E_{tot}^{\underline{E}}(\alpha^i)_{i=1:Nf} (d^{i*}, x^{i*}) = E_{el}^{\underline{E}}(\alpha^i)_{i=1:Nf} (d^{i*}, x^{i*}) + E_{diss} (d^{i*}, x^{i*}) \quad (3)$$

Where the loading parameters appear as parameters for the energy functions, while the damage parameters appear as arguments. The values of these energy functions will be specified later. Let us now explain the energy criterion used to compute the crack state of the sample for a given loading.

### 3.2.2. Postulate: law of evolution of the damage (Francfort-Marigo criterion)

The damage law we choose is as follows: the real damage parameters are found by minimizing the total energy over all possible virtual damage parameters, as proposed
by Francfort and Marigo, Fedelich and Ehrlacher, or Mielke [28, 26, 47]:

$$ (d^i, x^i)^{\underline{\underline{\underline{\underline{\mathcal{E}}}}}^{(\alpha^i)}_{i=1:Nf}} = \underset{(d^{i*}, x^{i*}) \in \mathcal{A}}{\operatorname{argmin}} \left\{ E_{tot}^{\underline{\underline{\underline{\underline{\mathcal{E}}}}}^{(\alpha^i)}_{i=1:Nf}} (d^{i*}, x^{i*}) \right\} \quad (4) $$

Where $\mathcal{A}$ is the set of admissible damage states. Without entering too much details about this, it requires that damage is increasing, and that cracks surfaces are not interpenetrating.

These definitions of the energy and evolution law form problem ($\mathcal{P}^*$). We now need to specify the shape of the total energy, which means mainly of the elastic energy. However it is a difficult problem since the attack degrees influences the energy in a complex way. Let us introduce a simpler problem where we assume pressure is controlled.

### 3.3. Auxiliary problem: loading with imposed strain and pressures

We do so by means of an auxiliary problem ($\mathcal{P}^{aux}$) to the virtual problem ($\mathcal{P}^*$). In this problem, the loading variables are changed. In problem ($\mathcal{P}^*$), the loading parameters are ($\underline{\underline{\underline{\underline{\mathcal{E}}}}}^{i}$, $\alpha^i$), corresponding to the physical parameters imposed in our representation of alkali-silica reaction. While we apply a strain and the chemical attack progresses, all regions filled with gel (pores in the attacked zone, pores in the ITZ, cracked aggregates interfaces, cracks in the cement paste) in family i are under pressure $p^i$.

#### 3.3.1. New loading parameters and resulting strain

For the purpose of writing explicitly the elastic energy, we decide to use ($\underline{\underline{\underline{\underline{\mathcal{E}}}}}^{i}$, $p^i$) as loading parameters for the mechanical auxiliary problem. Following the description we have proposed of the morphology of the attacked concrete, there are two type of sites we need to distinguish, which we will call Type I and Type II, as shown on Fig. 2. The sites

<figure><img src="image_3.png" /><figcaption>Figure 2: Attacked concrete: auxiliary problem</figcaption></figure>

of Type I are those where decohesion has not occurred yet. Therefore the porosity is only that of the *attacked aggregate* and of the *ITZ*. The sites of Type II have undergone decohesion. Therefore, from the point of view of the solid

skeleton, the aggregate is not visible anymore. Only the pressure on the cavity and the crack lips remain. The $N_f$ families are therefore divided into $N_I$ families of Type I and $N_{II}$ families of Type II.

In this new mechanical problem the displacement and hence the strain is, everywhere in the structure, proportional to the loading parameters because all materials are linear elastic. Let us define two sorts of localization tensors. First, a strain localization tensor $\underline{\underline{\underline{\underline{\mathcal{A}}}}}^{i}(\underline{x})$ such that if all pressures ($p^i$) are equal to zero, the strain in the structure writes:

$$ \underline{\underline{\underline{\underline{\mathcal{E}}}}}^{i}(\underline{x}) = \underline{\underline{\underline{\underline{\mathcal{A}}}}}^{i}(\underline{x}) : \underline{\underline{\underline{\underline{\mathcal{E}}}}}^{i} \quad (5) $$

It is a fourth order tensor, with the minor symmetry but not necessarily major symmetry. Second, $N_f$ pressure localization tensors $\underline{\underline{\underline{\underline{\mathcal{A}}}}}^{i}$, $i = 1..N_f$. If only $p^j \neq 0$, while $p^i = 0, i \neq j$ and $\underline{\underline{\underline{\underline{\mathcal{E}}}}} = 0$,

$$ \underline{\underline{\underline{\underline{\mathcal{E}}}}}^{i}(\underline{x}) = \underline{\underline{\underline{\underline{\mathcal{A}}}}}^{j}_{p}(\underline{x}) p^j \quad (6) $$

As a result our strain writes, for any set of loading parameters in the auxiliary problem:

$$ \underline{\underline{\underline{\underline{\mathcal{E}}}}}^{i}(\underline{x}) = \underline{\underline{\underline{\underline{\mathcal{A}}}}}^{i}(\underline{x}) : \underline{\underline{\underline{\underline{\mathcal{E}}}}}^{i} + \sum_{j=1}^{N_f} \underline{\underline{\underline{\underline{\mathcal{A}}}}}^{j}_{p}(\underline{x}) p^j \quad (7) $$

#### 3.3.2. Elastic energy of the solid skeleton and constitutive law

We are now able to write, at least formally, the expression of the elastic energy of the skeleton by unit volume of the porous medium in the auxiliary problem for a virtual damage state, using a poromechanical description at the scale of the REV. We make the strong assumption that the matrix of Biot compliances is diagonal, that is when pressure is applied in a given porous zone, the average strain on this zone is much larger than anywhere else in the porous medium. It makes the expressions simpler and avoids keeping terms that cannot be approximated easily by micromechanics tools. Let us illustrate this assumption by showing a numerical evaluation of the value of the crossed Biot compliances in a very simple case: the case where we have 2d plane strain elastic medium with circular pores, gathered in two families, each with volume fraction $f$. These results are unpublished but were established by the authors while comparing finite element results with micromechanical estimated. This work is already published in [13], where all details concerning the simulations can be found. The materials parameters chosen were $E = 1$ for Young's modulus and $\nu = 0.25$ for Poisson's ratio. As one can see on Figure 3, the proportion between diagonal Biot compliances $M^{11} = M^{22}$ and the extra-diagonal ones $M^{12} = M^{21}$ greatly varies with volume fraction. When reaching $f = 0.2$ for both families, which means a total porosity of 0.4 and is close to what can be reached with the random sequential addition algorithm we used for these simulations, crossed and diagonal Biot compliances are close, with opposite signs. This
means that at zero macroscopic strain, the pressure in one family induces an increase of its volume of just a little larger amount than the decrease of volume of the other family. This fact is surely to be considered to precisely study the coupled evolution of damage around different reactive sites. We underline that the use of the strain localization tensor to determine the Biot compliance which is made in this paper induces the replacement of $M^{ii}$ by $\sum_j M^{ij}$ and $M^{ij}$ by zero (§ 4.1). However we stick to this assumption which is exact if the various families undergo exactly the same evolution, even if obviously we lose a coupling phenomenon by doing this. We keep in mind that this could be an idea for an improvement of the model, which will be however hard to implement due to the lack of simple micromechanical formula.

<figure><img src="image_4.png" /><figcaption>Figure 3: Evolution of two dimensional plane strain Biot compliances of a homogeneous medium with circular holes gathered in two families of the same volume fraction $f$</figcaption></figure>

All quantities except the loading parameters depend in the following equations on the virtual damage state $(d^{i*}, x^{i*})$. The elastic energy of the skeleton writes:

$$E_{el,ske}^{\underline{E}(p^i)}_{i=1:Nf} = \frac{1}{2}\underline{E} : \mathcal{C}^{hom} : \underline{E} + \frac{1}{2}\sum_{i=1}^{N_f}(p^i)^2 M^{ii} \quad (8)$$

While the constitutive equations of the overall porous medium write:

$$\left\{ \begin{array}{l} \underline{\Sigma} = \mathcal{C}^{hom} : \underline{E} - \sum_{i=1}^{N_f} p^i \underline{B}^i \\ \phi^i - \tilde{f}^i = \underline{B}^i : \underline{E} + M^{ii} p^i \quad \text{for } i \in \{1, N_f\} \end{array} \right. \quad (9)$$

Where $\underline{B}^i$ is the overall Biot coefficient of family i, describing how a pressure in a pore creates a macroscopic stress $\underline{\Sigma}$ in a macroscopically restrained sample, or equivalently how a macroscopic strain induces a pore volume change in a pore free of pressure, $M^{ii}$ is the overall Biot compliance of family i, linking the pressure in a pore family to its

volume change when the macroscopic strain is zero, and $\mathcal{C}^{hom}$ is the macroscopic stiffness tensor. The computation of these poroelastic properties will be detailed in section 4.

The dual quantities of the imposed deformation $\underline{E}$ and the imposed pressures $p^i$ are the average stress on the porous medium $\underline{\Sigma}$ and $\phi^i - \tilde{f}^i$ which is the difference between the volume fraction of pores in site i relatively to the volume of porous medium in deformed configuration $\phi^i$ and the volume fraction of pores in site i relatively to the volume of porous medium in undeformed configuration $\tilde{f}^i$, that means the volume variation of the porous space of site i relatively to the volume of porous medium.

We have completed the presentation of the auxiliary problem, let us explain how it is related to the virtual problem.

### 3.4. Back to $\mathcal{P}^*$, introduction of the gel

In this section we wish to compute the pressures in problem $\mathcal{P}^{aux}$, and introduce them in the energies given in the previous paragraph to yield energies of problem $\mathcal{P}^*$.

#### 3.4.1. Undeformed and deformed gel volume

First let us write the undeformed gel volume, from the attack degree $\alpha^i$, of site i. We define the parameter $\delta$ that we call the expansion factor of the gel. It is the volume of gel created by unit volume of dissolved aggregate. The volume fraction of new porosity due to the attack at site i per unit volume of concrete writes:

$$F^i = f^i \rho^i \left[ 1 - (1 - \alpha^i)^3 \right] \quad (10)$$

While the undeformed volume fraction of the ITZ of family i writes:

$$V_{itz}^i = 3f^i \frac{l_c \rho_{itz}}{R^i} \quad (11)$$

Since we assume that the gel has a volume $\delta$ times that of the dissolved aggregate, the total undeformed volume of gel of site i writes:

$$V_0^i = \delta F^i \quad (12)$$

Then assuming our gel behaves elastically and is characterized by a bulk modulus $K_g$, the deformed volume of gel in sites i per unit volume of porous medium $V^i$ writes:

$$V^i = \delta F^i \left( 1 - \frac{p^i}{K_g} \right) \quad (13)$$

This gel volume has to be compatible with the volume available.

#### 3.4.2. Available volume

This volume of gel has to be compared to the available volume in deformed configuration for site i by unit volume of porous medium. For families of Type I, it simply is equal to the current volume fraction of the pores $\phi^i$, which according to the constitutive law at the macroscopic scale writes:

$$\phi^i = \tilde{f}^i + \underline{B}^i : \underline{E} + M^{ii} p^i \quad (14)$$
Where the undeformed volume $\tilde{f}^i = V_0^i + V_{itz}^i$ is the sum of the newly created porosity and the porosity of the ITZ and the deformed porosity is occupied by gel if there is enough gel $\phi^i = V^i$. Therefore for families of Type I, the equation we obtain is:

$$ (F^i + \underline{\underline{B}}^i : \underline{\underline{E}} + M^{ii}p^i) + V_{itz}^i = \delta F^i \left(1 - \frac{p^i}{K_g}\right) \quad (15) $$

Which can be solved for the pressure as:

$$ p^i(\underline{\underline{E}}, \alpha^i) = \frac{(\delta - 1)F^i - \underline{\underline{B}}^i : \underline{\underline{E}} - 3f^i \rho_{itz} \frac{l_c}{R^i}}{M^{ii} + \frac{\delta F^i}{K_g}} \quad (16) $$

In families of Type II, the undeformed porosity in the constitutive equation still writes $\tilde{f}^i = V_0^i + V_{itz}^i$. But the deformed porosity is occupied both by gel and a deformed aggregate of volume $V_a^i = (f^i - F^i)(1 - \frac{p^i}{K_a})$. Hence, $\phi^i = V^i + V_a^i$, and the pressure writes:

$$ p^i(\underline{\underline{E}}, \alpha^i) = \frac{(\delta - 1)F^i - \underline{\underline{B}}^i : \underline{\underline{E}} - 3f^i \rho_{itz} \frac{l_c}{R^i}}{M^{ii} + \frac{\delta F^i}{K_g} + \frac{f^i - F^i}{K_a}} \quad (17) $$

Where $K_a$ is the bulk modulus of the sound aggregate. In Eqs. 16 and 17, $F^i$, $\underline{\underline{B}}^i$ and $M^{ii}$ depend on the attack degree $\alpha^i$, and $\underline{\underline{B}}^i$ and $M^{ii}$ depend on the damage state. If Eqs. 16 and 17 predict negative pressures, they are taken as zero instead.

### 3.5. Total energy in problem $P^*$

The combination of energies of the problem $P^{aux}$ which are expressed in terms of $(\underline{\underline{E}}, (p^i)_{i=1:N_f})$ for a damage state $(d^{i*}, x^{i*})$, combined with the expression of the pressure under the loading $(\underline{\underline{E}}, (\alpha^i)_{i=1:N_f})$ at the same damage state, will now yield the energy in terms of imposed deformation and attack degrees, which corresponds to the situation of problem $P^*$.

#### 3.5.1. Contribution of the solid skeleton

Now that we know how to compute the elastic energy of the skeleton at given macroscopic loading and pressure (Eq. 8) and how to compute the pressure for a given macroscopic loading and given attack degrees (Eqs. 16 and 17), we can write the elastic energy of the skeleton as a function of our real loading parameters, for a given virtual damage state $(d^{i*}, x^{i*})$:

$$ E_{el,ske}^{\underline{\underline{E}}, (\alpha^i)_{i=1:N_f}} = E_{el,ske}^{\underline{\underline{E}}, (p_i(\underline{\underline{E}}, \alpha^i))}_{(i=1:N_f)} \quad (18) $$

#### 3.5.2. Contribution of the gel

The elastic energy of the gel is written as follows for site $i$ per unit volume of porous material:

$$ E_{el,gel}^{\underline{\underline{E}}, (\alpha^i)_{i=1:N_f}} = \frac{1}{2} \frac{\delta F^i}{K_g} (p^i)^2 (\underline{\underline{E}}, \alpha^i) \quad (19) $$

#### 3.5.3. Contribution of the aggregates

While for Type II families, in which the aggregates are no more accounted for as a part of the skeleton, the energy stored in the aggregates of site $i$ writes:

$$ E_{el,a}^{\underline{\underline{E}}, (\alpha^i)_{i=1:N_f}} = \frac{1}{2} \frac{f^i - F^i}{K_a} (p^i)^2 (\underline{\underline{E}}, \alpha^i) \quad (20) $$

#### 3.5.4. Dissipated energy

It is the sum of crack surfaces weighted by surface fracture energies, which we take equal to $G^{dec}$ for decohesion and $G^{fiss}$ for cracks in the cement paste. Therefore it writes:

$$ E^{fiss} = \sum_{i=1}^{N^s} \frac{f^i}{\frac{4\pi}{3}(R^i)^3} \left[ d^{i*} 4\pi (R^i)^2 G^{dec} + \sum_{k=1}^{3} \pi \left( (x_k^{i*})^2 - (R^i)^2 \right) G^{fiss} \right] \quad (21) $$

We have shown the way we want to solve this evolution problem. One information is missing: how are the homogenized properties of the concrete related to the material properties of the various regions for a given damage state. We address this question in the next section.

## 4. Resolution technique: microporomechanics

We have all the tools to write the elastic energy of our attacked concrete for any virtual damage state $(d^{i*}, x^{i*})$ and any imposed strain and attack degrees $(\underline{\underline{E}}, (\alpha^i)_{i=1:N_f})$. We need to compute the poromechanical properties of the concrete appearing in the overall constitutive equations (Eq. 9). We will not detail this part but give final expressions. Details about this derivation can be found in [19, 11, 55].

### 4.1. Poroelastic coefficients

The expressions write: Introducing the volume fractions of each zone, $\alpha$ referring to Type I families, $\beta$ to Type II families: $f_a^\alpha = \Omega_a^\alpha/\Omega$, $f_p^\alpha = \Omega_p^\alpha/\Omega$, $f_t^\alpha = \Omega_t^\alpha/\Omega$, $f_v^\beta = \Omega_v^\beta/\Omega$ (see Fig. 2). The homogenized stiffness tensor writes:

$$ \begin{aligned} \mathcal{C}^{hom} = \mathcal{C}_c & - \sum_{\beta=1}^{N^{II}} f_{\beta}^v \mathcal{C}_c : \langle \mathcal{A} \rangle_{\Omega_{\beta}^v} \\ & + \sum_{\alpha=1}^{N^I} \left[ f_a^\alpha (\mathcal{C}_a - \mathcal{C}_c) : \langle \mathcal{A} \rangle_{\Omega_a^\alpha} + f_p^\alpha (\mathcal{C}_p - \mathcal{C}_c) : \langle \mathcal{A} \rangle_{\Omega_p^\alpha} \right. \\ & \left. \qquad + f_t^\alpha (\mathcal{C}_t - \mathcal{C}_c) : \langle \mathcal{A} \rangle_{\Omega_t^\alpha} \right] \end{aligned} \quad (22) $$

Where, as shown on Fig. 2, $\mathcal{C}_c$ is the elasticity tensor of the cement paste, $\mathcal{C}_a$ is the elasticity tensor of the sound aggregate on domain $\Omega_a^\alpha$, $(\mathcal{C}_p, b_p, M_p)$ and $(\mathcal{C}_t, b_t, M_t)$ are the elasticity tensor, Biot coefficient and Biot compliance
of the attacked aggregate and ITZ, respectively on domains $\Omega_p^\alpha$ and $\Omega_t^\alpha$. $< A >_D$ denotes the volume average of the field A on domain D. The Biot coefficients of the various zones containing fluid write:

$$ \left\{ \begin{array}{l} \underline{B}^{\alpha} = \underline{B}_{\underline{p}}^{\alpha} + \underline{B}_{\underline{t}}^{\alpha} = f_{\underline{p}}^{\alpha} b_{\underline{p}\underline{1}}^{\underline{1}} : < \underline{A} >_{\underline{\Omega}_{\underline{p}}^{\alpha}} + f_{\underline{t}}^{\alpha} b_{\underline{t}\underline{1}}^{\underline{1}} : < \underline{A} >_{\underline{\Omega}_{\underline{t}}^{\alpha}} \\ \underline{B}_{\underline{v}}^{\beta} = f_{\underline{v}}^{\beta}\underline{1} : < \underline{A} >_{\underline{\Omega}_{\underline{v}}^{\beta}} \end{array} \right. \qquad (23) $$

And finally the Biot compliances (neglecting crossed terms, as explained at § 3.3.2):

$$ \left\{ \begin{array}{l} M^{\alpha\alpha} = f_{\underline{p}}^{\alpha} [b_{\underline{p}\underline{1}}^{\underline{1}} : < \underline{A}^{\alpha} >_{\underline{\Omega}_{\underline{p}}^{\alpha}} + M_{\underline{p}}] \\ \qquad + f_{\underline{\alpha}}^{t} [b_{\underline{t}\underline{1}}^{\underline{1}} : < \underline{A}^{\alpha} >_{\underline{\Omega}_{\underline{\alpha}}^{\alpha}} + M_{\underline{t}}] \\ M^{\beta\beta} = f_{\underline{v}}^{\beta}\underline{1} : < \underline{A}^{\beta} >_{\underline{\Omega}_{\underline{v}}^{\beta}} = \underline{1} : \mathcal{C}_{\underline{c}}^{-1} : [\underline{B}_{\underline{v}}^{\beta} - f_{\underline{v}}^{\beta}\underline{1}] \end{array} \right. \qquad (24) $$

Let us now use micromechanics to give approximations of the averages of the various localization tensors.

## 4.2. Approximation of the localization tensors

In our approach, one of the main objectives is to keep the model very simple. Hence, we have decided to approximate the localization tensors thanks to classical micromechanics schemes. We have decided to use a scheme called the interaction direct derivative (IDD) scheme, which was introduced by Zheng and Du [68, 20]. Let us first present a dilute estimate for our problem, since it is the easiest way to write and also serves as a basis for the IDD scheme.

### 4.2.1. Dilute approach

In this approach, we assume that the heterogeneities only occupy a very small volume fraction of concrete. Hence, to estimate the localization tensors averages on the various zones, we can use simpler mechanical problems, assuming that each grain lies in an infinite cement paste. Let us first show our approach for Type I families.

Type I sites. In these families which have not undergone decohesion, we have a number of spherical shells. At the center, the sound part of the aggregate in a sphere of radius $(1 - \alpha^i)R^i$ and modulus $C_a$. Then, the attacked zones between radii $(1 - \alpha^i)R^i$ and $R^i$, and poromechanical properties $(C_p, b_p, M_p)$. Next, the ITZ between $R^i$ and $R^i + l_c$ with properties $(C_t, b_t, M_t)$. Finally, an infinite cement paste of elasticity tensor $C_c$.

As explained in § 4.1, we need to compute the localization tensors under two types of loading: strain imposed at infinity on the one hand and pressure in the porous zones on the other hand. This is conveniently achieved using Love's solution [46] which is also used by Christensen and Lo in [15]. This solution yields the displacement field under an imposed shear strain. For the pressure and the spherical part of the imposed strain, it is much easier due to the spherical symmetry of the solution. We won't give the details of this solution which can be found in the author's PhD thesis and are classical [11].

The dilute estimate of the strain localization tensor is an isotropic tensor which can hence be written for the three regions as:

$$ < \underline{A}_{\text{dil}} >_{\underline{\Omega}_{\underline{z}}^{\alpha}} = r_{\underline{z}}^{\alpha} \underline{J} + s_{\underline{z}}^{\alpha} \underline{K} \qquad (25) $$

Where $J$ and $K$ are the classical basis tensors for isotropic fourth order tensors and the index $z$ takes values $a, p$ and $t$. Second, the space averages of the pressure localization tensors can be determined simply using a spherical symmetry solution where the porous zones are under pressure. We hence obtain the averages $< \underline{A}_{\text{dil}}^{\alpha} >_{\underline{\Omega}_{\underline{p}}^{\alpha}}$ and $< \underline{A}_{\text{dil}}^{\alpha} >_{\underline{\Omega}_{\underline{t}}^{\alpha}}$. These averages allow computing the poromechanical coefficients thanks to Eqs. 22, 23, and 24.

Type II sites. For these families we haven't yet specified the shape of the cracked cavities. We need a very simple geometrical description so we can write the solution inspiring ourselves of Eshelby's solution. These sites are made of a spherical cavity on which the pressure of the gel surrounding the aggregate which has undergone decohesion is applied, and cracks in the cement paste which appear due to the high pressures remaining even after decohesion.

First, the cracks are assumed to be penny-shaped (PS) and centered on the center of the cavity. In fact, they have an annular shape. Each crack is only described by its orientation and external diameter ($x$). With this description it is possible to have $N$ cracks of sizes $x_i$, $i \in 1:N$ around each cavity, and orientations $\theta_i, \phi_i$. We want to reduce the number of cracks. Thanks to finite element calculations in 2d, we have shown that the overall properties of a medium containing randomly oriented cracks are very close to those of a medium containing cracks oriented along 2 orthogonal directions only [11, 13]. Hence, we assume that anisotropic crack distributions can be represented by orthogonal cracks of adequate sizes. Therefore, we assume that in 3d, only 3 orthogonal cracks can grow around each cavity (see Fig. 4).

This geometrical simplification leads us to try to compute the desired averages of strain localization tensors using a superposition of various Eshelby solutions. Our complex cracked cavity is hence considered as a superposition of a spherical cavity of Eshelby tensor $S_c^{sph}$ and volume fraction $f^i$, and 3 PS cracks identical to the cracks around the cavity, of aspect ratio $X_f$ and sizes $x_1, x_2, x_3$ in directions $\epsilon_1, \epsilon_2, \epsilon_3$. The associated Eshelby tensors only depend on the orientation, the aspect ratio and the outside medium (here the cement paste). We call these Eshelby tensors $S_c^{f1}, S_c^{f2}, S_c^{f3}$. This decomposition was studied in detail in 2d in [11].

An important point in the choice of the volume fraction affected to the cracks (here we use volume fraction instead of the more common crack density parameter introduced by Budiansky and O'Connell [8]). If there is one crack in each direction $\tau$ for each family $i$ of size $x_\tau^i$ and aspect ratio $X^f$, the volume fraction for each direction and family writes:

$$ f \tau^i = f^i X_f \left( \frac{x_\tau^i}{R^i} \right)^3 \qquad (26) $$
<figure><img src="image_5.png" /><figcaption>Figure 4: Spherical cavity surrounded by three orthogonal cavities</figcaption></figure>

However, we have noticed that better results are obtained when correcting this volume fraction to account for the fact that the cracks intersect the sphere. Thanks to 2d computations detailed in [11], we came up with a correction volume fraction equal to:

$$f'_{\tau} = f_{\tau}^{i} \left(\frac{R^{i}}{x_{\tau}^{i}}\right)^{2} = f^{i} X_{f} \frac{x_{\tau}^{i}}{R^{i}} \quad (27)$$

The final volume fraction we use for the crack in direction $\tau$ of family i is then:

$$f_{\tau}^{i} - f_{\tau}^{i} = f^{i} X_{f} \left(\frac{x_{\tau}^{i}}{R^{i}}\right)^{3} \left(1 - \left(\frac{R^{i}}{x_{\tau}^{i}}\right)^{2}\right) \quad (28)$$

All the Eshelby tensors can be found in Mura's book [52]. Finally the dilute localization tensors corresponding to these families can be written:

$$< \mathbf{A}^{\text{dil}} >_{\Omega_v^\beta} = [\mathbb{I} - \mathbb{S}_c^{\text{sph}}]^{-1} + \sum_{\tau=1}^{3} \frac{f_{\tau}^{i} - f_{\tau}^{i}}{f^{i}} [\mathbb{I} - \mathbb{S}_c^{f_{\tau}}]^{-1} \quad (29)$$

One could think that since the volume fractions of the cracks are much smaller than those of the cavity, the corresponding terms are negligible, but in fact some terms of the Eshelby tensors diverge when the aspect ratio goes to zero, so that their product by the volume fraction converges to a constant value [19].

We now have a dilute approximation of the averages of the localization tensors of Eqs. 22, 23, and 24. Let us explain how we obtain a better estimate with the IDD scheme.

### 4.2.2. Use of the IDD estimate

The IDD estimate was proposed by Zheng and Du [68, 20]. From our point of view, it gives good results with cracks and is easy to use, so it makes it the best candidate to improve our estimate of the localization tensors [11, 13]. With this estimate, it is easy to account for the diversity of inclusion shapes and space distribution. In many cases it is identical to classical estimates: the two-phase estimate [65] when we only have one type of ellipsoidal inclusion, Mori-Tanaka [5] when the inclusions are spherical. We have shown its efficiency thanks to 2d finite element computations on randomly generated microstructures, in cases where the interpenetration of cracks was forbidden [13]. Therefore, we don't expect our model to remain valid at high cracking states where percolation occurs. Also, this explains why we have considered the IDD estimate to perform much better than self-consistent estimates. Let us now build the IDD estimate from the dilute estimate.

**Strain localization tensor.** The IDD localization tensors are obtained from the dilute ones by multiplication by a correction term accounting for the interaction between the heterogeneities. It makes use of a Hill tensor describing the shape of the ellipsoidal inclusions $\mathbb{P}_c^k$ and another one describing the ellipsoidal spatial distribution of the inclusions $\mathbb{P}_c^{Dk}$. Let us recall that the Hill tensor of a given inclusion $k$ in medium $c$ is related to the Eshelby tensor of the same inclusion in the same medium by the relation $\mathbb{P}_c^k = \mathbb{S}_c^k : \mathbb{C}_c^{-1}$.

$$< \mathbf{A}^{\text{IDD}} >_{\Omega^i} = < \mathbf{A}^{\text{dil}} >_{\Omega^i} : \left[ \mathbb{I} - \sum_k f^k \mathbb{P}_c^{Dk} : \left( [\mathbb{C}_k - \mathbb{C}_c]^{-1} + \mathbb{P}_c^k \right)^{-1} \right]^{-1} \quad (30)$$

The correction term is the same for all inclusions. We call it $\mathbb{T}$. More precisely, the following sum must be computed:

$$\mathbb{W} = \sum_k f^k \mathbb{P}_c^{Dk} : \left( [\mathbb{C}_k - \mathbb{C}_c]^{-1} + \mathbb{P}_c^k \right)^{-1} \quad (31)$$

However, this description only works for ellipsoidal inclusions, which is not the case for our Type I families. These cannot be fully described by Hill tensors, that is why we had to use the Love solution. If we go back to the idea behind this correction term, we can write it as a function of the dilute localization tensor:

$$ \begin{aligned} < \mathbf{A}^{\text{dil}} >_{\Omega^i} &= (\mathbb{I} + \mathbb{P}_c^i : [\mathbb{C}_k - \mathbb{C}_c])^{-1} \\ &= [\mathbb{C}_k - \mathbb{C}_c]^{-1} : \left( [\mathbb{C}_k - \mathbb{C}_c]^{-1} + \mathbb{P}_c^i \right)^{-1} \end{aligned} \quad (32) $$

It writes:

$$\mathbb{W} = \sum_k f^k \mathbb{P}_c^{Dk} : [\mathbb{C}_k - \mathbb{C}_c] : < \mathbf{A}^{\text{dil}} >_{\Omega^k} \quad (33)$$
Then, if the inclusions are not ellipsoids but we know their dilute localization tensor (we do thanks to Love solution) and a Hill tensor approaching their space distribution, we know all the terms needed. Concerning the spatial distribution of the aggregates and the cracks, we keep things simple assuming that it is spherical (which might be wrong if the aggregates are flat, but then a lot of things would have to be adapted). We get:

$$W = \mathbb{S}_c^{sph} : \sum_k f^k \mathbb{C}_c^{-1} : [\mathbb{C}_k - \mathbb{C}_c] : < \mathbf{A}^{\text{dil}} >_{\Omega^k} \quad (34)$$

Therefore, we know how to write the *IDD* estimate for strain localization tensors.

Concerning pressure localization tensors, we have not been able to make the same kind of modification and have decided to keep the same expression as for the dilute estimates.

## 5. Simple examples

We have explained in the previous sections the theoretical basis of our model. Let us now demonstrate how it works by very simple examples of concrete under attack. When using our model, the most time-consuming operation is the minimization of the total energy, which requires many evaluations of the energies for different crack configurations. Hence, to make things simpler, we need to decrease the number of crack sizes over which the energy is minimized. There are two simple ways to do that. First, reduce the number of aggregate families. Second, assume the symmetry of the crack state not to have 3 independent crack sizes for each aggregate family, but just one or two. In this section and next section we present our results, which were created using a second implementation of our program (compared to [11]), in Python, in the Materials Ageing Platform developed at EDF R&D [40].

### 5.1. Imposed strain

Our theoretical derivation has been performed in the case of an imposed macroscopic strain. In this example we assume that the strain is imposed to the value:

$$\underline{E} = 10^{-4} \underline{\epsilon}_1 \otimes \underline{\epsilon}_1 \qquad (35)$$

Which corresponds to a moderate traction along axis $\underline{\epsilon}_1$. The result of this traction is that when the pressure develops around the aggregates, once decohesion has occurred due to high pressures, an anisotropic cracking takes place due to the combination of the pressure and the macroscopic loading, which corresponds precisely to the situation on real structures which we are modeling.

Let us detail on this example how the degradation of the attacked concrete develops. We will choose simple values for the parameters of the model:

<table><tr><td>E<sub>a</sub></td><td>60 GPa</td><td>K<sub>gel</sub></td><td>0.3 GPa</td></tr><tr><td>E<sub>c</sub></td><td>20 GPa</td><td>ρ<sup>1</sup></td><td>0.1</td></tr><tr><td>ν<sub>a</sub> = ν<sub>c</sub></td><td>0.25</td><td>δ</td><td>1.5</td></tr><tr><td>G<sub>c</sub><sup>dec</sup></td><td>40 J.m<sup>-2</sup></td><td>l<sub>c</sub></td><td>10 μm</td></tr><tr><td>G<sub>c</sub><sup>fiss</sup></td><td>80 J.m<sup>-2</sup></td><td>ρ<sup>itz</sup></td><td>0.3</td></tr></table> (36)

And, in this example, there will only be one aggregate family of size $R^1 = 1$ mm and volume fraction $f^1 = 0.1$.

Let us start with the observation of the evolution of the pressure. One can see on Fig. 5 that there is first a short

<figure><img src="image_6.png" /><figcaption>Figure 5: Evolution of the pressure under uniaxial traction strain. $R^1 = 1$ mm, $f^1 = 0.1$</figcaption></figure>

time during which no pressure increase is seen. This period corresponds to the filling of the *ITZ* by the gel. Then the pressure increases very fast to high values (60 MPa) in this case, but here the parameters have not been identified properly, so this value is not to be trusted. Following this increase, there are two pressure drops. These correspond to the decohesion of the aggregates and the beginning of cracking, with a slight crack size jump in this case. This interpretation of the evolution of pressure is confirmed by the observation of the crack size (Fig. 6). We see the crack jump at the same time as pressure drop. The crack size plateau which occurs around attack depths of $4.10^{-4}$ m is purely due to our model: using the Full Range *IDD* scheme, at some point the spatial distribution of the cracks is assumed to flatten. We haven't discussed this issue in detail here, but it is explained extensively in [57, 3, 13, 11]. This evolution of the shape of the crack distribution induces this temporary crack size stabilization. It could be avoided by smoothing the relationship between the distribution and the crack volume fraction. Due to this increase in crack size in direction $\underline{\epsilon}_1$, there is a drop in mechanical properties in this direction, as can be seen on Fig. 7. The first property loss is due to decohesion, which is why it affects all three directions identically. The second drop is only in direction 1 and is due to cracking. The Young's
<figure><img src="image_7.png" /><figcaption>Figure 6: Evolution of the crack size under uniaxial traction strain. $R^1 = 1$ mm, $f^1 = 0.1$</figcaption></figure>

<figure><img src="image_8.png" /><figcaption>Figure 7: Evolution of the solid skeleton Young's modulus under uniaxial traction strain. $R^1 = 1$ mm, $f^1 = 0.1$</figcaption></figure>

modulus of the skeleton diminishes to a very low value, while the undrained modulus, which takes the value of the gel and aggregate in the cavities into account diminishes much less (Fig. 8). We also display the Biot coefficients (Fig. 9) and Biot compliance (Fig. 10). Their evolution are similar, but the Biot coefficient becomes anisotropic as the stiffness tensor, which is due to the crack orientation. The values of the anisotropic Biot coefficient might seem large (when it reaches 1, it means a stress inside the pore network is fully transmitted at the macroscopic scale, if strain is prevented), but are consistent with the Young's modulus values. When the Biot coefficient reaches 0.7 in direction 1, the drained Young's modulus is divided by almost a factor of 4. This large increase of the Biot coefficient (Fig. 9) and decrease of the Young's modulus (Fig. 7) are due for the isotropic part to decohesion (spherical cavities of volume fraction 0.1 here yield a Biot coefficient of

<figure><img src="image_9.png" /><figcaption>Figure 8: Evolution of the undrained Young's modulus under uniaxial traction strain. $R^1 = 1$ mm, $f^1 = 0.1$</figcaption></figure>

0.2 with the IDD estimate), and for the anisotropic part, to a high concentration of cracks normal to direction 1, not in terms of volume fractions which are always negligible for cracks, but in terms of crack density parameter which is not computed here (see Fig. 6), [8].

<figure><img src="image_10.png" /><figcaption>Figure 9: Evolution of Biot coefficients under uniaxial traction strain. $R^1 = 1$ mm, $f^1 = 0.1$</figcaption></figure>

Finally, during this pressure build-up and crack creation under fixed displacement, stresses develop (Fig. 11). While at the beginning of attack the largest stress is in the direction of the imposed traction, cracking creates a high compressive stress which is most important in the cracking direction. It is also interesting to take a look at the evolution of energies during the evolution of the concrete (Fig. 12). The total energy is continuous, even at moments where the crack state is discontinuous. The dissipated energy jumps, and the elastic energy drops exactly of the same amount at decohesion and at the crack jump. The potential and elastic energies are identical be-
<figure><img src="image_11.png" /><figcaption>Figure 10: Evolution of the Biot compliances under uniaxial traction strain. R¹ = 1 mm, f¹ = 0.1</figcaption></figure>

<figure><img src="image_12.png" /><figcaption>Figure 11: Evolution of stresses under uniaxial traction strain. R¹ = 1 mm, f¹ = 0.1</figcaption></figure>

cause there are no imposed forces in this example. We have then shown that our model is able to determine the evolution of cracking depending on an external loading of the imposed strain type. Let us now give an example with imposed stresses.

## 5.2. Imposed stress

We choose to impose an unidirectional compression stress.

$$ \Sigma = -5 \text{ MPa} \quad \varepsilon_1 \otimes \varepsilon_1 \qquad (37) $$

All parameters remain the same as given in Eq. 36, but the aggregate volume fraction is higher in this example (R¹ = 1 mm, f¹ = 0.3). In this example, a compression is imposed to the sample along the first direction. It is therefore more likely to see cracks appear along axes 2 and 3. As expected, this is what happens, as can be seen on

<figure><img src="image_13.png" /><figcaption>Figure 12: Evolution of energies under uniaxial traction strain. R¹ = 1 mm, f¹ = 0.1</figcaption></figure>

<figure><img src="image_14.png" /><figcaption>Figure 13: Evolution of crack sizes under uniaxial compression stress. R¹ = 1 mm, f¹ = 0.3</figcaption></figure>

Fig. 13. This cracking under imposed stress induces important strains in the uncompressed directions (Fig. 14). We have hence shown that our model is also able to deal with imposed stresses. The equations relative to this case have not been detailed because they are very close to those of the imposed strain case. The only difference is that when computing the total energy, the potential energy includes a term taking account of the work of the imposed stresses.

Finally let us present a more complicated test which was used by Multon in his thesis [50].

## 5.3. Multon's test

In this test, cylindrical concrete samples are submitted to various compressive loads along their longitudinal axis, while their radial expansion is restrained by steel rings [50]. This set-up is interesting because of the complex three-dimensional stress state which is induced. The idea
<figure><img src="image_15.png" /><figcaption>Figure 14: Evolution of strain under uniaxial compression stress. $R^1 = 1$ mm, $f^1 = 0.3$</figcaption></figure>

is that because of the longitudinal compressive load, cracking starts in the radial directions. This induces tension in the rings and eventually high stresses. Therefore, at some point radial cracking can be stopped and cracking can appear in the longitudinal direction where compression is still applied. This phenomenon of cracking reorientation is very interesting and a good way to test the qualitative behavior of our model. This behavior is shown here by using the same set of parameters as previously (Eq. 36), again with only one aggregate size ($R^1 = 1$ mm and $f^1 = 0.3$). As in Multon's test, the sample's radius is 6.5 cm and the steel ring's Young's modulus is 193 GPa. The rings cover all the height of the sample but are not connected. Hence they only impact the radial stress in the sample. In the present case the steel rings are of radius 3 mm and the compressive load is 10 MPa.

Let us take a look at the development of cracking in the sample. It starts with cracks in directions 2 and 3, which are the directions impacted by the rings (Fig. 15). At some point, compressive stresses become much more important in the radial direction that in the longitudinal direction (Fig. 16). Cracking changes orientation and stresses in the longitudinal and radial directions progressively equilibrate.

### 5.4. Conclusion on the qualitative behavior of the model

In this section we have shown some examples of the ability of the model to predict cracking of an attacked concrete under various loadings. However, no comparison with experiments has been shown yet. This will be the focus of next section.

## 6. Comparison with experiments

In order to test our model, we need to compare its results to experimental data. We have chosen to use Multon's

<figure><img src="image_16.png" /><figcaption>Figure 15: Evolution of crack sizes in Multon's test with 3 mm rings and 10 MPa compressive load. $R^1 = 1$ mm, $f^1 = 0.3$</figcaption></figure>

<figure><img src="image_17.png" /><figcaption>Figure 16: Evolution of stresses in Multon's test with 3 mm rings and 10 MPa compressive load. $R^1 = 1$ mm, $f^1 = 0.3$</figcaption></figure>

data [50, 51]. Until more data is available, particularly on triaxial tests in which the radial and axial loads are independently imposed (such as in the testing device developed by Giorla [31]), it seems to us it is the most complete data set for ASR under loading. Let us give some details about the test.

### 6.1. Multon's test

Multon performed 9 different loading cases on reactive samples, and 3 loading on non reactive sample, in order to correct the strains of the reactive sample from shrinkage and creep. Even if this correction is not perfect, we follow Multon's approach and use the corrected expansion curves. A better solution would be to simulate shrinkage and creep, which is done for example in Grimal's model [33, 34] in his phenomenological model and which was developed by Giorla in a microscale model [30, 31]. Since
creep of cementitious materials is faster at small scales, we think that creep plays a role very locally on crack propagation. This role is neglected here, and could be investigated in a further version of the model where the cement paste would be considered as viscoelastic, as is done by Giorla through finite element simulations.

The samples are vertical cylinders of 13 cm diameter and 24 cm length. They are sealed with aluminum adhesive foils. A creep frame is used to submit them to 10 MPa and 20 MPa loads. The steel rings are either 3 mm or 5 mm thick and 1 cm high, covering all the lateral surface of the samples but not touching each other. The temperature is kept at 38 °C. The 9 tests combine 3 radial restraints: rings of 0 mm, 3 mm, and 5 mm, and three axial loadings: 0 MPa, 10 MPa, and 20 MPa. In our model we assume that the volume fraction of reactive silica in the reactive aggregates is $\rho = 0.1$ [49].

The particle size distribution is shown on Fig. 17. There is 23 % of non reactive aggregates (represented with a radius of 1 mm, but this radius does not matter since they don't undergo degradation). The rest of the aggregates are reactive. They are limestone crushed aggregates with silica veins, which can react. We think the created gels can easily migrate to the ITZ because of the veins. Therefore, our model seems well suited for this aggregate since most fracture phenomena will take place in the cement paste. All reactive aggregates are rather large, starting from 2.5 mm to 2.8 cm. We here underline that our model is able to deal with such variety of particle sizes, very large computations would be required to deal with a representative elementary volume of such a material.

<figure><img src="image_18.png" /><figcaption>Figure 17: Aggregate size distribution.</figcaption></figure>

## 6.2. Comments on the test results

The strain and stress computed from the measurements are given on Figs. 18 and 19. Looking at the strains (Fig. 18), we see that the free expansion is quite anisotropic, which was also mentioned by Comi and Grimal who worked

on these tests results as well [34, 17]. The axial strain ($E_{zz}^{multon}$) decreases with increasing axial load, except when looking at the couples (3 mm, 10 MPa; 3 mm, 20 MPa) and (5 mm, 10 MPa; 5 mm, 20 MPa) and tends to increase with increasing rings thickness, but it is not very pronounced when looking at couples (3 mm, 10 MPa; 5 mm, 10 MPa) and (3 mm, 20 MPa; 5 mm, 20 MPa). The radial strain slightly increases with the axial load and decreases with the rings thickness. It is difficult to interpret these strains, possibly because of measurement errors, as mentioned by Multon concerning the radial strain, or of the correction for creep and shrinkage from experiments with non reactive aggregates and therefore, different stress states. The computation of the stress is more strainghforward (Fig. 19) and it seems easier to interpret: the radial stress increases with increasing axial load and rings thickness.

Let us use these experiments to identify the parameters of our model.

## 6.3. Parameter identification

We have performed the identification of some of the parameters of our model on Multon's data. First, we have considered that some parameters could be fixed because their values can be measured:

$$ 
\begin{array}{|l|l|l|l|}
\hline
E_a & 60 \text{ GPa} & \rho^{itz} & 0.3 \\
E_c & 20 \text{ GPa} & \rho^i & 0.1 \\
\nu_a = \nu_c & 0.25 & X & 10^{-6} \\
E_s & 193 \text{ GPa} & \nu_s & 0.25 \\
\hline
\end{array}
\qquad (38) $$

We have fixed the aspect ratio of the cracks since the model is insensitive to this parameter as long as it is very small. We then assume that the fracture energy of the cement paste is twice that of the interface $G^{fiss} = 2G^{dec}$, which is consistent with values given in [1]. We also introduce a time scale parameter $\eta$, which links the attack depth to the physical time:

$$ d(t) = \eta \sqrt{t} \qquad (39) $$

We assume the attack depth is the same for all aggregates and scales with the square root of time, which was also our choice in [12], following [48, 29]. The reason for that it that the progress of the attack is linked to diffusion processes. Therefore, the parameters of our model which need to be identified are $K_g, l_c, \delta, \eta, G^{dec}$.

The optimization of the parameters is done only on the expansion curves of two experiments: (0 mm, 0 MPa) and (3 mm, 10MPa). Only the strain curves were used to identify the parameters. Then all experiments were simulated, and stress were also computed for comparison with the experimental results. The obtained set of parameters is the following:

$$ 
\begin{array}{|l|l|l|l|}
\hline
K_g & 0.303 \text{ GPa} & \eta & 2.67 \ 10^{-5} \text{ m}/\sqrt{d} \\
l_c & 3.88 \ \mu\text{m} & G^{dec} & 338 \text{ J.m}^{-2} \\
\delta & 1.58 & G^{fiss} & 2G^{dec} \\
\hline
\end{array}
\qquad (40) $$

Where the time scale parameters $\eta$ is written in meters per square root of days. Let us first comment on this value of
$\eta$. With this value, 400 days of experiment correspond to an attack depth of 0.53 mm. Compared to the size of the aggregates (2.5 mm to 2.8 cm, see Fig. 17), we see that the smallest aggregates are significantly attacked, while for the larger ones, the attack is rather superficial.

The fit of the strains on experiments (0 mm, 0 MPa) and (3 mm, 0 MPa) give satisfactory results (Fig. 18). Our model cannot reproduce the anisotropy of the free expansion. One way to do it would be to use non spherical aggregates, as shown by Giorla [31]. This improvement could be introduced in a more advanced version of the model but was not attempted yet. The other strain curves are also well reproduced for loads of 0 and 10 MPa, but the model obviously fails at 20 MPa. The decohesion does not occur at such loads, because our model currently does not allow for partial decohesion, and also forbids interpenetration of cracked surfaces. In the case of high loads (20 MPa here), the pressure is not sufficient to prevent interpenetration of the aggregate and its surrounding cement paste after decohesion. Therefore decohesion is prevented. This could be improved by modeling partial decohesion of the aggregate but would require important improvements of its micromechanical description. This discrepancy of the model is also seen on the 10 MPa tests where the brutality of decohesion induces important strain jumps. This is also due to the low number of aggregate size and the total absence of variability of the material parameters. One could point out that the expansion seems to slow down at the end of the test period (450 days), which is less visible in the predictions than in the tests results. Indeed, the depletion of reactants such as alkali is not modeled here, while it might have occurred in the tests. Some models at the same scale consider this point [48].

The results are also satisfactory in terms of stress for most cases (Fig. 19), even if stresses are usually overestimated. This can be linked to the fact that creep is not taken into account in our model, but only by correcting the expansion curves using experiments on non-reactive samples.

Finally, let us comment on the identified parameters (Eq. 40). The gel bulk modulus is roughly a tenth of that of water, and two orders of magnitude below that of C-S-H. It is also much lower that the modulus chosen by Giorla in his simulations [31] and those measured by Leeman and Lura on samples from a 45 years old dam [41]. We think that the kind of gel that forms outside the aggregates is very different from the swelling pockets found in more slow-reactive concrete as modeled by Giorla. Considering the gel as a porous solid, this bulk modulus corresponds to its drained stiffness, since the experiment is long compared to the time of water transport in gels. Therefore we think the bulk modulus is in the right order of magnitude. We think that the expansion factor $\delta$ and the ITZ thickness are also realistic.

In contrast to that, the fracture energies $G^{dec}$ and $G^{fiss}$ are an order of magnitude above the values reported in the literature, for example by Wittmann or Alexander [1, 67].

We think that the reason for that is that our model creates fracture patterns that are more oriented than in real concretes, where crack directions are highly influenced by the relative position of aggregates. To compensate for this intrinsic lack, the identification lead to increased values of the fracture energies, allowing cracking to occur at higher pressures and distributing cracks more isotropically.

We think that the model is well suited for loads below 10 MPa. The drawbacks of our model, such as the brutality of decohesion and the high fracture energies, do not invalidate it, but give hints to improve the model. Another drawback which was not underlined is the overestimation of the decrease of the Young's moduli with our model. We haven't shown these results, but these moduli are roughly divided by four, which is too much compared to experimental data. We think the main reason is that in our model, the decohesion removes completely the shear stiffness of the aggregates, while in real concrete this stiffness would contribute to the overall stiffness in compression. Another interesting output of the model would be to compute the resistance of the concrete at different attack stages, as done by Esposito and Hendriks [24]. This would require that the failure mechanism during a resistance test be close to that modeled here for alkali-silica reaction, which needs to be discussed.

## 7. Conclusion

In this article our main goal was to show that it is possible to build simple models for alkali-silica reaction based on micromechanics. To do so, we have chosen to represent the behavior of a concrete containing aggregates such that most of the damage occurs at the interface between the aggregates and the cement paste and in the cement paste. We have presented the micromechanical framework as well as the energy criterion that allow computing cracking around the aggregates under macroscopic loading. We have been able to identify the parameters of our model on triaxial experiments by Multon. The model reproduces the interesting phenomenon of cracking reorientation under loading and behaves well up to compressive loads of 10 MPa. The model was implemented in the Materials Ageing Platform at EDF R&D.

Some improvements are needed for the model. For example, the description of decohesion is not refined enough, inducing brutal evolutions of strain in some cases, and forbidding swelling in any direction in cases where the concrete is highly compressed.

To make the model more practical, now that we are confident that that kind of model can be built and that parameters can be identified on experiments, we plan to build a version that would describe degradation for slow-reactive aggregates. That would require computing damage inside the aggregates rather than outside. On the one hand, from a micromechanical point of view, such a model seems easier to build since interaction between aggregates
<figure><img src="image_19.png" /><figcaption>Figure 18: Identification of parameters on Multon's experiments. Strain. Grey plots: experiment used for the identification.</figcaption></figure>

will be diminished if damaged is located inside the aggregates. On the other hand, on larger time-scales, the neglect of creep would become a bigger approximation. This path surely needs to be explored.

## Acknowledgements

This work was funded by the 'Chaire pour la durabilité des matériaux et des structures pour l'énergie' at ENPC.

## References

[1] M.G. Alexander. The effects of ageing on the interfacial transition zone in concrete. In J. C. Maso, editor, *Interfacial transition zone in concrete (RILEM report 11)*, pages 150–174. E & FN SPON, 1996. p. 14 and 15

[2] F. Bangert, D. Kuhl, and G. Meschke. Chemo-hydro-mechanical modelling and numerical simulation of concrete deterioration caused by alkali-silice reaction. *International journal for numerical and analytical methods in geomechanics*, 28:689 – 714, 2004. p. 1

[3] B. Bary. Estimation of poromechanical and thermal conductivity properties of unsaturated isotropically microcracked cement pastes. *International Journal for Numerical and Analytical Methods in Geomechanics*, 35(14):1560–1586, 2011. p. 10

[4] Z. P. Bažant, G. Zi, and C. Meyer. Fracture mechanics of ASR in concretes with waste glass particles of different sizes. *Journal of Engineering Mechanics*, 126-3:226–232, 2000. p. 1

[5] Y. Benveniste. A new approach to the application of MoriTanaka's theory in composite materials. *Mechanics of Materials*, 6:147–157, 1987. p. 3 and 9

[6] M. Berra, G. Faggiani, T. Mangialardi, and A.E. Paolini. Influence of stress restraint on the expansive behaviour of concrete affected by alkali-silica reaction. *Cement and Concrete Research*, 40(9):1403–1409, 2010. p. 2

[7] A. Binal. Investigation of mixing ratio and restraint condition effects on alkali-silica gel pressure in mortar by using a uni-axial alkali-silica gel pressure measuring device (ASGPM-D). *Construction and Building Materials*, 21(6):1218 – 1228, 2007. p. 2

[8] B. Budiansky and R. J. O'Connell. Elastic moduli of a cracked solid. *International Journal of Solids and Structures*, 12(2):81–97, 1976. p. 8 and 11

[9] B. Capra and J.-P. Bournazel. Modeling of induced mechanical effects of alkali-aggregate reaction. *Cement and Concrete Research*, 28:251–260, 1998. p. 2

[10] S. Caré and E. Hervé. Application of a n-phase model to the diffusion coefficient of chloride in mortar. *Transport in Porous Media*, 56:119–135, 2004. p. 3

[11] L. Charpin. *Modèle micromécanique pour l'étude de la réaction alcali-silice*. PhD thesis, Université Paris Est, France, 2013. p. 7, 8, 9, and 10

[12] L. Charpin and A. Ehrlacher. A computational linear elastic fracture mechanics-based model for alkali-silica reaction. *Cement and Concrete Research*, 42(4):613 – 625, 2012. p. 14

[13] L. Charpin and A. Ehrlacher. Estimating the poroelastic properties of cracked materials. *Acta Mechanica*, pages 1–19, 2014. p. 5, 8, 9, and 10

[14] L. Charpin and A. Ehrlacher. Simplified model for the transport of alkali-silica reaction gel in concrete porosity. *Journal of Advanced Concrete Technology*, 12(1):1–6, 2014. p. 4

[15] R.M. Christensen and K.H. Lo. Solutions for effective shear properties in three phase sphere and cylinder models. *Journal of the Mechanics and Physics of Solids*, 27(4):315 – 330, 1979. p. 3 and 8
<figure><img src="image_20.png" /><figcaption>Figure 19: Identification of parameters on Multon's experiments. Computed stresses (not used for the identification of the parameters of the model).</figcaption></figure>

[16] I. Comby-Peyrot. *Development and validation of a 3D computational tool to describe damage and fracture due to alkali-silica reaction in concrete structures*. PhD thesis, École des Mines de Paris, France, 2006. p. 1

[17] C. Comi, R. Fedele, and U. Perego. A chemo-thermo-damage model for the analysis of concrete dams affected by alkali-silica reaction. *Mechanics of Materials*, 41:210 – 230, 2009. p. 1 and 14

[18] O. Coussy. *Poromechanics*. Wiley, 2004. p. 3

[19] L. Dormieux, D. Kondo, and F.-J. Ulm. *Microporomechanics*. Wiley, 2006. p. 7 and 9

[20] D.-X. Du and Q.-S. Zheng. A further exploration of the interaction direct derivative (IDD) estimate for the effective properties of multiphase composites taking into account inclusion distribution. *Acta Mechanica*, 157:61 – 80, 2002. p. 8 and 9

[21] C. F. Dunant. *Experimental and Modelling Study of the alkali-silica reaction in concretes*. PhD thesis, EPFL, Switzerland, 2009. p. 2

[22] C. F. Dunant and K. L. Scrivener. Micro-mechanical modelling of alkali-silica-reaction-induced degradation using the AMIE framework. *Cement and Concrete Research*, 40:517 – 525, 2010. p. 1 and 2

[23] C. F. Dunant and K. L. Scrivener. Effects of uniaxial stress on alkali-silica reaction induced expansion of concrete. *Cement and Concrete Research*, 42(3):567 – 576, 2012. p. 2

[24] R. Esposito and M.A.N. Hendriks. Multiscale material model for asr-affected concrete structures. In *in Onate, E. and Owen, D.R.J. and Peric D. and Suarez, B. (Eds.). XII International Conference on Computational Plasticity. Fundamentals and Applications (COMPLASXII)*, Barcelona, Spain, 2013. p. 15

[25] J. Farran. Contribution minéralogique à l'étude de l'adhérence entre les consituants hydratés des ciments et les matériaux

enronés. *Revue des Matériaux de Construction*, 490-491-492, 1956. p. 3

[26] B. Fedelich and A. Ehrlacher. Sur un principe de minimum concernant des matériaux à comportement indépendant du temps physique. *Comptes Rendus de l'Académie des Sciences II*, 302:1311 – 1394, 1993. p. 5

[27] C.F. Ferraris, J.R. Clifton, E.J. Garboczi, and F.L. Davis. Stress due to alkali-silica reactions in mortars. *Mechanisms of chemical degradation of cement-based systems*, pages 75–82, 1997. p. 2

[28] G.A. Francfort and J.-J. Marigo. Revisiting brittle fracture as an energy minimization problem. *Journal of the Mechanics and Physics of Solids*, 46(8):1319 – 1342, 1998. p. 5

[29] Y. Furusawa, H. Ohga, and T. Uomoto. An analytical study concerning prediction of concrete expansion due to alkali-silica reaction. *Durability of concrete, Proceedings of the third international conference*, SP-145:757–780, 1994. p. 14

[30] A. B. Giorla, K.L. Scrivener, and C.F. Dunant. Finite elements in space and time for the analysis of generalized visco-elastic materials. *International Journal for Numerical Methods in Engineering*, 6:454 – 472, 2014. p. 2 and 13

[31] A. B. Giorla. *Modelling of Alkali-Silica Reaction under Multi-Axial Load*. PhD thesis, EPFL, Switzerland, 2013. p. 2, 13, and 15

[32] P.E. Grattan-Bellew. *Alkali-silica reaction - Canadian experience, in: the alkali aggregate reaction in concrete*. Ed. R.N Swamy, Blackie, London, UK, 1992. p. 2

[33] E. Grimal, A. Sellier, Y. Le Pape, and E. Bourdarot. Creep, shrinkage, and anisotropic damage in alkali-aggregate reaction swelling mechanism - Part I: Constitutive model. *ACI Materials Journal*, 105:3-227 – 231, 2008. p. 1 and 13

[34] E. Grimal, A. Sellier, Y. Le Pape, and E. Bourdarot. Creep, shrinkage, and anisotropic damage in alkali-aggregate reaction
swelling mechanism - Part II: Identification of model parameters and application. *ACI Materials Journal*, 105:3-236 – 242, 2008. p. 13 and 14

[35] Z. Hashin and P.J.M. Monteiro. An inverse method to determine the elastic properties of the interphase between the aggregate and the cement paste. *Cement and Concrete Research*, 32(8):1291 – 1300, 2002. p. 3

[36] E. Hervé, S. Caré, and J.P. Seguin. Influence of the porosity gradient in cement paste matrix on the mechanical behavior of mortar. *Cement and Concrete Research*, 40(7):1060 – 1071, 2010. p. 3

[37] D.W. Hobbs. *Alkali-silica reaction in concrete*. Thomas Telford, London, 1988. p. 2

[38] M. Kawamura and K. Iwahori. ASR gel composition and expansive pressure in mortars under restraint. *Cement and Concrete Composites*, 26(1):47 – 56, 2004. p. 2

[39] C. Larive. *Apports combinés de l'expérimentation et de la modélisation à la compréhension de l'alcali-réaction et de ses effets mécaniques*. PhD thesis, École Nationale des Ponts et Chaussées, France, 1997. p. 2

[40] F. Latourte, C. Toulemonde, J.-F. Rit, J. Sanahuja, N. Rupin, J. Ferrari, H. Perron, E. Bosso, A. Guery, J.-M. Proix, and B. Wattrisse. The materials ageing plateform: towards a toolbox to perform a wide range of research studies on the behavior of industrial materials. In *PhotoMechanics 2013*, Montpellier, France, 2013. p. 10

[41] A. Leemann and P. Lura. E-modulus of the alkali-silica-reaction product determined by micro-indentation. *Construction and Building Materials*, 44(0):221 – 227, 2013. p. 15

[42] P. Léger, P. Côté, and R. Tinawi. Finite element analysis of concrete swelling due to alkali-aggregate reactions in dams. *Computers and Structures*, 60(4):601 – 611, 1996. p. 2

[43] E. Lemarchand, L. Dormieux, E. Fairbairn, and F. Ribiero. A micromechanical approach to ASR-induced damage in concrete. *Third Biot Conference, Norman, Oklahoma, USA*, 2005. p. 2

[44] E. Lemarchand, L. Dormieux, and F.-J. Ulm. Micromechanics investigation of expansive reactions in chemoelastic concrete. *Philosophical Transactions of the Royal Society A*, 363:2581 – 2602, 2005. p. 1

[45] K. Li. *Modélisation chimico-mécanique du comportement des bétons affectés par la réaction d'alcali-silice et expertise numérique des ouvrages d'art dégradés*. PhD thesis, École Nationale des Ponts et Chaussées, France, 2002. p. 1

[46] A. E. H. Love. *A Treatise on the Mathematical Theory of Elasticity, Fourth Edition*. Cambridge University Press Warehouse, 1927. p. 8

[47] A. Mielke. Energetic formulation of multiplicative elastoplasticity using dissipation distances. *Continuum Mechanics and Thermodynamics*, 15(4):351–382, 2003. p. 5

[48] S. Multon, M. Cyr, A. Sellier, P. Diederich, and L. Petit. Effects of aggregate size and alkali content on ASR expansion. *Cement and Concrete Research*, 40:508 – 516, 2010. p. 14 and 15

[49] S. Multon, A. Sellier, and M. Cyr. Chemo-mechanical modeling for prediction of alkali silica reaction (ASR) expansion. *Cement and Concrete Research*, 39:490 – 500, 2009. p. 1 and 14

[50] S. Multon. *Evaluation expérimentale et théorique des effets mécaniques de l'alcali-réaction sur des structures modèles*. Laboratoire Central des Ponts et Chaussés, 2004. p. 2, 12, and 13

[51] S. Multon and F. Toutlemonde. Effect of applied stresses on alkali-silica reaction-induced expansions. *Cement and Concrete Research*, 36(5):912 – 920, 2006. p. 2 and 13

[52] T. Mura. *Micromechanics of defects in solids*. Martinus Nijhoff Publishers, 1991. p. 9

[53] J.C. Nadeau. Water-cement ratio gradients in mortars and corresponding effective elastic properties. *Cement and Concrete Research*, 32(3):481 – 490, 2002. p. 3

[54] A. Nielsen, F. Gottfredsen, and F. Thøgersen. Development of stresses in concrete structures with alkali-silica reactions. *Materials and structures*, 26:152–158, 1993. p. 3

[55] B. Pichler and L. Dormieux. Cracking risk of partially saturated porous media - Part I: Microporoelasticity model. *International journal for numerical and analytical methods in geomechanics*, 34:135 – 157, 2009. p. 7

[56] J.M. Ponce and O.R. Batic. Different manifestations of the alkali-silica reaction in concrete according to the reaction kinetics of the reactive aggregates. *Cement and Concrete Research*, 36:1148 – 1156, 2006. p. 2

[57] P. Ponte Castañeda and J.R. Willis. The effect of spatial distribution on the effective behavior of composite materials and cracked media. *J. Mech. Phys. Solids*, 43 - 12:1919–1951, 1995. p. 10

[58] S. Poyet, A. Sellier, B. Capra, G. Foray, J.-M. Torrenti, H. Cognon, and E. Bourdarot. Chemical modelling of alkali-silica reaction: influence of the reactive aggregate size distribution. *Materials and Structures*, 40:229 – 239, 2007. p. 1

[59] S. Poyet. *Étude de la dégradation des ouvrages en béton atteints par la réaction alcali-silice : Approche expérimentale et modélisation numérique multi-échelles des dégradations dans un environnement hydro-chemo-mécanique variable*. PhD thesis, Université de Marne-La Vallée, France, 2003. p. 1

[60] E.D. Ramesh, G. and Sotelino and W.F. Chen. Effect of transition zone on elastic moduli of concrete materials. *Cement and Concrete Research*, 26(4):611 – 622, 1996. p. 3

[61] H. W. Reinhardt and O. Mielich. Mechanical properties of concretes with slowly reacting alkali sentitive aggregates. *Proceedings of the 14th International Conference on Alkali-Aggregate Reaction (ICAAR)*. Austin, Texas, 2012. p. 2

[62] H.W. Reinhardt and O. Mielich. A fracture mechanics approach to the crack formation in alkali-sensitive grains. *Cement and Concrete Research*, 41:255 – 262, 2011. p. 2

[63] B. A. Schrefler, F. Pesavento, L. Simoni, D. Gawin, and M. Wyrzykowski. Modelling alkali-silica reaction under non-isothermal conditions in partially saturated cementitious materials. *Proceedings of the Mechanics and Physics of Porous Solids conference, Champs-sur-Marne, France*, 2011. p. 1

[64] A. Sellier. *Modélisations probabilistes du comportement de matériaux et de structures en génie civil*. PhD thesis, Ecole Normale Supérieure de Cachan, France, 1995. p. 1 and 2

[65] Z. Shige and W. Tzuchiang. Effective elastic moduli of inhomogeneous solids by embedded cell model. *The Chinese society of theoretical and applied mechanics*, vol. 15, no. 4:334 – 343, 1999. p. 9

[66] T.E. Stanton. Expansion of concrete through reaction between cement and aggregate. *Proceedings of the American Society of Civil Engineering*, 66-10:1781 – 1811, 1940. p. 1

[67] F.H. Wittmann. Crack formation and fracture energy of normal and high strength concrete. *Sādhanā*, 27 - 4:413 – 423, 2002. p. 15

[68] Q.-S. Zheng and D.-X. Du. An explicit and universally applicable estimate for the effective properties of multiphase composites which accounts for inclusion distribution. *Journal of the Mechanics and Physics of Solids*, 49:2765 – 2788, 2001. p. 8 and 9