Mechanical damage, chemical damage and permeability in quasi-brittle cementitious materials

ABSTRACT The serviceability of concrete structures is a coupled problem in which fracture and damage are coupled with several environmental attacks. In this paper, we start with the description of chemo-mechanical damage and the case of calcium leaching which is relevant to waste containment vessels. The second example of implementation of continuum damage models discussed deals with coupled damage permeability effects. In the case of diffuse damage, the material permeability is controlled by the decrease of average stiffness due to micro-cracking. After a macro-crack has formed, permeability is controlled by a power function of the crack opening (Poiseuille flow). A relationship between permeability and damage, consistent with the two above asymptotic cases, is defined.


Introduction
The serviceability concrete structures used for instance in waste containment vessels or electric production facilities is a coupled problem. These concrete structures are subjected to various loads (mechanical, hydraulic…), and at the same time to chemical aggressions. In order to determine the service life of such structures, prediction of the material response and structural failure with time is required. Typically, coupled chemical-mechanical damage involves two sources of damage: one due to mechanical loads, and another one due to environmental actions. In the following, we are going to give an example related to the safety of nuclear waste disposal systems. Then, we will discuss issues related to tightness of confinement vessels involving a coupled effect between cracking and damage on one hand, and permeability and mass transfer on the other hand.

Basics of continuum damage modelling
Fracture and damage mechanics are two correlated theories which describe the same degradation process in materials (Mazars and Pijaudier-Cabot, 1996). In a computational setting, it is more convenient to deal with continuum damage instead of having to describe the shape, length and orientation of cracks (including their interaction). For this reason, we are going to focus on a continuum damage approach, the scalar damage model developed by Mazars (1984), enriched by a nonlocal approach. The same developments could be performed on a basis of a more complex description of damage without major difficulties.
Damage is assumed to be isotropic and it produces a degradation of the elastic stiffness of the material through a variation of the Young's modulus: where V and H are the Cauchy stress tensor and the strain tensor respectively. C is the fourth order tensor of elastic moduli. The damage variable d ranges from 0 for virgin material to 1 for completely damaged material with zero stiffness and it depends on state variable Y: The state variable Y reaches a maximal value during loading history between the damage threshold Y D0 and the equivalent strain eq : H H [4] where i H is the ith positive principal strain. Damage follows a damage evolution law which distinguishes tensile damage d t and compressive damage d c : where . t and . c are the weights computed from the strain tensor. In the present applications, damage due to tension will be considered only and . t =1 and . c =0. The damage model used hereafter is based on the following evolution law of damage: where A t and B t are the model parameters.
In this formulation damage is defined locally and strain softening is observed. This local formulation exhibits spurious strain localization (like any other strain softening local formulation). Consequently, numerical simulations yield a pathological mesh dependency and physically unrealistic results are obtained (Bazant and Planas, 1998). One possible remedy consists in reformulating the constitutive model in a non local approach, with damage at each material point depending on the strain not only at this point, but also in its neighbourhood. This non local model was developed by Pijaudier-Cabot and Bazant (1987), (see also (Bazant and Pijaudier-Cabot, 1988)). In this model, a non local equivalent strain is the weighted average of the local strains over the representative volume surrounding each point x in the material: ) (¸¸¹ l c , the internal length of the material, is related to the material heterogeneity. Model parameters, including the internal length may be obtained from optimal fits of experimental size effect test data (Le Bellégo et al., 2003a).

Chemical damage
The long-term behaviour of concrete structures subjected to environmental damage is a topic of the utmost importance for waste containment disposals. In the case of long-life wastes (such as radioactive wastes), the models cannot rely on experimental data over the entire service life of the disposal, which spans over several hundred years. Assessing such structures and comparing several possible designs prior to their construction must rely on robust and accurate computational tools, in which the knowledge of the various degradation mechanisms at the material level must be incorporated. Here, the mechanism of degradation considered is leaching due to water, and particularly calcium leaching as calcium is the most important element in the hydrated cement paste. The consequences of calcium leaching are an increase of porosity, a decrease of the elastic stiffness and a decrease of strength (Carde, 1996;Gérard, 1996;Carde and Francois, 1997). It is this decrease of the mechanical properties that is addressed more specifically.
As leaching is a very slow process, accelerated experimental procedures are needed in order to calibrate models, and to validate them subsequently. Goncalves and Rodrigues (1991) and Chen (1998, 1999) among others, have studied the influence of calcium leaching on the residual strength of beams for different cements subjected to different aggressive solutions. Carde (1996) used ammonium nitrate as an aggressive solution for accelerated leaching. Le Bellégo et al. (2000) have carried out size effect tests on leached specimens for different degradation rates, and Ulm and co-workers (2001) used a similar method in order to accelerate leaching.
At the same time, various models have been developed to describe calcium leaching, its effect on the mechanical properties of the leached material, and mechanical damage as well. co-workers (1999, 2001) used a chemo-plastic approach based on the mechanics of porous materials. Kuhl et al. (2000) devised a continuum damage model on the same basis. Gérard (1996) proposed a non-local damage model coupled with leaching effects. It is based on the definition of two damage variables: one is due to mechanical loads and the other one, which could be called a chemical ageing variable, accounts for chemical damage. These models can provide predictions of the residual strength and of the lifetime of concrete structures subjected to aggressive water but their range of validity in structural computations needs to be explored still.
The constitutive relations due to Gérard (1996) will be considered as an illustration. Among the models mentioned above, it is the only one that contains an internal length and is therefore capable of handling strain softening properly (Gérard et al., 1998).

Constitutive relations
Concrete is a reactive porous material. Hydrates and water are the main components. Solid phases are in thermodynamic equilibrium with the surrounding pore solution chemistry and calcium is the main chemical component of the binder. The hydrates leaching process is caused by the difference of composition and chemical activity between water in contact with concrete and the pore solution inside concrete. The difference between the pore solution and the external solution causes the motion of ions out of the cement paste. This overall motion of ions breaks the thermodynamic equilibrium established between the pore solution and the hydration products, and favours the dissolution of solid phases. The present simplified approach avoids considering all the elementary chemical mechanisms and focuses on the evolution of one variable only, the amount of calcium ions in the liquid phase of the hardened cement paste (Gérard et al., 1999). It is obtained by solving a non-linear diffusion equation: where solid C denotes the amount of calcium in the solid phase and C denotes the amount of calcium ions in the liquid phase (not to be confused with the elastic stiffness C in the previous section). A mathematical relationship for the calcium solid-liquid equilibrium can be found in Gérard (1996). It is a function of the cement chemistry and of the water to cement ratio. The diffusion coefficient D is not constant. It is related to the volume fraction of the different solid phases initially. As we will see further, it changes guring the leaching process.
The influence of chemical damage on the mechanical properties of the material may be introduced in the stress -strain relations following at least two methods. In the first one, which is inspired from the works of Ulm and Coussy (1996), the chemical degradation enters in the plastic yield function and, at the same time, changes the elastic stiffness of the solid skeleton of the material. A second method consists in adding a new chemical damage variable V on top of the mechanical one. This chemically induced damage variable is an independent contribution to the weakening of the material. It is a new internal variable.
This second method will be followed. It has been chosen by Gérard (1996), Gérard et al. (1998), Saetta et al. (1999, Kuhl et al. (2000), Meftah et al. (2000), and Stabler and Baker (2000) for environmental (chemical and thermal) damage modelling. The corresponding stress -strain relation reads: V is the chemical damage, a function of the concentration C of calcium ions in the pore solution: . g is an experimentally determined function shown in Figure 1.  (Gérard, 1996)) It may be also related to the porosity of the material in the case of a cement paste (Ulm et al., 1999). Figure 2 shows the material responses in tension and in compression at different levels of chemical damage. The mechanical evolution of damage affects the diffusion of the ions in the material too. In the absence of experimental data, the variation of the diffusion Strain Stress (MPa) coefficient due to damage can be neglected. In more advanced computations, it can be arbitrarily fixed as a S-shaped function ranging in between D 0 which is the value of the diffusion coefficient in the absence of chemical attack and D max which is the value for a cracked material which cannot exceed the diffusion coefficient in free water.

Example -residual strength of leached beams
In order to illustrate the coupled chemical and mechanical effect, we are going to compare the model predictions to experimental results. Experimental data on size effect tests have been taken from Le Bellégo et al. (2000) and Le Bellégo (2001).
The experimental program performed consists in residual strength tests: the residual behaviour of mortar beams was obtained after leaching in an aggressive solution (amonium nitrate). The specimens are mortar beams of rectangular cross section. In order to study size effect, geometrically similar specimens of various height D = 80, 160, and 320 mm, of length L = 4D, and of thickness b = 40 mm kept constant for all the specimens have been tested. The length-to-height ratio is L/D = 4 and the span-to-height ratio is l/D = 3. The beams were immersed into the aggressive solution during different periods of time (28, 56 and 98 days). Only the two lateral sides were exposed in order to have a unidirectional leaching front. The leaching depth f X was measured at the end of each time period. It helped at defining the leaching rate L R = b X f / 2 which is the ratio of the leached cross section to the initial cross section of the beam. For 28, 56 and 98 days of leaching, the leaching rates are 45%, 62%, and 84% respectively. After the chemical degradation tests, the residual mechanical behaviour of the beams were obtained from three point bending tests. A notch of depth D/10 and thickness 3 mm was sawed just before the mechanical test. The notch thickness was kept the same for all specimens. Its size is approximately that of the sand particles and the notch could be seen as a local defect in the beam. The response of bending beams entirely leached was also extrapolated from these data following the method described in Le Bellégo et al. (2000). The parameters controlling chemical damage have been fitted independently, including the diffusion coefficient in Equation [9] whose variation depends on the type of accelerated leaching used in the experiments (complete data are available in (Le Bellégo, 2001)). The parameters entering in the mechanical model are provided in Table 1. These parameters result from a manual (trial and error) fit of the size effect test data on sound beams (Le Bellégo, 2003a). It is now possible to compute the responses of beams at different leaching rates and to compare with the experiments. A one dimensional diffusion problem is solved first within the beam thickness in order to compute the calcium concentration. Then, the distribution of calcium ions is converted into a distribution of chemical damage according to the function plotted in Figure 1. It is averaged over the thickness of the beam. Then, the average chemical damage is inserted into the mechanical model (decrease of Young's modulus) and the load deflection curves of the beams are computed. In the present computations, the finite elements placed in the fracture process zone (constant strain triangles) have a size of 1/3 of the internal length approximately. Figure 3 shows the numerical predictions for leaching rates of 45% and 64%. Due to chemical damage, the average Young's moduli of the beams have been decreased by 25% and 40% in the mechanical computations respectively. These comparisons can be considered to be relatively satisfactory. The softening response for large beams is not very well described but the peak loads are predicted almost correctly. Some more accurate fit would request an adaptation of the constitutive relations with a variable internal length, as shown by Haidar et al. (2005) and Le Bellégo et al. (2003b).

Transfer properties and permeability
Transport properties of concrete, like permeability or diffusivity, are particularly important in case of structures for which tightness is important, for instance for prestressed concrete containment vessels in nuclear power plants. For such sensitive vessels, tightness to gas is critical during their service life, where concrete remains at most micro-cracked, but also upon minor accidents, when macro-cracks may appear locally. Hence, it is important to provide relationships between the amount of cracking and damage in concrete and its intrinsic permeability.

Figure 4. Evolution of the permeability as a function of the applied strain in a compression test on a hollow cylinder with a radial gas flow
Experimental test data (Choinska et al., 2007), performed on hollow concrete cylinders subjected to compressive loading and at the same time to gas flow through their thickness exhibit three regimes of growth of permeability. The first regime ranges almost up to the peak stress. It exhibits a slow increase of permeability due to the increase of the density of diffuse micro-cracks. Due to strain localization (transition between diffuse micro-cracking and macro-crack formation), a second regime of growth of permeability with a fast increase is observed experimentally, nearby the peak stress/peak strain. Finally, and for an applied strain which is larger than several times the peak strain, a third regime is reached. It is characterized by a slower rate of growth of permeability that is consistent with Poiseuille flow. Figure 4 shows these regimes as measured by Choinska et al. (2007) in compression.
Relationships between permeability and diffuse micro-cracking, described by a damage variable, have already been derived theoretically (Dormieux and Kondo, 2004;Chatzigeorgiou et al., 2005) and investigated experimentally (Picandet et al., 2001). This case corresponds to the first regime of the above figure. When macrocracks have completely formed, it seems consistent to assume that the apparent permeability of the specimen should be governed by the crack opening. We expect that this case corresponds to regime 3 on Figure 4.
Between these two extreme cases, there is a transition regime that needs to be properly captured in order to propose a law governing the growth of intrinsic permeability of a quasi-brittle material (concrete) as a function of damage. For that, we are going to define a consistent matching relationship between the two above extreme regimes. Further details on the following developments may be found in ).
We are going to focus on the intrinsic permeability of the material free of any stress and we will disregard the reversible stress effect on the permeability. Before describing how the permeability can be defined over the entire range of variation of damage, let us first discuss the two extreme cases of diffuse and localised damage.

Case of diffuse damage
In order to represent the interaction between diffuse damage and permeability at material level, the phenomenological relation established by Picandet et al., (2001) for damage lower than 0.15 is used. The intrinsic material permeability of the stress free material is an exponential function of damage: where k D and k 0 are respectively the current and the initial material permeability.
D and are the parameters fitted by Picandet (11.3 and 1.64 respectively) in the case of compression damage. It is assumed that it holds similarly in the case of tension damage too. As observed by Equation [11] with the reported values of the material parameters holds for a variety of concretes (ordinary, high performance, fiber reinforced) and it has been selected so that for standard concrete mixes, it can be implemented without any determination of material parameters, as a first approximation.

Localised damage. Crack opening versus permeability
At complete rupture, a macro-crack (or several macro-cracks) is expected. Hence, fluid flow will be mainly governed by these cracks and Poiseuille's law may be considered. For a fluid flowing between two parallel (rough) plates, the permeability is provided by the expression: [12] [u] is the crack opening and ] is the tortuosity factor related to the roughness of the crack (Boussa et al., 2001). Accordingly, and if we consider a crack of length L in a specimen of cross section S exposed to fluid flow, the total flux per unit thickness of specimen provides the apparent permeability ap k , by adding the fluid flow through the crack and the fluid flow through the undamaged material around the crack:

]
[13] Outside the crack, the permeability is assumed to remain equal to the initial one. It is a standard assumption in double porosity models. The "crack permeability" called k f is: [14] This permeability is the material permeability that should be reached when damage is close to 1, at material failure. It is the second asymptotic case, beside the case where damage is zero. The difficulty is that it is a geometry dependent parameter. A more convenient formulation of the contribution of the crack to the apparent permeability can be derived assuming that the crack is replaced by a band of intense damage.
Let us call c l O the width of this band. It is proportional to the internal length since the width of the fracture process zone in a non local damage model is indeed a linear function of the internal length. The apparent permeability derives now from the flux in the band with constant material permeability l k (damage is assumed to be homogeneous over the band for simplicity): We may now equate Equations [13] and [15] and the "crack" equivalent permeability l k is now: The local -point wise -permeability in Equation [16] does not depend on the geometry of the cracked solid. It may readily be used in a continuum model. If we consider a crack of length L in a specimen of cross section S exposed to fluid flow, Equation [13] will be recovered at the structural level (apparent average permeability). The remaining difficulty is the calculation of the crack opening. This point will be discussed in Sub-Section 4.4.

Matching law
The two expressions of the permeability, k D (Equation [11]) and k l (Equation [16]), may now be combined as they represent two asymptotic cases. For this purpose a simple matching law based on logarithm of the permeability is used: The derivation of a crack opening as a result from a continuum damage mechanics calculation is a major challenge. We discuss next a simplified approach to this aim.

Computation of the crack opening displacement
In order to compute a crack opening from a damage field, several possibilities exist. The first one consists in the determination of the crack that is equivalent to the damage zone. This equivalence may be performed following the technique due to Mazars and Pijaudier-Cabot (1996), based on the energy dissipation due to cracking. Once the crack equivalent to the damage zone is obtained, fracture mechanics provides the distribution of crack opening depending on the applied loads. A second technique is the extraction of the crack opening by comparison of non local strain profiles: the non local profile that corresponds to a crack and the profile that is computed from a non linear FE calculation. This technique is due to Dufour et al. (2008).
Another possibility is the derivation of the crack opening from the constitutive relations.
On the left hand-side, the stress corresponds to the effective stress in the material: Substitution this equation in the previous one provides the inelastic strain: Once this inelastic strain is computed, integration of the component normal to the crack over the finite element which contains this crack is performed. The result is the crack opening displacement. The crack opening displacement may also be included directly in the constitutive relations.  devised such a constitutive relation recently.
Here, we are going to derive the equations corresponding to the case where the displacement is obtained from the integration of the strain over the damage band. This strain will follow directly from the law governing the growth of damage instead of following from the computation of an inelastic strain as shown above. In the case of monotonic loading, and close to failure where the elastic strain is negligible, this calculation should be close to the above derivation, more restrictive because unloading and closure of cracks cannot be accounted for. a) b)

Figure 5. a) cracked domain (discrete case); b) damaged domain (continuous case)
We assume again that the crack, according to the damage model, is described as a band of intense damage of width c l O . An equivalence between two domains is then considered: a cracked one with a discontinuous crack opening [u] (see Figure 5a) and a damaged one with a damaged zone of width Ol c (see Figure 5b). This equivalence states that the displacement across the crack is described as the integral of the strain across the damage band. Subsequently, we assume that the loading is monotonically increasing. According to Mazars' model, the two quantities are similar and the crack opening is now expressed as: Assuming that the damage distribution in the damaged zone is uniform, one obtains the crack opening as: To substitute this crack opening with a damage field, we relate the state variable with the damage through the inverse of the damage evolution law (see Equation [2]): Finally, by substitution of Equation [24] in Equation [23], the crack opening may be represented as a function of the damage: We may now plug this expression into the equations defining the growth of permeability.

Examples of implementation
For the sake of illustration, let us consider the material parameters given in Table 2 and let us compute the evolution of material permeability as described by the matching law. The internal length l c is arbitrary chosen equal to 0.02 m, while O, which influences the width of a damaged band, is chosen equal to 2 (Bazant and Pijaudier-Cabot, 1988). The initial permeability 0 k considered in simulations is taken equal to Evolutions of permeability with damage are shown in Figure 6. The matching law provides a representation of the permeability for small damage (close to Picandet's model), as well as for intense damage where it tends towards the permeability given by Poiseuille's law.
We consider now experiments on discs subjected to uniaxial tension according to the "Bipede" configuration (Gérard, 1996). A disc made of mortar is glued to steel plates (Figure 7). Tension is applied to the plates and it is transmitted to the disc in which one or two cracks appear. At the same time, water flows through holes in the steel plates and through the disc. The apparent permeability (to water) is measured in the course of loading and cracking of the disc. Bipede test geometry (after (Gérard, 1996) For the purpose of comparing with experimental data, a local permeability is computed first at each stage of damage in the mortar disc. This computation is performed according to the matching law. The expression of the crack opening displacement is taken from Equation [25] under the assumption that strain and damage are homogeneous within the damage band. Then, a structural permeability is determined by a standard averaging of the local permeability.
We consider first the results reported by Jason et al. (2007). In Figure 8, the comparison between the predictions of average permeability computed by considering that the permeability follows the diffuse damage equation (Picandet's relationship) locally and the test data is presented. This comparison is very weak, which justifies the need for a permeability law that matches properly the extreme cases of diffuse and localised damage.
Eq. (14) x Eq. (15) K/K 0 Damage Figure 8. Evolution of the permeability with damage. Comparison of Picandet's relationship (denoted as Equation [14]) and a power law for damage (denoted as Equation [15]) with experimental data The second one cannot be obtained as a consequence of the first one. It is not possible to capture the evolution of the apparent permeability with the diffuse damage case alone. This relationship overestimates the variation of permeability locally and it is not because such an overestimation occurs in a very small region of the structure (damage localisation band) that it can be neglected. Figure 9 shows the comparison between the same test data and the present matching law. In the comparison, we have entered in our formula the same physical constants as those in Gérard (1996): the parameters of the damage model, the crack roughness and the dynamic viscosity of water. We have kept the same coordinate system which relates the applied strain to the water permeability of the cracked specimens. The comparison with experimental data is now quite consistent. Figure 9. Evolution of the permeability to water with the strain in the "Bipede" test

Conclusions
Continuum damage mechanics offers a consistent framework for the description of cracking and failure of concrete and quasi-brittle geomaterials. We have illustrated our presentation with a simple isotropic damage model, but many refined damage models exist in the literature (see e.g. (Carol et al., 2001;Desmorat et al., 2006)).
Damage models may be complemented with chemo-mechanical coupled effects and the relationship with permeability has been also discussed. In the case of calcium leaching, an additional damage variable describing the mechanical effect of this chemical phenomenon is added. Solution strategies involve the resolution of a mechanical non linear problem coupled with a diffusion problem that is also non linear. Some satisfactory comparisons with accelerated tests can be achieved, including on size effect tests. Test data show that there is also a correlation between the internal length in non local damage models and the evolution of microstructure due to ageing which was not discussed here but should be accounted for in structural failure analyses (Le Bellégo et al., 2003b).
A formula describing the evolution of permeability with damage has been discussed. It matches consistently two extreme cases: in the first one, permeability is a function of distributed damage, while in the second one it is governed by crack opening. Several techniques for computing crack opening displacements from continuum damage modelling have been briefly recalled. We have used one of these, where the crack opening is related to the variable that governs damage in a non local integral damage model. Afterwards we have related this state variable with damage in order to arrive to an expression where the permeability is controlled by the variation of material damage only. Preliminary comparisons with experimental data provide qualitative and quantitative results that are quite consistent.