T H Èse

Cristian Julián

Nieto Gamboa

Grade De Docteur

M B Cambou

Président M C Mariotti

Cea-Dam Rapporteur

M A Daouadji

Université Paul

Verlaine Rapporteur

J-J Fry

Edf-Cih Examinateur

Mme A Modaressi

M F Lopez-Caballero

MECHANICAL BEHAVIOR OF ROCKFILL MATERIALS APPLICATION TO CONCRETE FACE ROCKFILL DAMS Soutenue le 28 Mars 2011 devant la commission d'examen composée de

ECAP0018

Après trois ans de thèse partagés entre le

Dedicada a Maria Alejandra, a mis padres y a mi hermano. i ii

Context

The concrete face rockfill dams (CFRD) are structures frequently chosen in a hydropower project due to their different advantageous characteristics. During the last century, the increase of electricity demand has led to construct continuously higher and higher CFRDs for hydroelectrical projects. Moreover, the imposition of cost reductions has led to reduce the volume of the embankments by imposing stiffer slopes up and downstream, affecting in this way the factor of safety. However, no significant efforts have been made to enhance the current state of practice, which remains totally empirical. The limits of such approach have been shown by recent damages, specially evidenced during first impoundment in high CFRDs. These damages showed how the extrapolation of current practice to larger projects requires a deeper knowledge, and the identification of the reasons that triggered the problems recently experienced. Two main issues have been identified as important to make a step forward in the design of CFRD, i) the understanding of the soil-structure interaction between the embankment and the concrete face slabs, and ii) the correct description of the mechanical behavior of rockfill materials. This thesis focus on the later issue, the mechanical behavior of rockfills, in this sense we look to contribute to its understanding. Traditionally the knowledge in soils mechanics have advanced through well controlled experiences in laboratory. However, for rockfill materials a restriction is imposed due to the particles size reaching 1m in diameter. Therefore, to elaborate a "representative elementary volume (REV)" it would lead to samples of 8-10m diameter and 16-20m height (following the traditional rules of triaxial tests). The dimensions of current soil laboratory facilities are not adapted to test such a representative elementary volume of material, and the construction of new facilities would be not only highly expensive but also technically impossible. Therefore, for practical purposes, it is required to define reliable techniques to estimate the mechanical parameters of rockfill materials.

As testing coarse granular materials as rockfills is very difficult, the few laboratory experiences in rockfill-like materials are priceless, and an effort has been made to recover different laboratory experiences characterizing the mechanical behavior of coarse granular materials. This compilation will show different mechanical characteristics, and especially the influence of particles breakage on the mechanical behavior of crushable soils.

Different approaches are found to model soils behavior including particles breakage. One trend regroups those based on the theory of plasticity, where the influence of particles breakage is represented as an increase of plastic strains linked to a form of total work. This type of approach does not take into account any mechanical characteristic concerning the particles and it is not able to give information about the amount of crushed particles. Other trend looks to relate the variation of solid surface due to breakage of particles with the mechanical behavior of materials. This relation is made soundly under a thermodynamic framework. This approach is explored here to understand the relation among elastic, friction and breakage phenomena in granular materials.

1 This thesis has been developed in the framework of the CIFRE2 program sponsored by the French national agency for the technical research (ANRT ) and Tractebel Engineering France in collaboration with Ecole Centrale Paris.

Introduction

A dam is an artificial obstacle to a naturally flow of water constructed for one or more of the following purposes: storing water for irrigation, raising the water level to provide a head which can be turbined and transformed in electrical power (hydro-electricity), supplying water to town or industry, prevention of flooding, forming an artificial lake for navigation, leisure, activities, etc. From a mechanical point of view the dam is a structure and should be able to support the hydrostatic pressure and transfers it to the foundation, at the same time the dam should be impervious enough to avoid important leakage of water that could either reduce the reservoir water level or compromise the safety of the structure.

In seismic areas, it should also resist to the earthquake loading.

There exist different kinds of dams that could be regrouped into two general families. The first family is composed of rigid dams constructed in general of masonry or concrete. In this family we find the gravity dam, the arch dam, the buttress dam and the variations of them. For these dams the impermeability condition is assured by the nature of materials used in construction. The second family is composed of flexible dams constructed either in earth or rockfill. In contrast to the first family the materials used here do not assure the low permeability needed to avoid leakage. Therefore an impervious structure should be added. This impervious structure can be a clay core, a bituminous core, a concrete or bituminous face, or a geomembrane.

The selection of one or another type of dam for a specific project is made in accordance to the most feasible option. It is known that the rigid dams require very sound foundations while the flexible dams adapt easier to foundations with lower characteristics. The economical feasibility is a key point during the dam's selection stage in a project. It is looked forward to reduce the costs associated to transport of materials due to the fact that dams are massif structures of important volume. Therefore there is a financial interest to construct the dam with materials that could be found not very far away from the site. In resume some of the factors affecting the selection of a type of dam for a particular project are: geological, technical, economical and environmental.

Nowadays one of the most popular kind of dams is the concrete face rockfill dam (CFRD). In this kind of dam the impervious organ is a concrete slab cast in place over the upstream face. The body is composed of granular materials of sandy type for filters and higher sizes for the rest of the structure. The height of CFRDs constructed during the last century has experienced a constant increase as shown in Figure 1.1. However, despite this evolution the design techniques have remained completely empirical.

Figure 1.1: Evolution of the concrete face rockfill dams' height during the XXth century, after ICOLD (2002).

The increasing tendency of CFRDs is explained by several reasons, among them: the relative commodity for construction, the low cost of materials, its overall stability and the short time of construction. Until not a long time ago, it was thought that there was no need to make any refined calculus (e.g., finite elements) for the CFRD. These features made of the CFRD very popular and encourage to its selection, but the requirements in size or height have increased over the years leading to some drawbacks.

The first part of this chapter aims to describe the components of a CFRD, in order to establish a framework to analyze the problems encountered in recent projects around the world. Afterwards, the problems analysis points two main topics of interest: i) the mechanical behavior of rockfill materials at high stresses, and ii) the complex interaction between the concrete face and the rockfill embankment. Unfortunately, only the mechanical behavior of rockfill materials will be treated here and it will be the main topic of interest during the rest of this thesis.

CRFD Components

A CRFD as a whole is composed of different sub-structures with different mechanical and physical characteristics. Figure 1.2 shows a schema of a typical CFRD section with its components, among them the plinth, the face slabs, the parapet wall, the grout curtain, the filters zone and the rockfill zone. Nowadays the general trend establishes upstream and downstream slopes with a relation of 1.4/1.0 (horizontal/vertical), however it depends on the rockfill and foundation characteristics. For sound compacted rockfill founded on a sound rock foundation outer slopes can be as steep as 1.3/1.0 (ICOLD, 2002). The characteristics of zones shown in Figure 1.2 are:

Zones 1A : upstream concrete face protection zone, composed of cohesionless silt or fine sand with isolated gravel and cobble sized rock particles up to 150mm. It can act as a joint or crack healer over the perimeter joint and the lower part of the face slab. This zone is placed in 200 to 300mm height Zone 1B : upstream concrete face protection zones; this zone provides support for Zone 1A and it is composed of a random mix of silts, clays, sands gravels and cobbles. As previously this zone is placed in 200 to 300mm height layers. Compaction of this zone is only performed by spreading equipment. It presents similar characteristics to zone 1A but with a higher maximum particle size.

Zone 2A : Is a processed fine filter (sand and gravel) with specific gradations to limit leakage in the event of waterstop failure. It spreads up to 3m from the perimeter joint. This zone prevents the movement of silt size particles through it and, thus, serves as secondary defense against leakage. Placement is done in layers of 200 to 400mm height. Compaction is required with vibratory compactors. The gradation requirements of this zone can be found in ICOLD (2002). This zone should be carefully designed to fulfill the following requirements: i) it will not segregate, ii) will not change in gradation, iii) will be cohesionless, iv) will be internally stable, v) will have sufficient discharge capacity and vi) will have the capacity to control and seal a concentrated leak.

Zone 2B : Provides support to the face slab and consists of sand and gravel-sized particles, placed in 400mm horizontal layers and normally compacted with 4 passes of a 10ton smooth-drum vibratory roller. The horizontal width of the zone varies from 2 to 4m depending on the height of the dam. This zone consists of material equal or nearly equal in quality to concrete aggregate. The gradation requirements of this zone can be found in ICOLD (2002). A special case nowadays is the construction of a concrete curb at the upstream face after every layer. This new structure allows to confine the zone in a more effective manner.

Zone 3A This zone is a transition between the zone 2B and the zone 3C. It consists of rockfill materials with a maximum particles size up to 400mm placed in layers of 400mm height and compacted with at least 4 passes of a 10ton or heavier smooth-drum vibratory roller. The horizontal width of the zone varies from 2 to 4m depending on the height of the dam. The gradation criteria are near those of the filters.

Zone 3B : This zone consists of rockfill material with maximum particle size up to 1000mm placed in 1000mm layers and normally compacted with 4 passes of a 10ton smooth-drum vibratory roller. The final decision concerning the appropriate number of passes can be estimated after conducting tests in which the average surface level of the rockfill layer is determined by surveying at intervals of 2 to 12 passes.

Zone 3C : This zone consists of rockfill material with maximum particle size up to 2000mm placed in 2000mm layers and normally compacted with 4 passes of a 10ton smooth-drum vibratory roller. As for the previous zone the compaction guidelines can be obtained from experimental results in place.

Zone 3D : This zone provides drainage conditions within the embankment. In well-drained rockfill embankments this zone is placed at the base of the dam within the valley section. In poorly-drained rockfills or in fills consisting of sands and gravels these drainage zones may consist of a continuous chimney or wall drain.

The above features are described in ICOLD (2002). However, the incidents experienced recently in CFRD have led the profession to move towards a more stringent specifications regarding the maximum particle size and layer thickness, which should not exceed 800 mm [START_REF] Frossard | On the structural safety of large rockfill dams[END_REF][START_REF] Frossard | Du sol à l'ouvrage: une vision multi-échelles de la géomécanique[END_REF]. The CFRDs are composed also of other structures as:

Concrete Parapet Wall : The main purpose of the parapet wall is to reduce the total volume of rockfill. The height of the parapet wall is function of the amount of rockfill saved, therefore to determine the height of the wall an economical analysis should be done.

Grout Curtain : It is performed under the toe slab (plinth) and outside of the embankment body. The main role is to reduce the high hydraulic gradient developed at the dam's toe by increasing the path of the flow. The technical specifications of the grout curtain depend on the geological conditions of the site. However, a depth within the range of 1 2 to 2 3 of the reservoir head at the location of the curtain is often specified. Details about water-cement ratios and injections pressure can be found in ICOLD (2002).

Toe Slab -Plinth : This structure connects the foundation with the face slabs, a schema is shown in Figure 1.3. The technical specifications depend on the geological characteristics of the foundation and the hydraulic gradient expected. The plinth is connected to the face slabs by the perimeter joints, completing then the upstream water barrier of the dam.

Perimeter joints : The main function of the perimeter joint is to maintain a watertight seal against full reservoir load while allowing for anticipated movements between the toe slab and face slabs. In order to fulfill this mission the perimeter joint includes deformable materials as shown schematically in Figure 1.3.

Face slabs : The concrete face is composed of concrete slabs and conforms the main watertight component of the dam. The face slabs are constructed over the concrete curb (if it exists) or the zone 2B. The design of the slabs comprehends its thickness and location of vertical and horizontal joints. Selection of face slab thickness is typically based on previous experience, height of the dam, and minimum conditions for reinforced concrete. However, current guidelines suggest a slabs' thickness increasing with water head for CFRDs over about 100m. The variation of the thickness according the current guidelines is given by:

T = 0.3 + x • H (1.1)
where T is the thickness of the concrete face slab (in meters) for a point located at a depth H from the maximum water level, and x is a parameters such that 0.002 ≤ x ≤ 0.004. The slab width ranges from 12 to 18m, with panels width of 15m being common. However special abutment conditions can lead to special designs as for instance at Mazar project in Ecuador (Frossard, 2006a).

Recent Problems in High CRFDs

Different CFRD projects around the world have experienced important problems during impoundment stage. The consequences are the cracking of the concrete face slabs and joints opening. Therefore the waterproof efficiency is reduced and the increase in leakage flow is observed. Several problems in high concrete face rockfill dams have been experienced during the last five years, leading to wonder about the limitations of the current state of practice.

Mohale CFRD located in Lesotho is one of these cases. Design description of the project is found in [START_REF] Gratwick | Mohale dam, lesotho[END_REF]. It concerns a 145m high CFRD constructed in accordance to the current guidelines. The first reservoir impounding was completed by March 2006 after an event of heavy rains. After full impoundment, concrete face slab cracking was observed (Figure 1.4). The slab cracking resulted in a significant jolt, which likely broke the bond between the concrete slab and the underlying extruded curb concrete. The final result was the increase of joint openings, both for vertical as well as for perimeter joints at flank sections [START_REF] Johannesson | Lessons learned from mohale[END_REF]. Subsequently seepage increased significantly, peaking at some 600l/s. A second recent case concerns Barra Grande CFRD in Brazil [START_REF] Pritchard | Taking the empirical approach[END_REF]. The dam of 185m high and crest elevation at 648m.a.s.l1 started the reservoir impounding in July 2005, and by September it presented an increase in leakage flow. The explanation to this fact remains on the face slab rupture discovered at vertical joint located near the midsection of the dam. The leakage flow at 92% of reservoir capacity reached 428l/s. Also in Brazil, the CFRD of Campos Novos, with a height of 202m and a crest length of 590m, suffered important damages during impoundment in 2005. When reservoir water level reached about 92% of the normal storage level, face slab rupture was discovered at the upper part of the vertical joint between slab 16 and 17, near the middle section. After some remedial works a problem in a derivation tunnel leaded to empty the reservoir. Afterward some other important cracks were revealed over the concrete slabs as shown in Figure 1.5. Figure 1.5: Failure of concrete slabs in Campos Novos rockfill dam, after [START_REF] Pritchard | Taking the empirical approach[END_REF].

Similar problems have been reported for Tianshengqiao No 1 project in China [START_REF] Xu | Research on concrete face slab rupture of high cfrd[END_REF], where similar crack patterns have been observed for the 178m high dam. However not all large CFRDs have suffered structural damages, there exist also successfully operating cases as Karanjukar project in Iceland.

Problem Analysis for High CFRDs

As shown before, several problems in high CFRDs have been observed. They can be classified chronologically as: i) Construction, ii) Impoundment, and iii) post impoundment problems. Some explanations to these problems are described hereafter.

Construction Stage

During construction, the concrete face slabs is casted in place at the same time as the embankment construction. During this stage there are no forces at the outer face of the slab. Once the concrete slab reaches some stiffness it does not follow anymore the deformation of the rockfill shell. Therefore, the construction deformations of the rockfill lead to an eventual opening between the slabs and the upstream face of the dam [START_REF] Shi | Characteristic and experience of the design, construction and performance of tsq-1 concrete faced rockfill dam[END_REF].

Impoundment Stage

At the filling stage, the stresses caused by the water pressure are transmitted to the rockfill body by the concrete slabs. The deformation is showed schematically in Figure 1.6 and a 3D view in Figure 1.7. The main issue during this stage concerns the cracks formation in the concrete face. According to the observations of the cases mentioned before three patterns of cracks can be identified. i) Vertical compression cracks at the upper middle section of the concrete face. Compression can lead eventually to jolt the middle section slabs. These cracks are induced by high compressive stresses. ii) Near the lower third part of the concrete face, horizontal cracks can appear as in the case of Campos Novos. Some other diagonal cracks have been observed near the abutments, more or less oriented to the same inclination of the abutments. Another important issue concerns the perimeter joint opening, schematically shown in Figure 1.6c. This opening is induced by the tensile stress concentrated at the rockfill-abutment interface.

An explanation for the vertical cracking phenomena can be as follows: In a first time the pressure of water-load acts over the concrete slabs then it is transferred to the rockfill, afterwards and as consequence of the load over the rockfill a settlement appears (i.e., vertical deformation). Subsequently the vertical movement produces a horizontal movement because of the presence of abutments. Therefore, and because of the contact between the concrete slabs and the rockfill, the horizontal movement in the vicinity of the abutments produces shear stresses at the rockfill-slabs interface that drag the slabs horizontally. Finally because of the same dragging effect acting at both abutments but in inverse sense the center of the concrete slabs is constrained to high compressive stresses. Concerning the horizontal tensile cracks they could be produced by the flexion of the slabs as shown in Figure 1.6a. In any case once the concrete face cracks, it looses the waterproof protection.

Post-Impoundment Stage

The effect of time is another important point and even today is a research topic of interest. The level of deformation reached in the long term behavior of the structure can be more or less important. Figure 1.8 shows the behavior of rockfill dams through time, and as mentioned before the reasons of this behavior are not yet well clarified. However from the mechanical point of view it can be treated as a viscous behavior.

Regarding Figure 1.8 the importance of compaction on the reduction of crest settlement can be noticed. Another point concerns the influence of water load on the deformation. By comparing settlement curves for El Infiernillo dam with a height of 146m and Alto Anchicaya dam of 140m, it is observed that the settlement at the crest level is higher in the former. One of the reasons of this difference could be the water saturation degree in the upstream part of the dam. For El Infiernillo the waterproof structure is a clay core at the middle of the cross section and consequently the upstream cell of the dam is fully saturated, while in contrast, for Alto Anchicaya the upstream cell is relatively dry.

A particular point concerns the influence of the abutments' geometry on the mechanical behavior of the dam. On this topic [START_REF] Zhang | Effects of abutments slopes on the performance of high rockfill dams[END_REF] searched the influence of abutments on rockfill dams with central clay core. Their research was performed using the finite elements method with a non linear elastic constitutive model. Their works concluded the following points : Figure 1.7: 3D schema of deformation of a CFRD due to impounding, after [START_REF] Harita | Stress and deformation of concrete slab in cfrd[END_REF] 1. The condition of a narrow valley restraints the movements of the earthdam and reduces the vertical stresses in the lower parts of the dam, which is called as the arching effect.

2. The slope of the banks has a significant consequence on the arching effect across the valley. The relationship between the vertical stress and the stress by self weight is a function of the height of the dam, the slope of the banks, the parameters of resistance of the interface between the rockfill and the bedrock and building materials. Significant changes in slopes of the banks will lead to a significant level of tension stresses at the contacts of the core and therefore must be avoided.

3. About the arching effect in the longitudinal direction the authors conclude that when the bank slope is steeper there is more probability of hydraulic fracturing.

Figure 1.8: Long-term settlement records in crest for different types of rockfill dams, after [START_REF] Oldecop | Theoretical investigation of the time-dependent behaviour of rockfill[END_REF].

Concerning the main problems of the CFRDs the main issue pointed out is the level of deformations. The rockfill materials present a lower stiffness modulus than the concrete, but they are both solidary in deformations. Therefore the concrete will reach its yield limits earlier, as shown by the different cases presented before. Theoretically the solution concerns the increase of stiffness modulus of the rockfill to obtain a reduction of the deformations. However, it is not as easy to do as to say, because there are different parameters affecting the mechanical behavior of rockfills. In this sense a review of the main characteristics of the mechanical behavior of rockfills is made hereafter.

Mechanical behavior of Rockfills

From an experimental point of view [START_REF] Leps | Review of shearing strength of rockfill[END_REF] presents a compilation of rockfill mechanical characteristics, mainly concerning maximum friction angle, confining pressure, particles size and sample density. Some other triaxial tests have been done and reported since Leps's compilation. Therefore to update this database, a consolidation of some tests reported initially by [START_REF] Leps | Review of shearing strength of rockfill[END_REF] and other experimental results from drained triaxial tests published in scientific papers or in technical reports is presented in Appendix A. It includes some additional information as uniformity coefficient, specific weight of solid particles and dilatancy characteristics at maximum stress ratio (σ 1 /σ 3 ).

Two main issues concern the mechanical behavior of rockfills, the ultimate shear strength which defines the stability of the structure, and the deformability characteristics which determine the admissible state of service. Along the time, different experimental studies have shown that mechanical behavior of rockfill materials is characterized by : i) a non-linear, inelastic and stress dependent stress-strain behavior, ii) an increase in shear strength, axial strain and volumetric strain at failure with confining stress, and iii) an increase in volumetric strain with the particles size for the same confining pressure. These three main characteristics are affected in a different degree by different factors, including : the initial stress state, the stress path, the initial compacity condition, the mineral origin of rockfill particles, the size and shape of particles, the water content condition and the grain size distribution. Hereafter some of the works revealing the influence of these factors are reviewed.

Stress State

It is well known that due to the non linear behavior of geo-materials, the response to a given mechanical solicitation is different depending on the initial stress state. For granular soils the effects are observed on the maximum deviatoric stress, the dilatancy, the deformation modulus and the friction angle. These conclusions have been drawn from experimental works that have been performed by several researchers looking at the mechanical response of sands. Nevertheless, on the rockfill side the experiences are less frequent, hopefully we can find some remarkable works by [START_REF] Marsal | Research on the behavior of granular materials and rockfill samples[END_REF], [START_REF] Fumagalli | Test on cohesionless materials for rockfill dams[END_REF], [START_REF] Marachi | Strenght and deformation characteristics of rocfill materials[END_REF], Nobari and Duncan (1972) [START_REF] Charles | The influence of confining pressure on the shear strength of compacted rockfill[END_REF], [START_REF] Indraratna | Large-scale triaxial testing of grey wacke rockfill[END_REF], [START_REF] Chávez | A constitutive model for crushed granular aggregates which incluyes suction effects[END_REF] [START_REF] Varadarajan | Testing and modelling two rockfill materials[END_REF], Varadarajan et al. (2006a), Kohgo et al. (2007a) and a few others.

From some of these works (compiled in appendix A) we can plot the peak friction angle (φ peak ) against the confining stress (σ c ). The result shown in Figure 1.9 suggests the reduction of the peak friction angle with increase of confining stress. 2004) [START_REF] Marsal | Large scale testing of rockfill materials[END_REF] Figure 1.9: Variation of peak friction angle with confining stress from data in Table A.1. The hatched area covers more than 90% of the data compiled.

Ultimate strength state is normally represented by the Mohr-Coulomb criteria, which is normally assumed as a linear function in the plane of normal and shear stresses. The slope of this linear function is related to a constant friction angle. Nevertheless, as shown in Figure 1.9 the friction angle can variate with confining pressure. This fact has led to formulate curved shear strength envelopes. For instance [START_REF] Charles | The influence of confining pressure on the shear strength of compacted rockfill[END_REF] propose a power criteria of the type:

τ = A • σ b ∴ b < 1 (1.2)
where A and b are rockfill constants. The curvature of the shear strength envelope seems to be related to the breakage of particles. Let us adopt the proposition of [START_REF] Terzagui | Soil Mechanics in Engineering Practice[END_REF], where the friction angle is the sum of two independent components. The first component is the inter-particles friction (φ µ ) given by the friction between two flat surfaces and attributed mainly to mineral composition of the grains surface. The second component is given by a geometric interference (φ g ) component, an elementary analysis showing this relation is presented in appendix A.1. Then the friction angle, associated to the slope of the failure envelope is represented by:

φ = φ µ + φ g (1.3)
When confining pressure increases the overlap and displacement of grains outside the shear plane becomes more difficult. The movement between particles is completed by crushing of flaws or fracture of the entire particles, resulting in a decrease of the angle of geometric interference. This effect has as consequence the reduction of the total friction angle. As breakage increases with the stress state, then it can explain the curvature of the Mohr-Coulomb envelope.

Influence of particles origin

The particles origin concerns not only the mineralogic characteristics but also the extraction conditions, e.g., materials quarried from a mine, extracted by explosives, or extracted by backhoe, such as alluvial materials. The origin of a given material influences in a double manner. First, the conditions of the particles surface affect the friction angle between particles, and second, the elastic and strength properties are dependent of the mineral characteristics. The mechanical characteristics of particles are then associated to the mechanical behavior of an equivalent continuum as a rockfill sample [START_REF] Biarez | Influence de la granulometrie et de son evolution par rupture de grains sur le comportement mecanique de materiaux granulaires[END_REF].

An example of the influence of mineral origin is found in [START_REF] Charles | The influence of confining pressure on the shear strength of compacted rockfill[END_REF]. The summary of its properties are shown in Table 1.1. Table 1.1: Rockfill properties tested by [START_REF] Charles | The influence of confining pressure on the shear strength of compacted rockfill[END_REF]. PS: Specific gravity of solids, Ws moisture content of particles. d1, d2, d3:particle's diameters following three orthogonal directions. I1, I2: strength indexes of particles cf., reference. ρ: sample's density, w sample's water content, n sample's porosity, Sr sample's saturation degree

Shape of Resistance Indexes

The results of triaxial tests carried out on the materials described in Table 1.1 are shown in Figure

Figure 1.10: Drained triaxial compression tests on heavily compacted rockfill with σ 3 = 100kP a, after [START_REF] Charles | The influence of confining pressure on the shear strength of compacted rockfill[END_REF]. Under same test conditions, samples with different origins show different mechanical characteristics 1.10. It is seen that the response of different materials with different origins is different even if the geometric characteristics, grain size distribution and test conditions remain the same.

Influence of grains shape

The influence of grains' shape is an important factor for compacity, frictional and crushing characteristics. For frictional characteristics, the effect of grain's shape has been studied by [START_REF] Frossard | Effect of sand grain shape on interparticles friction: indirect measurement by Rowe's stress dilatancy theory[END_REF] who concluded that interparticle friction angle seems to increase with angularity and that sphericity of particles has a notable effect on volumetric strains. For compacity (i.e., the volume fraction that is filled) because during compaction the longer angular particles are more difficult to put in a dense arrangement than rounded regular particles. During loading stages the flaws act as stress concentrators and will break easier than the rest of the particle. This will represent a higher deformability of the material. Because of the presence of flaws in angular particles and their breakage during loading, the angular materials could generate more finer particles than the non angular materials. Some experimental results reported by [START_REF] Varadarajan | Testing and modelling two rockfill materials[END_REF] with alluvial (rounded) and quarried (angular) materials show strong differences between the behavior of both materials concerning shear resistance and deformation response. From discrete element analysis Nouguier-Lehon and Frossard (2005) showed that the angular particles dissipate less energy by rolling than rounded particles. This is explained by the restriction of rotations due to face to face contacts. Moreover [START_REF] Guises | Granular packing: numerical simulations and the characterisation of the effect of particle shape[END_REF] show how the variation of the particles' shape influences the distribution of contact forces inside a sample.

Influence of water content

The water content, as the particles origin, have a double consequence. Initially the increase of water content will affect the surface of particles where the condensation of water occurs. By affecting the surface of particles it directly affects the inter-particle friction angle (φ µ in equation 1.3) and in consequence a reduction in the shear strength can be obtained.

The second consequence is related to the condensation of water in the particles micro-cracks, which influences the particles strength and facilitates the crushing. This topic will be exposed in more detail in Chapter 2.

Influence of grain-size distribution

The grain size distribution influences the compacity of a material. Uniform materials at the densest state, will be less dense than materials composed of different particles sizes. This is so because smaller particles can fill the interstitial spaces. It is known that for rockfill materials the denser it is, the stiffer it behaves. Then the grain size distribution plays a role on the deformability. In civil engineering a grain size distribution called Fuller distribution is often used. Its origin comes from the experimental research of aggregates for concrete. This theory proposed by Fuller and Thomsom, expresses that the grain size distribution which represents the densest material state is given by equation 1.4,

F = d D n (1.4)
where n is a power less than 1. In practice a value of 0.45 have been proposed for this parameter, D is the maximal diameter and d represents any particle diameter smaller than D. The grain size distribution plays another important role in particles crushing. The increase in compacity influences the distribution of forces between the particles, by transmitting forces to smaller particles that, as will be seen after, have a lower breakage probability.

Influence of grain breakage

Factors like particle's size, shape and gradation are hardly linked analytically to the mechanical behavior of materials. Therefore, their association is made empirically, as the work of [START_REF] Barton | Shear strength of rockfill[END_REF] who proposed as expression linking particle's size, shape and resistance to an equivalent friction angle of the media. The association made by Barton and Kjaernsli goes directly from particle's size and shape to the mechanical behavior in terms of friction angle.

The approach taken here is relatively different. It is believed that factors such as particle's size, shape and gradation affect directly the breakage of particles, and then, breakage of particles affects the mechanical behavior of the rockfill materials.

A detailed analysis of the influence of grain breakage on mechanical behavior of granular materials is the main subject of the next chapter. In the meanwhile Table 1.2 presents a schematic correlation of the influence of different factors on breakage.

Conclusions

Recently problems in concrete face rockfill dams have shown the limits of the currently empirical practice. The principal problem concerns the development of cracks at the concrete face, reducing in this way the impermeability condition of the structure. These cracks appear by the excessive deformation at the upstream face of the dam. This deformation is due to the compressibility of the rockfill material, i.e., to its mechanical behavior.

The compilation of experimental data has shown that shear resistance in terms of peak friction angle varies for different conditions. The range of peak friction angle varies from 32 o to 60 o . This is a wide range for rockfill materials showing that shear resistance is affected by different factors. Among them, the confining pressure, nature of particles, shape of particles, water content, grain size

Particles strength ⇓ Friction Angle a ⇓ Voids Ratio ⇑ b ⇑ ⇓ ⇓ ⇓ ⇓ Gradation ⇑ ⇑ Breakage ⇓ ⇓ ⇑ ⇑ ⇓ ⇑
a The friction angle between particles b A higher friction angle makes more difficult to compact the material, or requires more energy to reach a given void ratio. c A measure for the particles' shape is the ratio between two orthogonal directions

d i d j , with di ≥ dj
distribution and notably breakage of particles. This last factor, is studied in more detail during the next chapters, in order to identify its influence over the mechanical parameters of rockfill materials. Initially a compilation of some research studies is presented. This compilation allows to identify the different physical factors affecting particles breakage. Afterwards, from a modeling point of view, different constitutive models are studied in order to identify the different strategies proposed to take into account the particles breakage into the mechanical behavior of granular materials. Most of the constitutive models studied do not take into account the physical factors identified initially.

A particular model is found to present advances in this sense (Einav, 2007a). This model, developed in a thermodynamic framework, takes into account the grain-size distribution and microscopic quantities to model the influence of particles breakage on the mechanical behavior of granular materials. This model is studied through its application to a very crushable material. From this study two subjects are found to be important to focus on: i) the modeling of the grain-size distribution, and ii) the relation between frictional dissipation and particles breakage.

Concerning the first subject, an analytical model describing the variation of grain-size curve due to mechanical solicitations is proposed. Concerning the second, a framework to link frictional dissipation and particles breakage is developed.

In addition, the study of size-scale effect emerging from particles breakage is developed. It is shown that a size-scale rule allows, with a certain confidence, to estimate the mechanical behavior of crushable materials with large particles from a sample with smaller particles. This rule requires the existence of some defined conditions of similarity.

The final chapter of this thesis proposes a variation on the design of high CFRD and a methodology to estimate the mechanical parameters of the rockfill materials.

Chapter 2

Mechanical Behavior of Crushable Granular Materials Introduction

This research is looking for insights about the influence of particles breakage on the mechanical behavior of granular materials. The particles' breakage is evidenced by doing a set of granulometric tests, before and after solicitation. It is seen that the gradation curve varies from the initial grain size distribution to a new one in which the relative amount of bigger particles decreases and the amount of smaller particles increases. This subject has been treated with more interest in the domain of mining engineering than in geotechnics. In the former, the effect of crushing a particle into smaller particles is called the comminution process. This comminution process is led by the generation and spreading of cracks inside the particles. An explanation to this phenomena is given by the works of Griffith who proposed that the breakage of particles is due to the spread of internal micro cracks undergoing an increase in stress.

In this chapter the influence of different physical factors on particles' breakage is reviewed. Afterwards the influence of breakage on the mechanical behavior of rockfill materials is analyzed. Keeping in mind the practical applications through numerical models, different constitutive models integrating the particles breakage are studied. It is looked forward to identify the different strategies used in order to integrate particles' breakage phenomena into the mechanical behavior of a continuum media.

Description of grain-size distribution

The grain-size distribution is a statistical representation of the particles' size within a sample, its cumulative form is called the gradation curve. This curve has been employed currently in soil mechanics practice since its very beginnings for classification purposes. In this way and based on the percentage of mass passing some defined sieves a given material can be classified as clay, sand, gravel, cobbles, boulders, or rockfills. The general procedure is described by standards (see for instance -ASTM C136 ). The final result is a curve of retained material (Figure 2.1a), that in a cumulative form represents the gradation curve (Figure 2.1b).

The mathematical relationship between the retained and cumulative curves of a grain size distribution is given by equation 2.1

F (D) = D M Dm f (D)dD ∀D ∈ [D m , D M ] (2.1)
where f (D) corresponds to the relative mass of particles with diameter D in reference to the total mass

(f (D) = m D,D+dD m T
). In a practical way a mass of granular material is sieved through different mesh sizes with a decreasing opening length. Then for a given set of sieves the retained curve (f (D)) is constructed by plotting the retained mass in each size class (in percentage with respect to the total mass) against the opening of the sieve (size of particle class). It is usual to find in soil mechanics practice the description of particles' size through the gradation curve (F (D)), e.g., in filters design. In contrast, during this thesis the use of the grain-size distribution curve (f (D)) is done frequently.

Measures of grain crushing

As will be presented hereafter different factors have been defined to describe the amount of particles' breakage. The most employed is the so-called Marsal's factor. This factor is used in this chapter and also in chapter 4, 6 and 7. In contrast, in chapter 3, a different breakage factor is employed. At the author's knowledge there is no equivalence between the factors described on the following paragraphs.

Marsal's crushing factor

The methods for measuring the level of crushing in granular materials are based on the variation of the gradation curve or the retention curve. One of the first measures was proposed by [START_REF] Marsal | Research on the behavior of granular materials and rockfill samples[END_REF] using the retention curve. In this case the difference between the percentages of retained mass at initial (f 0 (D)) and final (f u (D)) stages in a test is calculated for each sieve. The amount of crushing is given by the sum of positive difference values and it is equal to parameter B, or:

B = 1 2 D M Dm |f 0 (D) -f f (D)|dD (2.2)
where | • | represents the absolute value of •. The procedure is illustrated in Figure 2.2. The factor B is a simple factor describing the amount of crushed material and has been currently used by engineers in soil mechanics practice.

Hardin's factors

Hardin (1985) proposes another factor to evaluate the amount of breakage. The calculus of this factor is based on the variation of the gradation curve (F (D)). The factors describing the grains breakage in [START_REF] Hardin | Crushing of soil particles[END_REF] are three, potential breakage, total breakage and relative breakage. Hardin takes into account that the required stress to break a particle is greater if its size is smaller. Then he sets the size of 0.074mm, corresponding to silt size, as the grain size at which particles are no more divisible. This breakage limit defines the first factor, called the potential breakage (B p ) which is defined by the area between the initial curve and a uniform curve corresponding to a size equal to 0.074mm (Figure 2.3a). The second factor relates the curves between the initial gradation curve and the final gradation curve and it is called the Total Breakage (B t ), see Figure 2.3b.

Finally the third factor represented by the relation 2.3 corresponds to the ratio between the two earlier factors, called relative breakage (B r ).

B r = B t B p (2.
3)

The magnitude of B r is always in the interval [0 1].

Other factors

Other factors have been defined, we can quote [START_REF] Lee | Compressibility and crushing of granular soil in anisotropic triaxial compression[END_REF] who defined a breakage factor corresponding to the relation between the D 151 before and after test. Einav (2007a) proposes factors quite close to Hardin's. However the theoretical final grain size distribution does not correspond to a distribution of silt size (0.074mm), instead of that, a final grain size distribution is assumed to be known. Einav defines also the potential breakage, total breakage and relative breakage as in Hardin's proposition. In Einav (2007a) the final grain size distribution is supposed to be fractal, based on several propositions described in literature. This factor is employed all along chapter 3 during the analysis of Einav's model.

Factors affecting Particles Breakage

The parameters affecting the particles breakage have been studied by [START_REF] Lee | Compressibility and crushing of granular soil in anisotropic triaxial compression[END_REF], [START_REF] Vesić | Behavior of granular materials under high stresses[END_REF], [START_REF] Marsal | Large scale testing of rockfill materials[END_REF], [START_REF] Marachi | Evaluation of properties of rockfill materials[END_REF], [START_REF] Hardin | Crushing of soil particles[END_REF], [START_REF] Bard | Comportement des Materiaux Granulaires Secs et a Liant Hydrocarboné[END_REF], [START_REF] Biarez | Influence de la granulometrie et de son evolution par rupture de grains sur le comportement mecanique de materiaux granulaires[END_REF], Nakata et al. (2001a) Oldecop and[START_REF] Oldecop | A model for rockfill compressibility[END_REF], [START_REF] Oldecop | Suction effects on rockfill compressibility[END_REF], among others.

As mentioned in the previous chapter, the breakage of particles is a function of different parameters. These parameters can be divided into three categories i) factors related to the particles, ii) factors related to the assembly, and iii) factors related to the mechanical solicitations. A summary of these factors is found in Table 2.1.

Mineral nature of particles

The influence of mineral nature of particles is logically related to the amount of crushing in granular materials. Materials with stronger particles will experience less breakage than those with weaker particles. Different works have been performed to identify the influence of mineral nature in particles breakage. e.g., [START_REF] Valdes | Influence of mineral composition on the simple shear response of sands: Experimental study[END_REF] tested under direct shear conditions mixtures of Ottawa sand (composed of rounded quartz grains) and calcareous sand (lower in resistance). The amount of [START_REF] Marachi | Strenght and deformation characteristics of rocfill materials[END_REF] for three different nature of rockfills, showed higher amount of breakage for Pyramid dam Material (originally angular and quarried from an argillite rock) than for a crushed basalt, other characteristics such as grain size distribution, confining stress and compacity being rather equal. [START_REF] Lo | Response of particulate materials at high pressures[END_REF] tested three different sands from different mineralogic origins : aluminium oxide, quartz and limestone. The initial gradations for the tests were nearly the same. The confining stresses were from 25psi (≈ 172kP a) to 1600psi (≈ 11M P a). The results showed that the amount of breakage increased with decrease in particles' strength. Consequently for the same confining pressure the dilatancy behavior decreased with decrease in particles strength.

Particles shape

Particles' shape is related to particles origin, alluvial deposits will be characterized by rounded particles while quarried or colluvial deposits are rather angular. According to [START_REF] Santamarina | Soil behaviour: The role of particle shape[END_REF] grain shape is established at three different scales: the global form, the scale of major surface features and the scale of surface roughness. Some studies on this subject have been carried out by [START_REF] Lee | Compressibility and crushing of granular soil in anisotropic triaxial compression[END_REF], [START_REF] Marsal | Embankment-Dam Engineering[END_REF] and [START_REF] Charles | The influence of confining pressure on the shear strength of compacted rockfill[END_REF]. In general trends, angular particles suffer more crushing than rounded ones. This is explained because the flaws of angular particles act like stress concentrators. According to [START_REF] Nakata | A probabilistic approach to sand particle crushing in the triaxial test[END_REF] testing singles particles of sand, (particles tested under compression between steel plates) the presence of flaws or asperities is observed in the stress-strain curve by a sequence of peaks lowers than the limit strength.

Particles Size

Variation of tensile strength with particles diameter represents the origin of scale effects treated in Chapter 6. This fact has been observed experimentally by [START_REF] Marsal | Embankment-Dam Engineering[END_REF], [START_REF] Barton | Shear strength of rockfill[END_REF], [START_REF] Lee | The angles of friction of granular fills[END_REF], Nakata et al. (2001b) among others. The explanation to this effect remains on the statistical distribution of defects inside the particles. Therefore, for the same mineral origin with an homogeneous distribution of flaws, a bigger particle will have more defects than a smaller one. Then for the same tensile stress it will have a higher probability of breakage. An experimental example is found in Figure 2.4 where the range of grading curve variation is showed for different size of particles.

It is remarked how the variation is smaller for the smaller particles, all other conditions being nearly equal. e.g., curves A, B and C.

Figure 2.4: Variation of grading curves, after [START_REF] Lee | Compressibility and crushing of granular soil in anisotropic triaxial compression[END_REF]. Materials with smaller maximum particles' size experience less breakage.

Gradation

The granular materials can have different gradation curves. A gradation curve is called uniform, when all the particles conforming the material have the same size. In contrast, when the material is composed of many different particles' size, it is called well graded. Two coefficients have been defined to describe the gradation of a material: i) the uniformity coefficient C u is equal to: D 60 /D 10 , wherein D 60 is the grain size corresponding to F (D) = 60%. In other words, 60% of the particles are finer and 40% coarser than this size. The second coefficient is the coefficient of curvature C c defined as

D 2 30 /D 10 • D 60 .
Gradation influence is related to the distribution of forces in the media. For instance, in the case of a regular assembly of uniform gradation under isotropic condition, all the particles will have almost the same stress state. Therefore, once the first particle breaks, there will be a chain effect, because the load that was carried by the first broken particle should be redistributed among the others, but as all the particles have nearly the same stress state, then the overload will generate again another particle breakage. This phenomena can explain the results of Figure 2.5 where a sample composed of identical glass spheres was tested under oedometric conditions, and after a "breakage pressure" an important increment in volumetric deformation was observed. Similar results with glass spheres assemblies have been reported by Nakata et al. (2001b), and for rounded silica gel by [START_REF] Feda | Notes on the effect of grain crushing on the granular soil behaviour[END_REF]. Nakata et al. (2001a) tested four different initial grain size distributions of silica sand, under oedometric conditions. All test samples had the same initial void ratio, and reached the same vertical stress. The crushing of particles was more important in nearly uniform gradations with higher maximum particles size. It seems that in the case of well graded materials the effect of particles' breakage is progressive, while in the case of uniform gradations the particles breakage happens more rapidly. This implies high deformation rates. From the data compilation in Appendix A, the breakage factor B (cf., section 2.2.1) can be drawn as a function of the uniformity coefficient C u = D 60 D 10 . Large uniformity coefficient indicates a material composed of a wide range of particle sizes. C u near to one indicates that the material is composed of particles of the same size. As mentioned before Figure 2.6 shows a lower breakage amount for high C u values.

Figure 2.5: Oedometric behavior of glass spheres, after [START_REF] Biarez | Influence de la granulometrie et de son evolution par rupture de grains sur le comportement mecanique de materiaux granulaires[END_REF]. e=void ratio, σ v =vertical effective stress.

Compacity

Compacity effects are observed experimentally in Figure 2.5. For instance, curves for particles' diameter of ∅ = 10mm with initial void ratios of e 0 = 0.72 and e 0 = 0.81 can be compared. The main difference is observed at the breakage stress, being around 6 MPa and 9 MPa respectively. However at high stresses, near to 20MPa the void ratio is almost the same. A similar observation has been signaled by [START_REF] Vesić | Behavior of granular materials under high stresses[END_REF] who tested the same material under two different initial void ratios. They found that after application of confining pressures higher than about 150 kg cm 2 (around 15MPa), both loose and dense samples reached the same void ratio. Therefore initial compacity has an effect in a lower-middle range of stresses, but at high stress conditions the same material reaches more or less the same void ratio.

Water content conditions

The experimental observation of water content effect in rockfills was first reported by [START_REF] Marsal | Large scale testing of rockfill materials[END_REF]. The experimental observations show that, in general, the increase of water content increases the particle crushing and hence the plastic deformations. A model taking into account the influence of water content is presented afterwards with the works of [START_REF] Oldecop | A model for rockfill compressibility[END_REF]. [START_REF] Oldecop | A model for rockfill compressibility[END_REF] explain the reduction of stiffness and increase in plastic deformations using the theory of stress corrosion. This theory states that the strained inter-atomic bonds at the tip of a crack are more vulnerable to the attack of a corrosive agent, such as water, than the unstressed material away from the tip. Then the corrosion reaction produces a weaker material, which is broken at a lower intensity of stresses than the un-corroded material. Frossard (2005a) establishes another theory called the capillary tightening process. In this process water content is present as menisci inside the cracks, then there exists a capillary force proportional to the inverse of the meniscus radius, i.e., a smaller radius leads to a high capillary tension and then to a reduction of tensional stress acting at the walls of the crack. When wetting occurs the radius of the meniscus increases and then the effective tensional stress on the walls increase. In consequence cracks propagate.

Stress level

The stress level in a soil sample is transmitted to the particles. Therefore with higher stress level, contact forces between particles are higher. In consequence there is a higher probability of particles breakage. Experimental results plotted in Figure 2.7 show how breakage factor B tends to increase with confining pressure.

Stress path

The stress path seems to affect the particles breakage. According to [START_REF] Lee | Compressibility and crushing of granular soil in anisotropic triaxial compression[END_REF] results, by increasing the stress ratio (k c = σ 1 σ 3 ) the breakage of particles increased also. This could suggest that the level of particle breakage depends on the magnitude of deviatoric stress. However, at the author's knowledge, no experimental works have been carried out to enlighten the effect of stress path on particles breakage. Breakage of particles have been evidenced under isotropic, oedometric, k c , and deviatoric conditions.

Time influence

Particles breakage seems not to happen instantaneously. The works of [START_REF] Vesić | Behavior of granular materials under high stresses[END_REF] reported an interesting plot showing the time taken for "consolidation" for different stress steps (Figure 2.8). They argued that as there is no reasonable interpretation for consolidation in such a material; (at least not from the usual point of view of hydrodynamic consolidation for clayey materials). Thus, this behavior is only explained by the character of deformation in question i.e., predominantly breakdown of particles. This suggests that particles breakage is a time dependent phenomena as suggested by [START_REF] Lee | Compressibility and crushing of granular soil in anisotropic triaxial compression[END_REF] who reported that each load increment was accompanied by loud cracking noises coming from the sample, and that the volume compression rate was decreasing. The time for zero compression rate was not determined, because the time between load increments was fixed to two hours.

Influence of particles breakage on mechanical behavior of granular materials

The factor influencing particles' breakage have been briefly described. Hereafter some studies concerning its influence on mechanical behavior of granular materials are reviewed. This review is classified Figure 2.8: Volumetric strain vs time during isotropic compression, after [START_REF] Vesić | Behavior of granular materials under high stresses[END_REF] according the stress path of the mechanical test involved for identification. The influence of particles breakage on the mechanical behavior of granular materials have been studied by [START_REF] Koerner | Effect of particle characteristics on soil strength[END_REF], [START_REF] Oldecop | A model for rockfill compressibility[END_REF], Varadarajan et al. (2006b), among others. In any case, the increase of particles breakage leads to higher deformability.

Oedometric and isotropic path

In oedometric or isotropic conditions the breakage of particles leads to an increase in compressibility. This increase can be more or less important according to different conditions, e.g., gradation, hydric conditions, particles' strength. The compressibility index or bulk modulus is reduced in presence of breakage. According to [START_REF] Biarez | Influence de la granulometrie et de son evolution par rupture de grains sur le comportement mecanique de materiaux granulaires[END_REF] for Hostun sand breakage settles beyond 1MPa of isotropic compression. Figure 2.9 shows a compressibility difference between the rounded and angular particles. Samples with angular particles experience more breakage than those with rounded ones, and consequently they present more compressibility. It is remarked that both type of particles were issued from the same mineral stock, and they had equal initial gradation curves.

Figure 2.9: Compressibility of granitic coarse sand tested under a k 0 = 1 path, after [START_REF] Lee | Compressibility and crushing of granular soil in anisotropic triaxial compression[END_REF]

Deviatoric path

The particles' breakage influences the strength resistance of soils. According to [START_REF] Vesić | Behavior of granular materials under high stresses[END_REF], at low confining pressures the strength envelope for the sand in dense state was curved. From tests at high confining stresses (up to 632 kg cm 2 ≈ 62M P a), where the influence of initial void ratio has no effect anymore, the friction angle seems to be constant and shear strength envelope is linear (results are resumed in Figure 2.10a). It is worth to mention that the final gradation curves obtained after shear stage are rather similar for the two highest confining pressures (Figure 2.10b). This suggests that breakage evolves in a non-linear form with the deviatoric stress. This evolution has been observed experimentally by [START_REF] Indraratna | Modelling of particle breakage of coarse aggregates incorporating strength and dilatancy[END_REF]. [START_REF] Vesić | Behavior of granular materials under high stresses[END_REF].

From data in Appendix A the influence of breakage in peak (φ peak ) friction angle can be drawn as shown in Figure 2.11. For each series of data the lower peak friction angles corresponds to the higher breakage factor, as defined by [START_REF] Marsal | Research on the behavior of granular materials and rockfill samples[END_REF].

Critical state

An important concept accepted nowadays in soils mechanics concerns the critical state2 [START_REF] Schofield | Critical State Soil Mechanics[END_REF] which is an idealization of the observed behavior of saturated remoulded clays in triaxial compression tests, and it is assumed to apply to normally consolidated soils. It states that soils and other granular materials, under continuous shear load will come into a well-defined critical state. At the onset of the critical state, shear strains occur without any further changes neither in mean effective stress nor deviatoric stress nor volumetric strain. However, as it is a concept extracted from idealized materials and applied into the framework of elasto-plastic behavior theory. It does not take into account particles breakage. Therefore, the parameters determining the critical state of a material can vary with breakage. [START_REF] Biarez | Influence de la granulometrie et de son evolution par rupture de grains sur le comportement mecanique de materiaux granulaires[END_REF] present the breakage effects on critical state, concluding than the critical state line on the plane of voids ratio vs logarithm of isotropic stress will shift toward lower values of voids ratio and will increase its slope. On the plane of isotropic vs deviatoric stress, their conclusions Many of the constitutive models looking for the description of soils' mechanics are based on the critical state concept. Especially those developed under the framework of the plasticity theory. Therefore, the variation of such a state will influence the response of the constitutive model.

Modeling mechanical behavior of crushable materials

Previously, different factors affecting particles breakage have been shown. Hereafter some constitutive models are studied in order to understand the strategies used to incorporate the aforementioned factors into the constitutive models of materials. The constitutive equations correspond to the mathematical description of the physical response of a material to different solicitations (e.g., loads, displacements). The constitutive equations are defined at the meso-scale (REV's scale) for a continuum material. They are supposed to take into account different phenomena (at meso and micro scales) that affect in a way or another the mechanical response of the media, e.g., hydraulic couplings for soils. There exists different frameworks to develop constitutive equations. Among them, the plasticity theory has been largely used in soil mechanics, and therefore it will be shortly reviewed in order to define a comparative framework supporting the analysis of different constitutive equations depicted afterwards.

General theory of plasticity

The elasto-plastic constitutive equations for modeling of materials were initially developed based on experimental results obtained on metal samples. Nowadays the application of this modeling technique goes beyond its initial domain, and it is commonly applied for structures or geo-structures analysis. It is pointed out that this framework does not take into account time effects or creep phenomena. The elasto-plastic constitutive laws are supported on two principal concepts.

• The yield surface 3 , which defines the onset of non reversible deformations, and

• The plastic flow rule, which defines the way of evolution for non reversible (plastic) deformations or any internal variables.

Yield surface

The initial elasticity domain is defined by a scalar function (F ) of the stresses (σ ij ), called yield function for which:

• F (σ ij ) < 0, the stress state is inside the space limited by the yield surface, and the response of the material is supposed to be purely elastic

• F (σ ij ) = 0, the stress state is on the limit of elasticity domain, and any other load increment will cause non reversible deformations

• F (σ ij ) > 0, the stress state is outside the space limited by the yield function. Theoretically is a state that can not be reached.

The yield function can evolve leading to the possibility of modeling hardening or even softening behaviors. In order to make the yield function evolutive it should be a function of other state variables, i.e., F (σ ij , α), where α is an evolutive hardening variable. According to the nature of the variable α the evolution of the yield function is different. According α the evolution of F can be isotropic or kinematic. The yield function determines when a stress state will generate non reversible deformations, but in order to know how they are going to happen the flow rule should be defined.

Flow rule

If the state variables of a material lead to F (σ ij , α) = 0 the stress state lies on the limit of elastic domain. In order to describe the behavior it should be identified if the variation of stress state will lead to a loading or unloading situation. Then: There is loading if:

F (σ ij , α) = 0 (2.4) ∂F ∂σ ij dσ ij > 0 (2.5) then dε ij = dε e ij + dε p ij (2.6)
There is unload if:

F (σ ij , α) = 0 (2.7) ∂F ∂σ ij dσ ij < 0 (2.8) then dε ij = dε e ij
(2.9)

3 critère de plasticité in French

The flow rule defines the non reversible deformation dε p ij as a function of σ ij , dσ ij and α. Using the principle of maximum plastic work given by Hill in 1950 [START_REF] Mestat | Lois de comportement des geomateriaux et modelisation par la methode des elements finis[END_REF] the increment of non elastic deformations are given by:

dε p ij = dλ ∂F ∂σ ij , dλ > 0 (2.10)
where dλ is called the plastic multiplier. The materials following the principle of maximum plastic work are called to be standard and their plastic potential is called associated. This is in general true for metals and other materials for which the plasticity criteria is independent of the mean stresses.

For soils and other materials called non standard, it is possible to define another plastic potential G different from the yield function and for which the non-elastic deformation is given by:

dε p ij = dλ ∂G ∂σ ij , dλ > 0 (2.11)
Therefore the flow rule is called non-associated to the yield surface and the plastic potential defines the direction of plastic deformations. It is important to note that the yield function F can have different forms in the space of stress. Notably, it can be an open or closed surface. In the case of a closed surface it means that the material can have non-reversible deformations even in isotropic or K 0 conditions. Modelling soil behavior requires to take into account particular phenomena as:

• Dilatancy

• Critical state concept

• Hardening and softening behaviors

The problem of particles breakage and its influence on mechanical behavior of geomaterials is not completely new, however the interest of its influence on rockfill behavior is relative recent and motivated by the several cases of high CFRDs showing problems related to important deformability. Nowadays several approaches looking for integrate the effects of particles breakage have been developed. Hereafter some of the constitutive models or theories incorporating particles breakage are described and analyzed. [START_REF] Daouadji | An elastoplastic model for granular materials taking into account grain breakage[END_REF] The proposition of [START_REF] Daouadji | An elastoplastic model for granular materials taking into account grain breakage[END_REF] is a modification of the Hujeux's model [START_REF] Hujeux | Une loi de comportement pour le chargement cyclique des sols. Genie Parasismique[END_REF] described completely in section B.1. This elasto-plastic constitutive model combines four mechanisms, three deviatoric and one isotropic, each mechanism is associated to a plane. Adopting the soil mechanics sign convention (compression positive), the deviatoric primary yield surface of the k plane is given by:

Model of

f k (σ , ε p v , r k ) = q k -sin φ pp • p k • F k • r k (2.12)
with:

F k = 1 -b ln p p c (2.13) p c = p c0 exp(β ε p v ) (2.14)
where, σ is the effective stress tensor, φ pp is the friction angle at the critical state, p is the mean effective stress and p c0 represents the critical state stress corresponding to the initial voids ratio. The internal variable r k , called degree of mobilized friction, is associated with the plastic deviatoric strain. This variable introduces the effect of shear hardening of the soil and permits the decomposition of the behavior domain into pseudo-elastic, hysteretic and mobilized.

Incorporation of breakage effect

The methodology to take into account the effect of particles breakage consists in modifying one parameter defining the yield surface. The variable identified as critical for the description of breakage influence corresponds to the initial consolidation pressure p c0 . Through different analyses [START_REF] Daouadji | An elastoplastic model for granular materials taking into account grain breakage[END_REF] concluded that the evolution of this parameter with plastic dissipation could reflect the influence of breakage. The formulation is of the form:

1 - p c0 p c0i = W * p B + W * p (2.15)
where

W * p = σ ij |dε ij | (2.16)
p c0i is the initial pressure without breakage and B is a parameter taking into account the effect of grain ruptures which depends on the initial properties (shape, size, distribution, etc.) of the grains. In the original version of Hujeux's model the parameter p c is the same for all the mechanisms, however [START_REF] Daouadji | Modélisation de l'influence de la rupture des grains sur le comportement des matériaux granulaires[END_REF] found that this conditions did not allow to represent correctly the influence of grains breakage. Therefore, it was necessary to express p c independently for isotropic and deviatoric mechanisms. In consequence equation 2.14 becomes :

p ck = p coi 1 - W * p k B k + W * p k exp(β ε p v ) (2.17)
where k represents any of the deviatoric mechanism, a similar formulation is also proposed for the isotropic mechanism. The inclusion of B k parameter in this formulation and the actualization of p c at each step will leads to a hyperbolic variation of the influence of breakage. According to [START_REF] Daouadji | Modélisation de l'influence de la rupture des grains sur le comportement des matériaux granulaires[END_REF] B k parameter takes into account the effect of grain ruptures including shape, size, distribution, etc. Therefore it would be expected that such a parameter changes if any of the factors influencing breakage change.

Analysis of the model

This approach for modeling crushable materials seems to obtain good agreements for the cases applied in [START_REF] Daouadji | An elastoplastic model for granular materials taking into account grain breakage[END_REF], where materials showed a hardening behavior without peak. Parameter B k should be estimated from the fitting of stress-strain and volumetric curves. Therefore it is an internal parameter that does not have a physical meaning, but which rather gathers the effect of different factors as mentioned by the authors. Unfortunately no relations between factors affecting breakage and B k parameter have been given. The model does not give estimations concerning the degree of breakage, or estimation of the new grading curve.

Barcelona's model for rockfills

What here is named as Barcelona's model corresponds to a series of models proposed initially for unsaturated soils. It is worth to mention that in this model the increase of moisture remains as the

p * p q s Cap Y ield
Y ield limits (saturated)

s i p(s)
Figure 2.12: Yield surfaces after [START_REF] Chávez | A constitutive model for crushed granular aggregates which incluyes suction effects[END_REF] main mechanism triggering breakage. This model has been performed to represent the mechanical behavior of rockfill-like materials. Initially under oedometric conditions notably after the paper of [START_REF] Oldecop | A model for rockfill compressibility[END_REF]. Afterwards it has been extended to triaxial paths with the work of [START_REF] Chávez | Estudio del Comportamiento triaxial de materiales granulares de tamaño medio con énfasis en la influencia de la succión[END_REF]. A review of these works is made in section B.3. In this model breakage of rock particles and fracture propagation are basic underlying mechanisms controlled by the relative humidity of the air filling the rockfill voids. From triaxial test [START_REF] Chávez | A constitutive model for crushed granular aggregates which incluyes suction effects[END_REF] and [START_REF] Chávez | Estudio del Comportamiento triaxial de materiales granulares de tamaño medio con énfasis en la influencia de la succión[END_REF] present the following experimental observations and conclusions:

• the shear strength envelope is curved,

• once the samples compacted, any subsequent loading immediately induces plastic strains and

• critical state relationship strongly depends on the prevailing suction Therefore their model, which is written in terms of total stresses, comprehends two yield functions as shown in Figure 2.12. The basic equations are resumed in Table 2.2. A more detailed description is given in Appendix B.4. As relative humidity is the main factor affecting the mechanical behavior of the materials, different mechanical parameters are expressed as a function of suction (s) (or its equivalent water content (w))4 . The suction variation will affect not only the parameters depending on it directly, but also to the yield functions, leading to a modification of elastic domain.

Analysis of the model

The model proposed by [START_REF] Oldecop | A model for rockfill compressibility[END_REF]; [START_REF] Chávez | A constitutive model for crushed granular aggregates which incluyes suction effects[END_REF] and [START_REF] Chávez | Estudio del Comportamiento triaxial de materiales granulares de tamaño medio con énfasis en la influencia de la succión[END_REF] takes into account the effect of water content on the mechanical behavior of rockfill-like materials.

The main idea supporting the effect of water content concerns the development of cracks inside the particles. Depending on the state of stress the development of cracks will lead to particles' breakage. This phenomenon, as seen before, causes the variation of the critical state and therefore an increase in plastic deformations. Despite the identification of breakage phenomenon as the main factor in rockfill plastic deformation, and its description using linear fracture mechanics, the model does not include parameters representing neither the strength of particles nor the concept of stress intensity factors.

f q = q p -η s f p = p(s) -p * ; ∀p < p y f p = p[λ i + λ d (s) -κ] -... ...p y λ d (s) -p * (λ i -κ); ∀p ≥ p y
p isotropic stress q deviatoric stress η s stress ratio s suction Plastic Potentials

g q = 1 2 q -1 6 q sin(ψ d ) g p = f p ψ d dilatancy
,s) = ε p q b(W P E )+ε p q dp * = v dε p p λ i -κ
The model requires near to 17 parameters, which are mainly internals. This is so because different parameters of the model are related phenomenologically to the variation of suction. The model does not give estimations concerning the degree of breakage, or estimation of the new grading curve.

2.5.4 Model of Kohgo et al. (2007a) The model proposed in Kohgo et al. (2007a) follows a similar way of development as the previous model. Initially proposed for unsaturated soils in [START_REF] Kohgo | Verification of the generalized elastoplastic model for unsaturated soils[END_REF] It is written in terms of effective stresses and adapted to the description of rockfill-like materials. As the model described in the previous section, the breakage of particles has been identified as the principal cause of plastic deformations, and the water content as the principal factor dealing with. The model is composed of two yield surfaces, two plastic potentials and the modification of the state surfaces by suction effects. Details about the constitutive equations are found in section B.6. This constitutive model is intended to reproduce two principal phenomena observed experimentally.

1. an increase in suction increases effective stresses, and 2. an increase in suction enhances yield stresses and affects the resistance to plastic deformations

The hardening behavior is isotropic and the hardening parameter concerns the volumetric strain. The definition of the actions of water depends on its configuration inside the granular material, in this sense [START_REF] Kohgo | Verification of the generalized elastoplastic model for unsaturated soils[END_REF] define the saturation condition as insular, pendular or fuzzy saturation. It is proposed that insular and pendular are the main conditions for rockfill materials. These conditions of saturation are reproduced by the soil-water retention curve. The main adaptation for rockfill materials concerns the modification of the state surface represented in the space of void ratio, logarithm of isotropic stress and suction.

Analysis of the model

The present model has been applied for the estimation of mechanical response of rockfill-like materials tested under oedometric and triaxial conditions by Kohgo et al. (2007b). The comparison between model and experimental results shows a good agreement, but this requires the definition of 29 parameters describing elasticity, plasticity, state surfaces and soil water retention curve. This model represents phenomenologically the effect of suction on the stress-strain and shear resistance of granular materials. The propositions consist mainly to variate the critical state surfaces as a function of suction. However, as in the previous models, this model does not take into account microscopical phenomena but rather a mean behavior of the REV. The model does not allow to estimate the degree of breakage, or estimation of the new grading curve.

Clastic yielding and hardening

A different trend that does not follow exactly the current formulation of elasto-plastic constitutive model has been developed by the works of [START_REF] Mcdowell | The fractal crushing of granular materials[END_REF]; [START_REF] Bolton | Clastic mechanics[END_REF] and [START_REF] Mcdowell | On the micromechanics of crushable aggregates[END_REF]. They present a conceptual relationship between the fragmentation and current concepts of soil mechanics leading to the proposition of two new concepts named : Clastic Yielding and Clastic Hardening. The work of [START_REF] Mcdowell | The fractal crushing of granular materials[END_REF] concerns a research on a twodimensional probabilistic model of fracture trying to incorporate the following characteristics: i) crushing strength of a particles is a function of its size, ii) high coordination number5 reduces the probability of fracture because the load is well distributed, and iii) particles are more likely to crush as stress on a sample of granular material increases. It is noted that this model does not deal with local equilibrium or kinematics, but simply with the probability of splitting of grains. The three characteristics described before are conjugated using the so-called Weibull probability. This probability defined by the equation 2.18 represents the probability of survival P s (V ) of a block of material with volume V , under a tension stress σ.

P s (V ) = exp - V V 0 σ σ 0 m (2.18)
In equation 2.18 V 0 represents a volume for which the probability of survival at σ = σ 0 is equal to 37% and m is called the Weibull modulus. In addition, the tensional stress in a particle within a granular assembly depends on the diameter of the particle (D) and the quantity of contacts with neighboring particles or coordination number (C). Taking these effects into account the equation 2.18 becomes (cf., McDowell et al. (1996)) :

P s (D) = exp    - D D 0 2 σ σ 0 m (C -2) a    where a = m (2.19)
where is a factor affecting the degree of importance of the coordination number on the tensile stress. Supposing that the amount of split particles is known as well as a given split pattern, then it is possible to estimate the amount of new particles created for a given increment in stress. Therefore the simulation of the grain crushing can then be linked to the energy required for the creation of new surfaces

ΓdS

Vs(1+e) , and then this energy could be included as a term in a plastic work equation as that of Cam Clay models (equation 2.20) [START_REF] Mcdowell | The fractal crushing of granular materials[END_REF]. From a work equilibrium equation we can link the increment of strain with the change in the gradation curve. In equations 2.206 and 2.21 dS is the change in surface area of a volume V s of solids distributed in a volume V s (1 + e) where e is the voids ratio, and Γ is the surface free-energy, and M is the slope of shear strength envelope in the plan pq and then a function of the friction angle φ.

qδ p q + p δ p v = M p δ p q
(2.20)

qδ p q + p δ p v = M p δ p q + ΓdS V s (1 + e)
(2.21) Equation 2.21 can be rewritten for the case of uni-dimensional compression, leading to a relation between variation of plastic void index (de p ) and creation of new surfaces (dS) due to particles crushing.

de p = - ΓdS (1 -µ)σV s where µ = 2 9 M (1 -2k 0 ) (2.22)
where σ is the macroscopic stress and k 0 is the lateral/axial effective stress ratio (k 0 ≈ 1sin(φ)).

Hence, the plastic deformations are related to the successive splitting of grains. The grains are allowed to split, probabilistically, depending on the applied stress and the co-ordination number and size of each particle. The application requires to set a stress level σ, and therefore to calculate the fracture probability (Equation 2.19) for each particle. This fracture probability is compared against a random number between 0 and 1. If the fracture probability is higher than the random number the particle is selected for fracture following a self similar pattern. An important hypothesis laying behind the model is the evolution of breakage (split pattern), which is self similar, this means that once a particle is broken it will give origin to two identical particles that keep the same proportions of the parent particle. This hypothesis is depicted graphically in Figure 2.13

D 0 i = 0 D 0 D 1 D 1 = D 0 √ 2 2 i = 1 D 2 D 2 = D 1 √ 2 2 i = 2 D 3 D 3 = D 2 √ 2 2 i = 3
Figure 2.13: Breakage pattern after [START_REF] Mcdowell | The fractal crushing of granular materials[END_REF] From their results (Figure 2.14) these researchers concluded that the gradation curve can be described by a fractal distribution with a fractal dimension7 that is constant at every time (slope of parallel curves in Figure 2.14c). The fractal dimension is then related to the number of particles in the sample.

The analysis just described has been developed for an initial uniform gradation (Figure 2.14a). The fractal dimension obtained is related to the hypothesis of breakage evolution, i.e., to the fact that particles split in a self similar pattern. It is interesting to note that the analytical relation found between the probability of fracture and the evolution of void index can be expressed in a similar form to that (experimental) employed regularly in soil mechanics, i.e., e = e λλ ln(p). It is also important to note that this probabilistic model starts from microscopical conditions to explain the behavior at the mesoscopic scale. In resume, by knowing the amount of particles surface generated for a given increment in stress the amount of variation in voids index can be estimated. However, this demands the knowledge of contacts between particles (coordination number), the stress level for each one, and the split pattern. For real materials as sands, gravels, etc. the coordination number is nearly impossible to estimate, it depends on the packing conditions, the grain size distribution, and the randomly distribution of particles in space. The stress state of each particle is equally difficult to estimate. However, the macroscopic results suggest an evolution of the grading curve of fractal type. This hypothesis can be reviewed based on real experiences on granular materials.

Analysis of the model

What has been presented here does not concern a complete constitutive model but rather a framework relating the grains crushing on the mechanical behavior of granular materials. The application directly to real cases is not straightforward, due to the requirement of definition of variables at the particles' level. However it has been shown that the use of a probability function could regroup the effect of particles' size and stress state dependency for the application to crushable granular materials. This framework relates the increase in particles' surface (due to crushing) to plastic deformations through a work equilibrium equation. The macroscopic variation of the grain size distribution proposed as a fractal pattern is a key point of the framework. It is remarked that this framework is able to approximate analytically the linear trend observed since many time ago on the plane e-ln(p). However, the model suggests that the smaller particles have a higher breakage probability due to the reduced co-ordination number. This is in general contradictory to the experimental results where the larger particles experience more breakage.

Fractal dimension for the materials of civil works

The fractal dimension (α) of the number of particles is an important parameter in the clastic theory. Therefore it is worth to check if this parameter can be verified for real materials. In this sense, we are interested in the verification of such a fractal dimension for gradations related to the regular works in civil engineering, and more specifically in the field of dams and embankments. This research is based on the experimental works of [START_REF] Marsal | Research on the behavior of granular materials and rockfill samples[END_REF] and [START_REF] Marachi | Strenght and deformation characteristics of rocfill materials[END_REF].

Initially we focus on the work presented by [START_REF] Marachi | Strenght and deformation characteristics of rocfill materials[END_REF], which consist of drained triaxial tests of rockfill-like materials. The work cited here, were focused on the study of the scale effect for rockfills. The scale effect study has been done by testing the same material with three parallel gradations with different maximum size particle (Figure 2.15a).

The material used in this experimental campaign comes from three different origins, the first material called Pyramid dam material corresponds to a quarried material using explosives with very angular particles. The rock itself was an argillite, the specific gravity was around 2.62 and did not vary for different particle sizes. The second material is called crushed basalt; this material is also quarry blasted but afterwards it is crushed to various smaller sizes in a crushing plant. The rock itself was a fine grained olivine basalt, containing some phenocrysts of olivine and plagioclase feldspar. The specific gravity was around 2.87. Finally the third material called Oroville dam material is a well-rounded to rounded material. The rock particles were mostly amphibolite, a meta-volcanic rock containing augite with small percentage of biotite and quartz.

From the experimental results of triaxial tests we obtain the stress-strain curves, volumetric deformation -axial deformation, and a final gradation curve. Therefore to estimate the fractal dimension the particles retained in a sieve are supposed approximately spherical and their diameter unique. This diameter corresponds to the average size between the mesh sieve size where the material was retained and that just before. As the density of each particle size is required the value of density reported by the researchers is used for all grain sizes. The procedure is described as follow:

1. estimate the grain-size distribution curve; 2. knowing the total mass of the sample, estimate the mass retained at each sieve; 3. calculate the mass of one single particle for each sieve supposing spherical particles; 4. determine the number of particles retained on each sieve; 5. plot the number of particles larger than a size D i against the log 2 ( D i Dmax );

6. find the linear relation that fits the best the data obtained, the slope of this relation corresponds to an approximation of the fractal dimension. [START_REF] Marachi | Strenght and deformation characteristics of rocfill materials[END_REF] The results of this procedure applied to the materials of Figure 2.15b are shown in Figure 2.16. Figure 2.16d shows how the estimated fractal dimension changes with confining stress for the particular case in analysis. The same observation is valid for all the materials tested in [START_REF] Marachi | Strenght and deformation characteristics of rocfill materials[END_REF] as seen in Figure 2.17. This particular observation does not agree with the conclusions of [START_REF] Mcdowell | The fractal crushing of granular materials[END_REF] (Figure 2.14c) where the fractal dimension is constant because of the self similarity pattern of crushing.

It is pointed out that the initial gradation used in this analysis is not similar to that of [START_REF] Mcdowell | The fractal crushing of granular materials[END_REF]. In order to study the influence of the initial gradation some experimental results of [START_REF] Marsal | Research on the behavior of granular materials and rockfill samples[END_REF] have been analysed. The results concern a material with an initial uniform grain size distribution. The material is referred to as "La Soledad" material. For this material all the mass is retained between the sieves of 2" and 3" it represents a practically uniform gradation. The results of the fractal analysis are shown in Figures 2. 18 and 2.19. Despite the similarity between the gradation curves of the model (Figure 2.14e) and experiments (Figure 2.18a), the fractal dimension is still dependent on the confining stress. For this particular analysis it is remarked that the grain-size distribution for the most important confining stress approaches very well a totally fractal type gradation. i.e., the linearity is verified.

It has bee shown that the evolution of gradation curves is not completely characterized by a fractal dimension. This is explained partially by the creation of finest particles since the beginning of a test. [START_REF] Marachi | Strenght and deformation characteristics of rocfill materials[END_REF] These finest particles can be the product of breakage of flaws or asperities. This phenomena precedes the breakage of the entire particle [START_REF] Nakata | A probabilistic approach to sand particle crushing in the triaxial test[END_REF], especially in angular particles. In this sense, the model proposed by [START_REF] Mcdowell | The fractal crushing of granular materials[END_REF] can not represent this phenomena because it deals only with major particles breakage, i.e., the breakage of the entire particle.

Continuum Breakage Mechanics

As described before the inclusion of the breakage phenomena into the description of the mechanical behavior of soils can be more suitably performed through a thermodynamic framework. In this sense Einav (2007a) proposed recently a thermodynamic framework for crushable materials that leaded to the proposition of a constitutive model. This model has been retained to a deeper analysis that will be done during the next chapter.

Conclusions

In this chapter it has been shown how several physical factors affect the level of crushing in granular materials, and how crushing influences the strength and deformability characteristics of rockfill-like materials. From a modeling point of view, the identification of particles crushing as an important mechanism in plastic strains generation has been remarked for all the models analysed. However there is no inclusion of physical or mechanical characteristics of particles. Therefore the inclusion of particles breakage is represented through internal hardening variables combined to a form of the total plastic work.

A different trend to classic plasticity theory applied to soil modeling has been analysed also. This Figure 2.19: Fractal dimension analysis for result of oedometric tests on dense material of "La Soledad", after [START_REF] Marsal | Research on the behavior of granular materials and rockfill samples[END_REF] trend given by the works of [START_REF] Mcdowell | The fractal crushing of granular materials[END_REF] includes real physical parameters as crushing resistance, size dependency and coordination number through a probabilistic approach. Although it is not a mechanical constitutive model per se, it shows under certain assumptions how to introduce the particles' crushing effect to the modeling of soil behavior. By taking into account real physical parameters at the particles level and to relate non-reversible deformations to the energy required in the creation of new surfaces, which have a real physical meaning. A key aspect concerns the evolution of the crushing pattern as self-similar, leading to obtain a constant fractal dimension. However this fractal dimension has not been verified for sixteen experimental results analysed. The calculus of the fractal dimension is mainly affected by the smallest particles.

The distribution of particles as a fractal distribution is a restrictive condition, especially taking into account that initial grain size distributions can be non-fractal. Therefore the evolution of materials can not be represented by a constant fractal dimension. Nevertheless, it seems that for initially uniform gradations a fractal form can shows up for high confining pressures.

The model just mentioned relates the creation of solids surface to the non-reversible deformations through a thermodynamic equilibrium equation. Hence, the estimation of the creation of smallest particles is very important. This is so because, for the same mass of material, when smaller are its particles higher is the total solid surface on it. Therefore, in the light of this model higher non-reversible deformations will be generated.

The relation between the creation of solid's surface and non-reversible deformations has been done through a thermodynamic equilibrium equation. This trend is also developed in a recently proposed framework named continuum breakage mechanics (Einav, 2007a). This framework is analysed during the next chapter to assess its general applicability.

Introduction

In cohesionless granular materials subjected to constraints the external forces acting over the boundaries are transmitted and distributed among the particles contacts. The particles configuration inside the assembly and their interactions (contacts) form chains of forces that are randomly distributed inside the material and that can evolve during loading. This fact has been verified by several researchers using discrete element analysis e.g., [START_REF] Radjai | Force distributions in dense two-dimensional granular systems[END_REF], [START_REF] Kruyt | Statistics of the elastic behaviour of granular materials[END_REF], [START_REF] Bagi | Statistical analysis of contact force components in random granular assemblies[END_REF], among others. If contact forces are weak with respect to the particles' strength the overall response of the material depends on the configuration of particles, particles' deformability and friction coefficient between particles. On the other hand, if contact forces are of the same magnitude than particles' strength the probability of breakage is high, and in consequence the mechanical response is affected by inclusion of particles breakage. Nowadays, several constitutive models developed mainly on plasticity theory have incorporated the effect of particles' breakage, e.g., [START_REF] Daouadji | An elastoplastic model for granular materials taking into account grain breakage[END_REF], [START_REF] Salim | A new elastoplastic constitutive model for coarse granular aggregates incorporating particle breakage[END_REF], [START_REF] Chávez | Estudio del Comportamiento triaxial de materiales granulares de tamaño medio con énfasis en la influencia de la succión[END_REF], Kohgo et al. (2007a). However the approach is still more phenomenological than physical. On a different trend, the author's belief is that a sound thermodynamic framework can lead to a closer approach to the physical phenomena. In this sense the aim of this chapter is to review the principal hypothesis developed in a relative recent approach called "Continuum Breakage Mechanics" (CBM ), initially proposed by Einav (2007a). The present review is done by implementing the aforementioned theory to a very crushable material. For the sake of clarity, in a first part some necessary concepts to fix the framework of the subject are defined. This is followed by the review of the principal hypothesis and formulation of CBM. Afterwards the application to a very crushable material is made, and from the simulation's results some internal variables are redefined so as to obtain a better response for the case in analysis.

The main motivation of analyzing this model concerns the inclusion of the grain size distribution as a factor influencing the mechanical response of the media, and the dependency of different quantities to grain's sizes.

Einav's CBM framework

The present study concerns the works of Einav (2007a,b), which, at the author's point of view, give an interesting framework to include the particle breakage into the constitutive models for granular materials. This framework is resumed along the following paragraphs, and therefore, no formal proofs

Breakage evolution

When breakage or crushing of particles happens the grain-size distribution curve changes. In consequence the gradation curve measured at different stages of a test can reveal information about the amount of material crushed during the test. As made by other authors, it is proposed that under a mechanical test the gradation curve is going to evolve from the initial grain-size distribution, F 0 , to a final grain-size distribution at ultimate state, F u . This final gradation curve is supposed to follow a fractal1 form as:

F u (D) = D D M 3-α (3.1)
with α being the fractal dimension and D M the maximum particle diameter. Taking this hypothesis the breakage factors initially proposed by [START_REF] Hardin | Crushing of soil particles[END_REF] are redefined in the form showed in Figure 3.1. In this case B t is the total breakage, B p is the potential breakage and B r is the relative breakage index. As seen in Figure 3.1, the indexes are related to the area between any two grain-size distribution curves or their ratio (B r ).

Statistical average

The breakage factor just defined is a function of the grain-size distribution curve, or more precisely, of its variation in a global form. Therefore, it is possible to define a fractional breakage B(D) to take into account the breakage of different particle fractions. It is supposed that this function is independent of the fraction, meaning that B(D) = B and is given by equation 3.2.

B = f 0 (D) -f (D) f 0 (D) -f u (D) (3.2)
where f 0 (D), f (D), f u (D) represents the initial, current and ultimate frequency distribution of the particles respectively. The statistical homogenization via statistical average is implemented using the grain-size distribution as a weighting average function over the microscopic variables. This allows to take into account the variation of the grain-size distribution into any variable depending on particle's size. For instance the statistical average of a microscopic variable A(D) is given by:

A ≡ A = D M Dm A(D)f (D)dD (3.3)
Equation 3.3 in combination with 3.2 can be reformulated as follow:

A ≡ A = A 0 (1 -B) + A u B (3.4)
where A ≡ A represents the average value of the function A(D), and

A 0 = D M Dm A(D)f 0 (D)dD (3.5) A u = D M Dm A(D)f u (D)dD (3.6)
A 0 and A u correspond respectively to the initial and ultimate average value of the variable A(D). Therefore the mean value of variable A (equation 3.4) at any time can be assessed by applying equations 3.5 and 3.6 and by knowing the amount of breakage B. Hereafter statistical average will be applied to different functions depending on particle's diameter.

Overview of the thermodynamic framework

The Einav's model is based on the first law of thermodynamics for continuous media:

W = δΨ + Φ (3.7) Φ ≥ 0 (3.8)
where Ψ and δΨ are the Helmholtz free energy and its increment, Φ is the non negative increment of energy dissipation and W = σ : δε is the increment of mechanical work done on the representative elementary volume (REV ). The use of the tilde symbol over W and Φ represents an increment, deliberately different from the proper notation of an increment δ. This is made specifically to highlight that while Ψ is a state function, which therefore has a proper differential, W and Φ are not and only their increments could be defined.

It is proposed that Ψ can be expressed in terms of its density distribution, as:

ψ(D, ε) = f ψ (D)ψ r (ε) (3.9)
where f ψ (D) is an energy split function equal to 1 for a reference diameter D r and ψ r (ε) is the free energy density depending only on strain. Concerning equation 3.9, by choosing f ψ (D) = D Dr n the initial and final split functions can be obtained using the statistical averages as follows:

m 0 = f ψ (D) 0 = D M Dm f ψ (D)f 0 (D)dD;
(3.10)

m u = f ψ (D) u = D M Dm f ψ (D)f u (D)dD (3.11)
where the subscript "0" refers to the initial state and "u" to the hypothetical ultimate fractal distribution. It is possible to select D r in a way that m 0 = 1 (see more details in Einav 2007a). Therefore, applying equation 3.4 to the energy split function leads to find the mean split energy function corresponding to a breakage index B as:

f ψ (D) = (1 -B) D M Dm D D r n f 0 (D)dD + B D M Dm D D r n f u (D)dD (3.12)
In consequence the increment of Helmholtz free energy (equation 3.8) can be associated to the breakage index B and its increment δB by :

δΨ ≡ ∂ψ r (ε) ∂ε [(1 -B)m 0 + Bm u ]δε + ψ r (ε)(m u -m 0 )δB (3.13)
It is found that the stored energy in the system is proportional to the energy function at the reference grain-size, and it may be expressed as:

ψ r (ε) = Ψ(B ≡ 0) m 0 = Ψ(B ≡ 0) (3.14)
The last term in equation 3.8 corresponds to the energy dissipation increment. By now only the effects of breakage are taken into account as dissipative mechanisms. Therefore the energy dissipation increment by breakage is given by: ΦB (B, δB) = E B δB ≥ 0 where δB is the increment of breakage.

The factor E B = ψ r (ε)(m 0m u ) is defined as the "breakage energy" describing the total stored energy that can be released from the system during fracture. Combining the equations just mentioned with equation 3.2 and after rearranging, the following expression is obtained:

ΦB = δE B (1 -B) -E B δB = δE * B (3.15)
which means that the energy dissipation from breakage ( ΦB ) is equal to the loss in residual breakage energy (δE * B ) (Einav, 2007a).

Mathematical representation of CBM

By now, the general framework has been described, but for application purposes the definition of several functions or variables should be done. Among them the Helmholtz free energy potential has to be defined. The application of the model given by [START_REF] Einav | Breakage mechanics -part II: Modelling granular materials[END_REF] has been restricted to the isotropic path. In this particular case the results of [START_REF] Walton | The effective elastic moduli of a random packing of spheres[END_REF] and [START_REF] Einav | Pressure-dependent elasticity and energy conservation in elastoplastic models for soils[END_REF] propose that the Helmholtz free energy potential can be given by:

ψ r (ε e v ) = p r (ξ(ε e v )/p r ) 2-m K(2 -m) (3.16) Ψ = p r (ξ(ε e v )/p r ) 2-m K(2 -m) (1 -ϑB) for 0 ≤ m < 1 (3.17) ϑ = 1 - J nu J n0 = 1 -m u (3.18)
where ξ(ε e v ) is the isotropic pressure for elastic behaviour, p r is a reference pressure, K and m are elastic constants of the particles and J n0 , J nu are the n th moments of the distributions f 0 (D) and f u (D). It is pointed out that m is a mechanical parameter and it should not be confused with statistical averages m 0 and m u . Equation 3.16 represents the potential of the elastic behaviour for isotropic loading, and in absence of breakage (only elasticity), the stress-strain relationship associated is given by:

ξ(ε e v ) = p r 1-m K(1 -m)(ε e v -ε v0 ) + 1 (3.19)
With the previous definitions the following relations are found for the isotropic pressure and the breakage energy:

p = ξ(ε e v )(1 -ϑB) (3.20) E B = ϑ p r (ξ(ε e v )/p r ) 2-m K(2 -m) (3.21)
For a loading path, at the initial stage the work of external forces is stored elastically in the different particles, and the total elastic energy over the range of particles is equivalent to the stored energy at the reference particle as suggested by equation 3.14. The model proposes a yield limit G B , which is a strain energy constant of the material. Therefore when the stored energy at the particles (E B ) reaches a value E B = G B breakage begins. The yield pressure corresponding to this energy is given by equation 3.22 and the breakage is proposed to evolve according to equation 3.23.

p y0 = p r 2-m G B K(2 -m) ϑp r (3.22) B = 1 - K(2 -m)G B ϑp m-1 r ξ(ε e v ) m-2 (3.23)
By now, only dissipation due to breakage has been developed, however the total dissipation is also composed of a frictional component. This frictional component is linked to the breakage dissipation and together establish a flow rule given by:

Φ * p Φ * B 2 = σ : δε p E B δB (3.24)
where

Φ * p = sin -1 (ω)Φ p (3.25) Φ * B = cos -1 (ω)Φ B (3.26)
Equations 3.25 and 3.26 represent the coupling relations between breakage dissipation and plastic dissipation. The parameter ω represents a coupling parameter. Then the increment of plastic dissipation is redefined to be linked with the elastic energy and for an isotropic loading is given by:

Φp = pδε p v = p G B E B (1 -B) 2 δε p v (3.27)

Modeling crushable materials with CBM

The model proposed in Einav (2007a) is applied in [START_REF] Einav | Breakage mechanics -part II: Modelling granular materials[END_REF] to simulate the response of a silica sand tested by Nakata et al. (2001a). The general stress-strain response showed good agreements between simulation results and experimental data. However the assessment of breakage evolution did not show the same degree of agreement. Taking into account that breakage evolution is an important parameter in the model it seemed important to evaluate its applicability in a broader sense, where the influence of the different hypothesis can be revealed.

In this section we applied the model to simulate the behaviour of the material tested by [START_REF] Bard | Comportement des Materiaux Granulaires Secs et a Liant Hydrocarboné[END_REF]. The material is a petroleum coke with an initial void ratio of e 0 = 2.34. It has been tested in oedometric path for vertical stresses in the range of 5 to 100MPa (Figure 3.2a). In consequence the variation of the grain-size distribution curve has been observed due to an important breakage effect (Figure 3.2b). [START_REF] Bard | Comportement des Materiaux Granulaires Secs et a Liant Hydrocarboné[END_REF] In order to apply the model given by CBM different parameters must be defined. These parameters are classified as elastic, breakage and coupling parameters. The procedure for the parameters' assessment is depicted in Figure 3.3a. We begin by the elastic parameters estimation ( K and m), afterwards the breakage parameters assessment (pressure of beginning of breakage (p y0 ) and ϑ) value, and finally the coupling parameter (ω). Once the parameters have been assessed the computing stage begins by defining a small increment in elastic deformation (Figure 3.3b), then the equivalent elastic pressure is calculated (using expression 3.19), this pressure is checked against the yield pressure for breakage p yo in equation 3.22. In the case of yielding, the breakage (equation 3.23) factor is calculated as well as the other variables depending on this factor, finally the plastic deformation is updated. A flowchart of this procedure is shown in Figure 3 Begin

ε e vi = ε e vi-1 + δε e vi δε e vi , e 0 , ϑ, K, m, ω ξ(ε e v ) i eq. 3.16 ξ(ε e v ) i > p y0 ? p i = ξ(ε e v ) i , B = 0 ε e vi ≥ ε e vlim ?
B eq. 3.23 p i eq. 3.20 To facilitate the parameters assessment task, an average loading curve have been defined (Figure 3.4). The estimation of an isotropic stress "p" from vertical stress is made using a k 0 value of 0.371. The value of k 0 is obtained by applying Jacky's formula with φ = 39 o (corresponding to a M = 1.6 reported by [START_REF] Bard | Comportement des Materiaux Granulaires Secs et a Liant Hydrocarboné[END_REF]). [START_REF] Bard | Comportement des Materiaux Granulaires Secs et a Liant Hydrocarboné[END_REF]

δB = B i -B i-1 ε vi = ε p vi + ε e vi δε p v eq. 3.24 plot p -ε v ε p vi = ε p vi-1 + δε p

Assessment of elastic behaviour

In order to assess the elastic behaviour equation 3.28 has been fitted to the data corresponding to the unload path of a test (Figure 3.4). The fitting procedure has been done by minimizing a least square error function in order to find the values of K and m. In elasticity the pressure is given by equation 3.19, which can be rewritten as:

ε e v -ε v0 = (p/p r ) 1-m -1 K(1 -m) (3.28)
The parameters describing the elastic behaviour correspond to: K = 5791 and m = 0.45. It is pointed out that the m value is near to 0.5 corresponding to the theory of elasticity and it lies well between the theoretical interval given by Einav (2007b) 1 3 -1 2 .

Assessment of breakage parameters

The ϑ factor is given by equation 3.18 and it is calculated based on the initial and final gradation curves. Initially, from the initial grain-size distribution curve the second moment of the distribution given by equation 3.29 is obtained, its value corresponds in this case to J 20 = 58.3mm 2 . Afterwards, the second moment of the final distribution (equation 3.30) is obtained analytically using the hypothesis of a final fractal distribution.

J 20 = D 3 M -D 3 m0 3(D M -D m0 ) (3.29) J 2u = 3 -α 5 -α D 2 M (3.30)
In equations 3.29 and 3.30 D M is the maximum grain diameter in the sample, D m0 is the smallest particle diameter in the initial grain-size distribution and α is the fractal dimension. In this case the second moment of the final distribution is equal to J 2u = 16.6mm 2 for a fractal dimension of α = 2.6. With this last parameter and using equation 3.18 a value of ϑ = 0.71 is obtained. Before going further it is interesting to observe how this hypothesis looks like. For doing so the estimated final grain-size distribution is drawn in Figure 3.5 together with the experimental curves. It is observed at first glance that the fractal hypothesis does not represent the real tendency of the final grain-size distribution for this particular material. [START_REF] Bard | Comportement des Materiaux Granulaires Secs et a Liant Hydrocarboné[END_REF] and hypothesis of final fractal distribution for α = 2.6 A discussion about this hypothesis will be made in section 3.4.2, by now the analysis is going to continue under the assumption that the final grain-size distribution curve is fractal.

The next step is to find the pressure at the beginning of breakage (p y0 ), in a consolidation curve it corresponds to the point where the experimental curve separates from the elastic model. By doing so a p y0 ≈ 36kP a is found. Consequently by using equation 3.22 G B ≈ 20 J m 3 is obtained.

Assessment of coupling parameter

Once the breakage and elastic parameters have been estimated, the following step is to define the coupling parameter ω, which determines the slope of the consolidation curve over the non elastic behaviour. The procedure employed for the estimation of this parameter is as previously, the minimization of an error function. Figure 3.6 shows three results. The first for the the best fit (B.F.) and two others corresponding to the model response for ±5% of variation on the B.F. of ω. The error estimator calculated (as shown in Appendix C) corresponds to 19.5% for the case with ω = 84 o . It is highlighted that hereafter the stress-strain plane corresponds to the volumetric strain (ε v ) -logarithm of isotropic pressure (p). The "most satisfactory result" from Figure 3.6 corresponds to ω = 84 o . However, it seems not to be enough representative of the mechanical behaviour of the material.

Hypotheses review

Up to now, the application of CBM to the representation of oedometric laboratory tests on petroleum coke has shown an unsatisfactory result. In consequence, during this section the principal hypotheses are reviewed and some propositions are made according to the experimental results.

Coupling parameter

The coupling parameter ω determines the slope for the post-breakage domain in the stress-strain curve. By now the proposition in Einav (2007a) consists of supposing this factor as constant. It has been seen that the ω parameter seems to be an important factor of the model. For instance, considering only qualitatively the differences among the B.F. and ±5% curves in Figure 3.6, it seems that the model is very sensible to the parameter ω. Therefore, its correct estimation is of crucial importance.

Breakage pressure

In CBM after reaching p y0 crushing settles and the strain-stress curve changes abruptly. However for the present case the variation of the strain-stress curve is more gentle than predicted by CBM. Therefore, in order to increase the representativity of CBM the p y0 can be increased. By doing so it is possible to find a pressure near to p y0 = 175kP a and corresponding to the intercept of the elastic curve with the extrapolation of the linear part of the consolidation curve. By doing so a new ω = 82 o is found. The result obtained is shown in Figure 3.7a showing as well two responses corresponding to ±5% of the previous value. The error estimator calculated as shown in annexes, corresponds to 7.1% for the case of ω = 82 o . As seen, the result obtained is more representative of the experiences and, vis-a-vis the response of the model, an implicit relation is inferred between the parameters p y0 and ω . The result showed here suggests that the progressive breakage can occurs less abruptly than predicted by the formulation. It seems that the model is not able to represent a gradual variation between the completely elastic domain and the domain including breakage. This characteristic is particulary observed for the material under analysis. This gradual variation can be issued from different reasons, for instance: i) in the case of a high initial void ratio, as in this case, because some frictional dissipation is achieved before beginning of breakage dissipation, or ii) the particles themselves can experience plastic deformation, which is not taken into account in the model. From a modeling point of view, this could be translated into an evolutive ω parameter or the inclusion of another plastic mechanism.

Breakage factor evolution

The model is based on the evolution of the breakage parameter B. Theoretically this parameter evolves according to equation 3.23. Then for each increment in deformation a new B is obtained. The experimental and theoretical evolution of B for the case in analysis is shown in Figure 3.7b. The calculus of B requires a final grain-size distribution. For the theoretical case it is supposed to be fractal, but for the experimental results it corresponds to the grain-size distribution curve after the test with the highest level of stress.

Partial conclusions and remarks

Modeling petroleum coke material with CBM has shown some deficiencies. The deviations of the model in reference to the experimental response could be due to one or several of the following reasons:

• The parameter ω would not be constant for different load steps, i.e., it is a function of pressure or deformation.

• The final grain-size distribution is not fractal.

• There might be a variation of some parameters considered as constants during the load, for instance : there exist an explicit relation between G B and p y0 , consequently as p y0 increases G B should vary.

Hypothesis improvement on CBM model

This section concerns a parametric study of the model. The model hypotheses are reviewed by including the experimental evolution for each parameter. It is pointed out that the comparison between model and experimental results are restricted to a domain of isotropic pressures lower than 100M P a.

Variation of coupling parameter

The first analysis concerns the variation of the parameter ω considering all the others parameters as constants. The analysis is made over the load stage only. In this opportunity the pressure at the beginning of breakage p y0 , is set as the pressure at the beginning of non elastic behaviour and corresponding to a value near to p yo = 36kP a as in section 3.2.3. The procedure is actually a back-analysis because the known data are the total deformation and the final pressure p. The elastic parameters are taken equals to those obtained from the fitting curve stage. It is possible to combine the expression for p defined by equation 3.20 with equation 3.23 in order to obtain:

p = ξ(ε e v ) -ϑ ξ(ε e v ) -ξ(ε e v ) m 2 K(2 -m)G B ϑp m-1 r (3.31)
which relates the pressure p with the elastic parameter ξ(ε v ). However by doing so an implicit formulation for ξ is obtained. But, ξ is the searched parameter. Therefore equation 3.31 is solved for each interval of pressures.

Then, for each consecutive pair of points obtained from the averaged curve a value of ξ(ε e v ) is calculated. Afterwards the estimation of the elastic deformations, the breakage and its increment, the plastic deformation and the breakage energy are done. Finally an explicit expression for ω can be established by combining equations 3.24, 3.27,3.25 and 3.26. It is reminded that ω is a coupling factor between the energy dissipation by both friction and breakage. This parameter determines the slope of the oedometric curve (Einav, 2007a).

The results obtained are plotted in figures 3.8a and 3.8b. As expected the value of ω angle increases until a maximum value corresponding to the linear branch of the curve. The results show a high agreement between the averaged experimental curve and the simulated one until the limit of comparison. The variation of angle ω suggests an evolution of the relationship between frictional dissipation and breakage dissipation.

Final gradation curve different from fractal

One of the main assumptions in the model is the existence of a final grain-size distribution of fractal type, but we have seen in Figure 3.5 that this assumption is not completely fulfilled. The influence of the final grain-size distribution is direct over the parameter J 2u and consequently on parameter ϑ. Therefore, it is proposed to examine the influence of a final grain-size distribution different from fractal on the model response. In a first time by performing a sensibility analysis concerning J 2u and afterwards by performing an analysis with the real values of the final grain-size distribution curve.

The sensibility to a variation in J 2u is modeled, by keeping constant the coupling parameter ω = 82 o and setting p yo = 175kP a. The simulation results are summarized in Figure 3.9a, and the error2 for each case in Figure 3.9b. The sensibility is not symmetric, this means that the influence of taking +5% or -5% (in the value of J 2u ) is not the same over the final response.

It is observed in Figure 3.9a that in the case that J 2u changes, a new value of ω is required to approximate the real response. Therefore the fact of taking a final distribution not equivalent to the real one will give an artificial value of ω. Now the J 2u is estimated (see Appendix C) for the real ultimate grain-size distribution shown in Figure 3.5. A value J 2u = 3.87mm 2 is found and a new value of ω = 75 o has been calculated fitting the experimental data. A new simulation is done using this last value and the response is observed in Figure 3.10a, one more time two more curves are plotted for ±5% of the fitted value.

The simulation results show significant improvements in comparison with the previous simulation with a fractal distribution (Figure 3.7a). However breakage evolution seems to get farther from the real evolution as seen in Figure 3.10b. This fact is explained because by reducing J 2u we are increasing ϑ and then B increases faster.

Including experimental B

It has been shown that the assumption of the final grain-size distribution different from the fractal one can lead to a different response, and then, it is an important parameter of the model. But unfortunately the evolution of the breakage parameter is not well assessed. However, from experimental results a reliable assessment of the breakage evolution can be measured and calculated as proposed by Einav (2007a). Therefore the experimental value of B could be included instead of calculating it from equation 3.23. It is highlighted that now the p y0 = 36kP a corresponding to the pressure where the experimental curve bifurcates from the curve of elastic behaviour. The first simulation is done by using the experimental value of J 2u and the results are plotted in Figure 3.11, when a value of ω = 71 o has been adopted. It is seen that a better response is obtained by comparing it with Figure 3.6 As mentioned in a previous section ω parameter could not be constant, and it has been shown that the use of a ω variable can improve the simulation results. ω variable means that the relationship between breakage dissipation and plastic dissipation is going to evolve over the test. By implementing this variation into the model the simulation results increase in agreement with experimental data as seen in Figure 3.12. Finally two values of ω were required to obtain a good representation of the model. A last verification simulation is made to observe the response of the modified model (including experimental final grain size distribution and breakage evolution) to the case that was studied by [START_REF] Einav | Breakage mechanics -part II: Modelling granular materials[END_REF]. The results are very satisfactory and the ω angle differs in only 1 o . (68 o for Einav (2007b) and 69 o here).

Conclusions

In this chapter the application of a thermodynamical model has been studied. The model proposes a number of hypotheses that, as shown, are not always fulfilled. In consequence its application is restricted. However, it has been shown that important improvements in the model response can be obtained by the inclusion of the real variation of some important parameters in accordance to experimental observations. Especially, by changing the evolution of the breakage factor B an important improvement has been obtained. The author's opinion is that an overestimation of the evolution of this parameter is given by the formulation of Einav (2007a) which considers a fractal breakage scheme. Thereafter, searching for a better explication of the breakage evolution could give important improvements to the model. The inclusion of the experimental variables (grain-size distribution and breakage evolution) reduce also the sensitivity of the model. It should be pointed out that the application showed here was done to a material showing a highly non linear behavior in the strain-stress response, and modeled in a large scale of stresses. These factors could make appear problems to any constitutive model.

One disadvantage presented in the model corresponds to the coupling between breakage and plastic deformations. Because in the case that the material dissipates only by friction, the response of the model is not realistic. In the model this means to chose a ω near to 90 o (just plastic dissipation), and then once the threshold of the elastic domains is reached, infinite deformations appear even for a small increment in stress, which is not physically possible taking into account the nature of the test simulated.

The authors opinion is that the model proposed by [START_REF] Einav | Breakage mechanics -part II: Modelling granular materials[END_REF] gives a good framework to develop a better understanding in the participation of breakage as a factor in deformation, but that frictional (or plastic) dissipation could be developed differently to its original proposition. As the inclusion of real evolution of the grain-size distribution can improve this kind of modeling approach, it is very important to propose a model able to represent the variation of the grain-size distribution curve under mechanical solicitations. Moreover, the coupling between frictional and breakage dissipation should be studied more deeply. Some insights in this subject will lead to a better modeling approach of the non reversible mechanisms.

Chapter 4

Micromechanical Analysis of Particles' Breakage Introduction

In chapter 2 we showed a framework to take into account particle breakage into the mechanical behavior of crushable granular materials. This framework [START_REF] Mcdowell | The fractal crushing of granular materials[END_REF] requires the definition of some physical factors at the grains scale. Moreover, in the previous chapter an interesting model relating the grains' crushing effect to the mechanical behaviour of granular materials given by Einav (2007a,b) was reviewed. Its capacity to represent the mechanical behaviour for a highly crushable material improved while the parameters included in the model were closer to the real (experimental) evolution. Among the different parameters studied, the evolution of the grain size distribution curve was highly important, therefore the correct assessment of gradation during loading steps could improve considerably this kind of approach. In this chapter we begin by a review of the factors affecting particles crushing from the micro-mechanical point of view. Afterwards, we propose a model that allows to estimate the evolution of the grain-size distribution for crushable materials under mechanical test.

Fracture phenomenon

Crushing modes

The crushing of particles can be schematically described by the three general processes [START_REF] Guyon | Du sac de billes au tas de sable[END_REF] depicted in Figure 4.1. The first one is the fracture schema corresponding to the case when high compression forces are acting around the grains (Figure 4.1a), the second is the attrition or piling, that appears when moderate compressive forces break the flaws of the main grain (Figure 4.1b). Finally, the third one corresponds to abrasion, that takes place when frictional forces are predominant in the interaction (Figure 4.1c).

However as a particle can experience different load cases, the final real pattern of fracture can be a combination of the three modes. For all the cases, the fracture arises after the creation and spreading of internal cracks. The issue concerning the creation and spreading of internal cracks in solids materials is studied in the field of fracture mechanics.

Crack propagation and fracture mechanics

In continuum mechanics, usually stresses and strains are used to describe the mechanical behavior of a solid (cf., [START_REF] Ling | Fracture Mechanics,With an Introduction to Micromechanics[END_REF], for instance). These quantities, likewise very important in fracture mechanics, can not always be directly applied for the characterization of fracture processes. One reason for this is that stresses become infinitely large at a crack tip. Another one follows from the simple fact that two cracks of different lengths behave differently when loaded by the same external stress. Under increasing load the longer crack will start to propagate at a lower stress than the shorter one. For these reasons additional quantities like stress intensity factors or the energy release rate have been introduced in fracture mechanics, because they are able to characterize the local stress state at a crack tip and the global behavior of the crack during propagation, respectively.

To simplify the problem of cracks propagation and fracture at the level of particles it is convenient to consider that particles experience brittle fracture, i.e., fracture happens without considerable plastic strains. The crack problem is a field discontinuity problem [START_REF] Maugin | The Thermomechanics of Plasticity and Fracture[END_REF]. Under the effect of traction the two faces of a crack are separated. It is admitted that no cohesion forces allow to bind the faces again. To simplify the problem it is commonly accepted that the discontinuity surface is plane and that a local frame can be attached to the endpoint of the crack called also crack front. Therefore corresponding to the direction of the displacement vector on the local plane there is possible to define three different modes of crack propagation, as illustrated in Figure 4.2.

Mode I denotes a symmetric crack opening with respect to the x 1x 2 -plane. Mode II is characterized by an antisymmetric separation of the crack surfaces due to relative displacements in x 2 -direction (normal to the crack front). Finally, mode III describes a separation due to relative displacements in x 1 -direction (tangential to the crack front).

From an experimental point of view [START_REF] Luong | Basic fracture modes and breakage resistance of granite, Chapter Mechanical Properties of Rocks and its analysis[END_REF] tested a granitic material under different configurations in order to estimate for each fracture mode its mechanical characteristics. One important conclusion is that the energy required for fracture as well as the strength were both lower in mode I than for the other modes. Therefore, it is assumed that for granular materials used in civil works (rockfill dams, embankments, roads, etc) the breakage of particles is dominated by the mode I.

x 2 x 3 x 1 (a) Mode I x 2 x 3 x 1 (b) Mode II x 2 x 3 x 1 (c) Mode III

Experimental results of grain breakage

When the stress level in a particle is higher than its strength the breakage happens. The significant factors determining the strength of a particle are the properties of the mineral components (rock's nature) and the physical weaknesses of the mass itself (fissures, cracks and cavities) [START_REF] Auvinet | Statistical model of grain breakage[END_REF]. As mentioned before an explanation to the breakage phenomenon of a single particle is done by the fracture mechanics. More precisely the mode I of fracture (Tensile forces acting orthogonal to the crack plane) is supposed to be the most probable mode for the cracks propagation within the grains. This approach is reflected in the work of [START_REF] Charles | The influence of confining pressure on the shear strength of compacted rockfill[END_REF] and also by [START_REF] Oldecop | A model for rockfill compressibility[END_REF]. The mechanics involved is supposed to be the same as in the so-called Brazilian test (Figure 4.3), where compression forces at two opposite bounds of the material generate tension constraints within the body and then the crack propagation occurs. The speed of cracks propagation is related to the intensity of external forces. Griffith, who associated the dissipation of energy during crack evolution with the surface energy of newly created crack area. With the works of C. Jaeger the role of the coordination number (number of contacts with neighboring particles) was revealed, increasing tensile resistance with the number of coordination.

There has been some theories trying to define a pattern of crushing in particles, but before going The first peak is related to the breakage of asperities and the second or major peak corresponds to breakage of the particle in several pieces. Afterwards, some triaxial tests have been done and some particles have been chosen and marked to trace the breakage pattern. As a conclusion, for a given material the probability of failure decreases as the type of failure increases. It means that for a given stress state there is more probability to survive to a failure of Type 3 than of Type 2 or 1. At the same time there is more probable to survive to failure of Type 2 than 1.

In a further article, Nakata et al. (2001a), deals with the crushing under uni-dimensional compression and concludes that:

• for the same material the yielding characteristics depend on the gradation of the material; for uniform gradation, yielding occurs in a more marked way than for well graded materials.

• in well graded silica, smaller particles experienced major fracture while that bigger particles experienced fracture of asperities.

• it was observed that the distribution curves for small size particles evolved towards a common line and that an initial well graded grain size distribution lying above this line suffered little particle crushing. Then, a hypothesis shows up which states that there exists a critical grain size distribution for which no particle crushing occurs in one dimensional compression.

However, the works of [START_REF] Nakata | A probabilistic approach to sand particle crushing in the triaxial test[END_REF] did not establish the way of crushing, i.e., the way of cracks development. In contrast, [START_REF] Schonert | Breakage of spheres and circular discs[END_REF] based on uniaxial test on glass discs and spheres concluded that breakage by compression follows the following scheme: i) releasing of cleavage cracks near one contact and crack propagations almost according to the original stress field, and ii) as the primary cracks in a compressed sphere or disc reach the opposite contact, then new cracks are created there and propagate in the opposite direction. These conclusions are drawn from high velocity photographs (frequency of 1M Hz) taken during laboratory tests. From the pictures and the schematic representation of the crack development, Figure 4.6 suggests that a kind of compression nail is formed at the contact and produces the beginning of cracks from the contact area, and according to the conclusions of [START_REF] Schonert | Breakage of spheres and circular discs[END_REF], once one crack reaches the other side of the sphere new cracks are propagated from this side. What is retained from this work is that once the crack is initiated there is a rapid evolution leading to the crushing of the entire particle. In a similar approach, Cheng et al. ( 2001) use a high resolution photography to follow the breakage of particles. Figure 4.7 shows the evolution of breakage for a particle under oedometric test for different vertical stresses. It is interesting to observe that the development of cracks in a particle is progressive and that all the particles are not created at the same time.

Through observation of civil works and experimental results, several researchers have showed that there are several factors affecting the level of crushing or fragmentation of grains (see section 1.4). In summary the speed of cracks propagation within the particles depends not only on the level of mechanical constraints but also on other factors, such as water content or particle shape. According to the works of Griffith the energy consumed to propagate a crack is proportional to the energy required to create a new surface, in this regard some theories have been developed to assess the way of crushing in granular materials.

As the distribution and nature of forces determine the way of crushing it is important to understand its distribution. In literature there exist nowadays several articles concerning the statistical description of forces and associated variables in a granular assembly. The breakage of a particle depends on the distribution of forces acting on it. For instance the breakage probability reduces if a particle is in contact with many particles around it. In micro-mechanics words it is translated as the increase in coordination number. In this sense it is important to identify correctly the variables affecting the micro-structure of granular materials. The next section intends to recall some of the works dealing with statistical description of granular matter, focusing on dense granular materials.

Statistics of dense granular systems

The mechanisms of force transmission in a granular assembly are due to particles' contacts, where in a simplified form only normal and tangential forces at the contact of grains are considered. This underlies the fact that electrical or Van Der Waals forces are neglected. Moreover the randomness of particles' distribution within a volume makes the distribution and orientation of forces random as well. Nowadays with the use of discrete element method (DEM) the interaction between particles can be seen more clearly, revealing the existence of strong and weak chains of force, which in absence of any cohesion are reduced to compression and shear forces at the particles contact.

In granular assemblies there are different variables that are better described from a statistical point of view due to their random distribution. Among the different variables we found the contact forces, coordination number and relative displacements.

Force distribution

As mentioned earlier the forces in a dense granular material are randomly distributed. [START_REF] Radjai | Force distributions in dense two-dimensional granular systems[END_REF], working with 2D discrete element simulations of circular particles under oedometric path, showed that the probability of the normalized normal and tangential contact forces could be described by two functions. A power function for forces lower than the mean value and an exponential decay function for contact forces higher than the mean. In the reference just cited, the presence of strong and weak force chains as shown in Figure 4.8 is also pointed out . Concerning the force values, they were found in a range of 6 orders of magnitude. It was also observed that 60% of the contacts carried a force lower than the mean. The DEM simulations mentioned are summarized in Table 4.1.

In the log-log plane the probability distributions for the normalized normal and tangential forces show similar characteristics as can be observed in Figure 4.9

From another standpoint [START_REF] Kruyt | Statistics of the elastic behaviour of granular materials[END_REF] studied the statistics of the elastic behaviour of granular materials through two-dimensional discrete elements models composed of 50.000 particles. The difference of this work is that the authors supposed bonded and non rotating particles, which leads to pure elasticity at the contacts. The simulations showed that the probability function of the normal component of displacement at contacts, resembled to the particles distribution, which corresponds to a log-normal distribution.

According to [START_REF] Bagi | Statistical analysis of contact force components in random granular assemblies[END_REF] based on previous works by others, the distribution of force magnitude is rather close to the exponential distribution in the range of values larger than the average force, while that for the range of small forces there is no a general agreement. Bagi's theoretical works leaded to find an exponential density function of the forces. This is done by minimizing the statistical entropy. The verification of this proposition was made by simulating isotropic assemblies under hydrostatic load in two-dimensional discrete elements models. Bagi found that the exponential distribution was a good approximation for assemblies with low friction coefficient or even frictionless but less representative in the case of important frictional coefficient between grains. The DEM simulations of [START_REF] Bagi | Statistical analysis of contact force components in random granular assemblies[END_REF] consisted of assemblies of 4800 circular particles with size ratio varying between 1 and 4. The porosity of all the numerical samples was of 14%.

In a further work, [START_REF] Radeke | On probability distributions of contact force magnitudes in loaded dense granular media[END_REF] continue the theoretical approach for the probability distribution of contact forces. The analysis is done in two and three dimensions, and again the exponential type seems to represent the best the forces distribution. The comparisons against numerical models of DEM type seem to reveal good agreement for the distribution.

The works of [START_REF] Ngan | Mechanical analog of temperature for the description of force distribution in static granular packings[END_REF][START_REF] Ngan | On distribution of contact forces in random granular packings[END_REF] and [START_REF] Chan | Statistical distribution of contact forces in packings of deformable spheres[END_REF] deduce a theoretical force distribution based on the entropy minimization of an assembly by supposing Hertzian contacts (i.e., only normal contact force). The results for two and three dimensional cases reflect also an exponential form.

In these works the expressions obtained analytically were compared against DEM simulations. As described by the authors the approximation of the theoretical and numerical probability density functions are in good agreement. However it should be reminded that there are two important hypotheses in these works, i). The granular media is composed by particles of almost the same diameter (only 10% in radius variation), and most important ii). the contact between particles is rigid and no tangential forces are produced (no friction). Moreover, [START_REF] Guises | Granular packing: numerical simulations and the characterisation of the effect of particle shape[END_REF] showed that the stress distribution in a granular assembly is affected by both the particles shape and friction. They remarked that at low forces the shape of the probability density function (pdf ) of contact forces (normal and tangential included) is affected by particle shape and diverges from a strictly exponential distribution. Whereas for high forces the shape of the pdf follows an exponential distribution form. The effect of frictional forces was not clearly identified. Nevertheless, the contact forces of frictionless systems seem to decrease slightly faster compared to the frictional ones.

There exists much literature about the estimation of a probability density function for contact forces. However there is still significant differences between them. It is worth to say that not all the simulations are performed under the same conditions. It seems generally accepted that forces higher than the mean force present an exponential distribution, but the form of the distribution for weaker forces seems to have no consensus. It is remarked that the most part of literature deals with a grain size distributions nearly uniform or with a small variation. In contrast the real grain size distribution employed in dams engineering is rather distributed along several orders of magnitude, i.e., from sand size to rock fragments.

Coordination number

The coordination number is defined as the number of contacts having a particle. It is one of the factors that controls the distribution of the stress within a granular medium. The coordination number is affected by factors like grain shape and friction coefficient between particles. The works of [START_REF] Guises | Granular packing: numerical simulations and the characterisation of the effect of particle shape[END_REF] focused on this issue. They tested numerically different configuration of particles. From rounded particles to elliptical particles with an elongation factor up to 5. The numerical simulations using two dimensional FEM/DEM analysis, were carried out to simulate the formation or deposition of a granular assembly. i.e., they started from a completely free configuration and let the grains to fall under the action of gravity. The results reveal that both, particles shape and friction coefficient affect the coordination number differently. For elongation factors lower than 1.8 the coordination number seems to be linearly related to particles shape (elongation factor) while for factors beyond 1.8 the coordination number seems to be nearly constant. It was found also that for frictional particles the coordination number attained was lower than that reached in the frictionless case.

The study of micromechanics allows to enhance the understanding of soils behavior at the mesoscopic scale. At the particles scale different local variables can be defined. In practical situations only mesoscopic mean values of stresses or displacements can be measured (for instance in laboratory tests) but particles breakage happens at the microscopic scale. Therefore, it is interesting to study the manners to define mesoscopic (REV) variables from the microscopic ones.

From micro to macro state variables

One of the research domains in micro-mechanics consists in finding the equivalent "macro" (equivalent continuum) variables from the micro variables, e.g., to find the Cauchy stress tensor knowing the forces distribution. As described by [START_REF] Bagi | Stress and strain in granular assemblies[END_REF] the mechanical state of the assembly and its evolution can exactly be described and predicted if the following characteristics are fully known:

• position and geometry of each grain

• displacements (translations and rotations) of each grain

• contact forces

• material properties of the individual grains However it seems that such a detailed description is not necessary, and is very complicated to define from a practical point of view. Instead, the aim of microstructural approach is to find macro-level state variables through the proper averaging of micro variables.

Stress variables

The study of microscopic forces has leaded to a rather well accepted definition of the cauchy stress tensor σ as a function of variables at the contact between particles (e.g., [START_REF] Weber | Recherches concernant les contraintes intergranulaires dans les milieux pulvérulents[END_REF][START_REF] Bagi | Stress and strain in granular assemblies[END_REF] as:

σ = 1 S c∈S f c ⊗ l c for 2D (4.1a) σ = 1 V c∈V f c ⊗ l c for 3D (4.1b)
where S or V represent the surface or volume domain respectively, f c represents the contact vector force c and l c the vector connecting the centers of the two grains forming the contact c as seen in Figure 4.10.

δ ω c 1 δ ω c 2 l c v c 1 v c 2 c 1 2 n c 1

Kinematic variables

The microstructural definition of the strain tensor is not as easy to find as that of stress, and there is not an expression generally accepted. There exist many different ideas to interpret the strain tensor from a microstructural point of view. Most of them belong to one of the following two approaches [START_REF] Bagi | Analysis of microstructural strain tensors for granular assemblies[END_REF]:

1. Strains based on an equivalent continuum: In this approach the assembly is replaced by a continuous domain formed generally by polygonal cells having as corners the particles centers.

The definition of which particles form the cell is different for each approach. Therefore the gradient of a suitable translation field defined by the translation of particles' centers is associated to the macroscopic strain. The propositions made by [START_REF] Bagi | Statistical analysis of contact force components in random granular assemblies[END_REF][START_REF] Bagi | Analysis of microstructural strain tensors for granular assemblies[END_REF], [START_REF] Kruyt | Statistics of the elastic behaviour of granular materials[END_REF], [START_REF] Cambou | Change of scale in granular materials[END_REF] among others belong to this approach.

2. Best-fit strains: a gradient tensor is found giving the smallest deviation from the characteristic displacements of the system. Depending on the kind of displacement to be approximated different versions of best-fit strains are gained.

An interesting work comparing 10 expressions for the strain tensor based on micromechanical considerations is exposed in [START_REF] Bagi | Analysis of microstructural strain tensors for granular assemblies[END_REF].

State variables

The geometrical configuration of particles inside a granular assembly represents a characteristic state of the media. Macroscopically variables as void ratio, or equivalents, can be defined. Microscopically, two kinds of variables can be identified [START_REF] Cambou | Change of scale in granular materials[END_REF]: a scalar N c , called the coordination number, defined as the number of contacts per particle, and a tensor, called the fabric tensor which allows the description of the anisotropy of the media. Different fabric tensors have been defined in the literature [START_REF] Cambou | Micromechanics of Granular Materials, Chapter Multiscale Techniques for Granular Materials[END_REF]. For instance :

• Contact fabric tensor

H ij = n i n j = 1 N N c=1 n c i n c j = Ω P (n)n i n j .dΩ
• Branch fabric tensor

H ij = l i l j = 1 N N c=1 l c i l c j = Ω P (n)l i l j .dΩ
• Combined fabric tensor

H ij = n i l j = 1 N N c=1 n c i l c j = Ω P (n)n i l j .dΩ
where N is the number of contacts in the considered sample, P (n) is a distribution function, n and l correspond respectively to the normal vector of the contact and the vector connecting the particles centers of the particles in contact, as seen in Figure 4.10. The usefulness of a fabric tensor remains on the description of anisotropic materials, for instance the elasticity tensor can be a function of a given fabric tensor as described by [START_REF] Oda | Mechanics of Granular Materials: An Introduction[END_REF].

The expressions described before for the stress, kinematic or state variables allow to express mesoscopic quantities based on microscopic variables. However, the inverse way seems to be much more difficult, i.e., to define microscopic quantities based on measured mesoscopic or macroscopic variables as stresses or displacements. This would mean, for instance, to invert equation 4.1. We could think that in order to know the amount of breakage in a granular assembly we should know if each single particle breaks or not. This would lead to determine the forces acting in each particle, and to compare them against a breakage criteria. However, this is not necessary because we are not interested in describing the breakage of each single particle, but the breakage of the entire population. i.e., a mean tendency of the assembly. Therefore, it is thought that the description of a given microscopic variable, as contact forces or elastic energy, can be done through a probability density function of the variable. This approach will be developed afterwards in section 4.5 during the development of a model describing the crushing of particles. In the following, a review of some models proposed in literature aimed the description of particles crushing is done. This literature review will let us to give a short "state of the art".

Models for grading evolution

As mentioned in the introduction to this chapter, the main interest concerns the assessment of grain size distribution evolution. In literature there exist a few models already proposed for this goal, here a short review and analysis is proposed for some of them. It seems that this topic has been underestimated in geotechnical engineering. The variation of grain size distribution has been treated as a merely qualitative observation instead of a meaningful factor affecting soils' mechanics. In a neighbor branch of engineering, in mining engineering the process of size reduction of grains is called the comminution theory. It is conceived as a process in which a feed material is transformed into a product material. Two main comminution theories exist in mining engineering, named Rittinger and Kick theories. In a broad sense, both theories relate the energy required for size reduction of a solid to the particles sizes of the feed and product. The parameters for each one of the theories are found experimentally. In mining engineering the interest is to find the energy required to change the whole material of size a to a size b, which is expressed by the aforementioned theories. However, our interest is to know how an amount of energy will change the material a into smaller pieces, not necessarily of the same size. [START_REF] Fukumoto | A grading equation for decomposed granite soil[END_REF][START_REF] Fukumoto | Particle breakage characteristics of granular soils[END_REF] The original paper of [START_REF] Fukumoto | A grading equation for decomposed granite soil[END_REF] dealt with the description of the grain size distribution curves for any granular material. This led him to propose the description of the grain size distribution curve as a geometric progression. This representation was called the grading equation. Afterwards [START_REF] Fukumoto | Particle breakage characteristics of granular soils[END_REF] showed the application of the grading equation to several mechanical tests. The variables describing the geometrical progression were related to the mechanical solicitations (loads, energy, etc) phenomenologically. Therefore, they depended completely on the nature of test from where they had been estimated. In consequence, the generalization of this kind of approach seems difficult. [START_REF] Liu | Modelling of interparticle breakage[END_REF] This approach is based on the Weibull probability distribution. This model is based on three main hypotheses as follows:

Model of

Model of

1. The breakage of particles depends only on the energy transmitted to them and on the particle size. Therefore the fraction of broken particles S(X, E) and the breakage function B(x, X, E) can be represented in terms of the energy absorption E and the particle size X. The energy absorption expresses the consumed energy per unit mass.

2. The breakage of particles obeys the S and B functions determined in experiments with monosized particle beds. The feed size (initial grain size distribution) does not influence directly these functions.

3. The interaction between different size fractions is describable by an energy split function k, which expresses the distribution of the energy to the size fractions in the bed. The k function depends on the feed distribution f (X) and on the stressing energy;

k = k(X, E, f (X)).
The formulation is made by two balance equations, one for mass and one for energy.

p i = (1 -S i )f i + i-1 j=1 B i,j S j f j (i = 1...n) (4.2) k j = E j E (4.3) n j=1 k j f j = 1 (4.4) S j = S j (E j ), B i,j = b i,j (E j ), k j = k j (X j , E, f 1 , f 2 , .., f n ) (4.5)
where f i , p i are the mass fraction of size i in feed or product respectively; k j : energy split factor, ratio of the energy absorption in the class j related to the energy absorption of the bed. The percentage of breakage particles S corresponds to the probability of breakage given by a Weibull distribution and corresponding to:

S = S ∞ 1 -exp - E E c β (4.6)
where S ∞ (X) is the boundary value of S(X); E c (X) is a characteristic absorption energy and β the curve shape parameter. [START_REF] Liu | Modelling of interparticle breakage[END_REF] estimate the probability of breakage with tests in mono-sized beds charged in uniaxial compression. The breakage function B is also estimated from experimental results with a statistical analysis, however this analysis makes appear several parameters to be determined. The energy split function k is calculated by back analysis and an expression like equation 4.7 is proposed.

ln k = -c ln X X e -d ln X X e 2 (4.7)
where X e is the particle size for k = 1, c is the slope at X e and d is the curve shape parameter.

The application of this model is no longer developed here. However, it is pointed out that the results of [START_REF] Liu | Modelling of interparticle breakage[END_REF] are very satisfactory, especially because they assessed very well the product curves for different feed distributions. Nevertheless the estimation of the energy split function is not easy and depends on several factors that are rather internal to the model than physical parameters of the solids.

Grain size evolution model

In this section a new model intended to represent the grain size evolution of granular material is proposed. This model combines a probabilistic analysis gathers to fundamental principles as mass and energy conservation. The evolution of the grain size distribution is based on the fracture probability of particles. It is convenient to depict the mechanical reasoning supporting this model, which is as follows: It is known that when a representative elementary volume (REV ) of granular material is loaded two main types of deformation occur. Reversible and non-reversible deformations. From a thermodynamic point of view the reversible deformations are associated to a recoverable (elastic) energy. In contrast non-recoverable are associated to energy dissipation. The elastic energy stored at the REV comprehends the contribution of each single particle composing it. The non-reversible deformation of the REV is associated to the relative displacements of particles, which dissipates energy by friction. Another component of non-reversible deformations is the plastic strains of particles. However, it is supposed that particles behave as a brittle material and therefore no plastic strains are developed at the particles. This hypothesis of brittle material has been indeed supported by different experimental works at least for rockfill-like materials.

As particles are supposed to have an elastic -brittle behavior only a quantity of energy can be stored before it is converted or dissipated by breakage. A given increment in elastic energy is split through the different particles sizes. However, different particles sizes have different fracture thresholds. Then, once this energy threshold is reached, breakage appears. In consequence, new particles are created, with different and smaller masses and higher breakage thresholds. The creation of these smaller particles allows to store the increment of elastic energy. So it is proposed that a given elastic energy should be in equilibrium with a given grain size distribution able to store the elastic energy without particles breakage.

The mathematical description of the model requires the definition of the following equations:

• Mass conservation.

• The fracture criteria for particles, which should take into account fracture size dependency and is obtained from a probabilistic fracture analysis. • The way that the energy at the REV scale is connected to the particles scale, which corresponds to an energy distribution function. • The way that particles split or comminution rule.

Mass conservation

When a particle breaks it is divided into smaller particles, i.e., smaller masses. The total broken mass is distributed over the range of smaller sizes. In that case a new grain size distribution is obtained. Therefore for a new energy state at a time t + 1 the grain size distribution is given by:

f t+1 (D i ) = (1 -S(U e p , D i , y 1 , ..y n )) • f t (D i ) + Dmax D>D i B(D) • S(U e p , D, y 1 , ..y n ) • f t (x)dx (4.8)
Previous expression is an intuitive relationship for mass conservation proposed by [START_REF] Liu | Modelling of interparticle breakage[END_REF]. It represents that for a time t + 1 the mass percentage of particles size D i is equal to the percentage presented at time t of the class i (f t (D i )), less the broken amount S(U e p , D i , y 1 , ..y n ) during the passage from t to t + 1, plus the contributions (B(D)) of the bigger particles (D > D i ) that have been also broken.

The function B(D) represents the pattern of crushing, i.e., how the mass of a broken particle is distributed into smaller particles. S(U e p , D i , y 1 , ..y n ) is a function describing the amount of broken particles for a size D i . It is a function of the energy at the particle level U e p , the size of the particles and possibly other parameters (y 1 , .., y n ). As mentioned earlier the function S(U e p , D i , y 1 , ..y n ) is obtained from a probabilistic fracture analysis.

Probabilistic fracture analysis

In contrast to a deterministic analysis the use of statistical distributions and probability principles in a mechanical problem allows to take into account the variability of loads and resistance characteristics.

Probability distributions

The frequency by which some quantity x occurs is described by the probability density function (pdf ) f (x). If we assume that the values of x are restraints to the interval [x m , ∞] the probability distribution is given by

F (x) = x xm f (x)dx (4.9)
Therefore the probability P that a random variable X lies in the interval x m ≤ X ≤ x is:

P (X ≤ x) = F (x) (4.10)
Since P can attain values between 0 and 1 the following relations hold:

P (X < ∞) = ∞ xm f (x)dx = 1, P (X ≥ x) = 1 -F (x), (4.11) P (a ≤ X ≤ b) = b a f (x)dx = F (b) -F (a).
The mean ( or expectation ) value X of a random variable and the variance var X are defined as

X = ∞ xm xf (x)dx (4.12) var X = ∞ xm [x -X] 2 f (x)dx (4.13)
The latter can also be described as the average square deviation from the mean value X. The square root of the variance is called the standard deviation: σ = √ var X. There exists several expressions for the probability density functions, depending if the random variable x is discrete or continuous, some examples for continuous variables are: normal distribution, lognormal distribution, gamma distribution, Weibull distribution, etc. Special attention is given to the latter distribution that is used widely in statistical description of brittle fracture.

1 f u (U e ) P f F B (U e ) f B (U e ) U e f /F

Fracture probability

In order to compute the fracture probability of a particle we need to know two probability distributions, the fracture criteria distribution and the constraint distribution. First, let us assume that a given elastic energy U e defines the rupture criteria distribution. Then the probability density function (pdf ) for rupture is given by f B (U e ), which defines the probability that the body fails for a value of U e . Second, the distribution of the variable U e follows another pdf f y (U e ). It represents the probability to find a body with a constraint of U e .

The probability that the fracture criteria is smaller than a certain value U e is given by:

P (U e ≤ U e B ) = F B (U e ) = U e B 0 f B (U e )dU e (4.14)
Correspondingly, f u (U e )dU e is the probability for a constraint over the interval U e < U e i < U e + dU e . The product

dP f = F B (U e )f u (U e )dU e (4.15)
then describes the probability that both applies, i.e., that the component fails. Integration over all possible levels of constraints finally yields the total failure probability:

P f = ∞ 0 F B (U e )f u (U e
)dU e (4.16) Therefore the pdf s f u (U e ) and f B (U e ) should be defined. The precedent procedure is illustrated in Figure 4.11.

Fracture criteria

distribution -f B (U e )
The fracture of an isotropic material subjected to a given stress (σ t ) is due to the presence of internal defects. The defects are assumed to be distributed in a statistically homogeneous manner, i.e., the probability of the occurrence of a characteristic defect is everywhere the same. In addition it is assumed that total failure of the material occurs if a single defect becomes critical and starts to grow. The probability that at a tensile stress σ t a certain volume V does not contain a critical defect is given by F * (V ). The respective probability for some other arbitrary volume V 1 (which does not contain V ) is F * (V 1 ). If events V and V 1 are independent of each other the probability that no critical defect is found in V + V 1 is given by:

F * (V + V 1 ) = F * (V )F * (V 1 ) (4.17)
Differentiation at fixed V 1 and subsequent division by 4.17 yields to:

dF * (V + V 1 ) dV = dF * (V ) dV F * (V 1 ), dF * (V +V 1 ) dV F * (V + V 1 ) = dF * (V ) dV F * (V ) (4.18) or d dV ln[F * (V + V 1 )] = d dV ln[F * (V )] = -c (4.19)
Here, c is a variable which depends only on stress: c = c(σ). Integration and knowing that F * (0) = 1 finally leads to the probability that no critical defect is contained in the volume V :

F * (V ) = e -cV (4.20)
In consequence, the probability that V contains a critical defect is F (V ) = 1e -c(σ)V . As only one critical defect leads to failure, then the probability of failure P f is given by:

P f = 1 -e -c(σ)V (4.21)
The relation c(σ) is an unknown function for which often the empirical representation according to [START_REF] Weibull | A statistical distribution of wide applicability[END_REF] is expressed by:

c(σ) = 1 V 0 σ-σu σ 0 m for σ > σ u 0 for σ ≤ σ u (4.22)
where V 0 and σ 0 are normalization parameters and σ u is the threshold stress below which the fracture probability is zero. For simplicity σ u is frequently set to zero. Then by combining equations 4.21 and 4.22 the probability of fracture is given by a Weibull distribution:

P f (V, σ) = 1 -F * (V, σ) = 1 -exp - V V 0 σ σ 0 m (4.23)
Instead of expressing c in equation 4.20 as a function of the stress, we can express c as a function of the elastic energy in the particle U e p . This is physically coherent, because the energy required for the spread of a crack in a particle should be present in the media. In this case we obtain:

P f (V, U e p ) = 1 -exp - V V 0 U e p U 0 mu (4.24)
where U 0 and m u are parameters of the pdf. It will be shown in section 4.6.1 that the experimental data are better represented by a more general expression given by:

P f (D, U e p |D 0 , U 0 ) = 1 -exp - D D 0 γ U e p U 0 mu (4.25)
where (D) is the diameter of particles and γ and D 0 are parameter of the pdf.

Energy distribution functionf u (U e )

A granular assembly can schematically be reduced to a group of particles interacting each other through contacts. Therefore, when a representative elementary volume (REV) of a cohesionless granular material is subjected to constraints in a quasi-static state, i.e., near equilibrium, the forces at the REV boundaries are distributed among the contacts between particles. Then, the overall response (e.g., deformation) of the REV depends on the redistribution of these forces at each increment of constraints. The redistribution of forces is translated in terms of particles movement and the deformation of particles themselves.

Deformation energy in particles assemblies

As mentioned before, the forces at the REV contact are transmitted from one particle to other through the contact. In a simplified form this will generate two types of forces at the particles contact: i) normal compression forces and ii) shear tangential forces, the relation between these two forces can be described by a Mohr-Coulomb relationship. As external forces increase, normal forces do the same, and in terms of energy, more energy (elastic energy) is stored in the particles. At the same time, in pseudo-static behaviour, an equilibrium of shear forces should also be reached. The elastic properties observed at the REV level are the direct consequence of the elastic energy stored in the particles.

Taking into account energy conservation it could be established that:

dW REV = dΦ REV + dU e REV (4.26)
which means that any given increment of the total work (dW REV ) at the REV's scale is related to an increment in both elastic energy (dU e REV ) and dissipation (dΦ REV ). These increments are related to the microscopic increments of the same nature at the particles scale. Which can be expressed as:

dΦ REV = c dΦ i (4.27) dU e REV = p dU e i (4.28)
where the sum over c represents the sum over all the contacts, and the sum over p represents the sum over all the particles composing the REV. At any stress state the elastic energy at the REV is related somehow to the elastic energy of the particles. Let us suppose that for a given increment in stress at the REV boundaries a new stress state is reached (σ), and in consequence an increment of strains (dε) is observed. This implies an increase in the mechanical elastic energy, dU e REV = σ • dε e , for the REV scale. But at the particles contact it is no that evident to find the increment of elastic energy. Therefore, hereafter a development intended to relate the elastic energy at the particles level to the elastic energy at the REV level is proposed.

The elastic strain energy in a representative elementary volume is of the form:

U e REV = σ ij dε e ij (4.29)
going deeper into the particles composing the REV we find a force network connecting the particles through contacts. If an elastic behaviour of the particles is supposed, the elastic energy stored in a single particle p is given in a similar way to 4.29 by: U e p = σ p ij dε e,p ij (4.30)

As both stress and strains are intensive variables (e.g., see [START_REF] Maugin | The Thermomechanics of Plasticity and Fracture[END_REF] then expressions in equations 4.29 and 4.30 are energy densities by unit volume. Then equation 4.26 is rewritten as:

U e REV = 1 V REV Np p=1 U e p • V p (4.31)
where V p is the particles volume, V REV the REV's volume and N p the total number of particles. The elastic energy U e p is supposed to be a random variable distributed according to a probability density function f u (U e p ). Moreover, equation 4.31 can be discretized for some particles subsets. For instance, the range of particles is divided into different groups having approximatively the same size. These subsets are called classes. Then the total volume is composed of N c classes, and for each class there are N i number of particles. In this sense, the probability density function of the elastic energy for a class i is given by f i u (U e j ), which corresponds to the relative amount of particles ( n j N i ) found in the class with an elastic energy of U e j (i.e., f i u (U e j ) =

n j N i ).
Thus the total energy by class is obtained by:

U e c i • V i = j U e j • V p(i) • n j (4.32) = j U e j • V p(i) • N i • f i u (U e j ) (4.33) = N i • V p(i) j U e j • f i u (U e j ) (4.34) U e c i • V i = N i • V p(i) • Ū e p(i) (4.35)
Equation 4.32 represents that the energy of a class i is equal to the sum over the range of energies j of all the particles (n j ) having an energy U e j proportional to the volume of the particles V p(i) , that for each class is taken as constant. In equations 4.33 to 4.35 N i is the total number of particles in class i and Ū e p(i) is the mean elastic energy of the particles in the class i. From equation 4.35 we finally obtain:

U e c i = N i V p(i) V i • Ū e p(i) (4.36) But N i • V p(i)
corresponds to the total volume of class i, then:

U e c i = Ū e p(i) (4.37)
Thus, the total elastic energy at the REV level is equal to the sum of elastic energy from each class, which is expressed by:

U e REV = 1 V REV N c i=1 V i • U e i = N c i=1 ρ t ρ p m i m t • Ū e p(i) (4.38)
where ρ t is the total density of the REV, ρ p is the density of particles, and m i mt is the relative amount of mass of class i to the total mass of the sample, i.e., the retained mass in a sieve for a grain size analysis. Under dry conditions, particles density is related to REV density by:

ρ rev = ρ p (1 -n) (4.39)
where n is the porosity. Then 4.38 becomes:

U e REV = (1 -n) N c i=1 m i m t • Ū e p(i) (4.40)
Equation 4.40 expresses the relationship between the elastic energy at the REV and particles scale.

Comminution rule

The comminution rule is the pattern of particle sizes obtained after breakage, i.e., once a group of particles of class i is broken, they form a new distribution of particles following the defined comminution rule. It is pointed out that it is not mandatory that every single particle in a class breaks exactly in the same way, but that the class follows the comminution rule. 1996). From an analysis of several grain-size distributions subjected to different stress states, it was found that the grain-size distributions do not follow completely a fractal distribution (see section 2.5.6). However, it seems that a fractal proposition is more acceptable for an initial uniform gradation curve i.e., only one class of particles.

Based on this case it is assumed that the comminution rule for a given particle class follows a fractal distribution. A fractal distribution is represented by an expression of the type:

N (D > D 0 ) = C • D -α 0 (4.41)
that expresses the number of particles N that have diameters D above the size D 0 , C being a constant of proportionality and α the fractal dimension. In a conventional grain size analysis (by sieving), the mass of all the particles that are bigger than sieve mesh D i and finer than the sieve mesh of size D j is normally expressed by:

M (D i < D ≤ D j ) = D j D i sρD 3 dN (D), (4.42)
where s is a shape factor, ρ is the density of particles of size D. The relative mass retained between sieve meshes D i and D j in relation to the total mass is obtained from:

f (D) = M (D i < D ≤ D j ) M d (D m < D ≤ D M ) = D j D i sρD 3 dN (D) D M Dm sρD 3 dN (D) , (4.43)
where D m and D M are the minimum and maximum diameter of the particles. For practical purposes the minimum diameter D m is chosen to the level of clay size ≈ 0.078mm. From 4.41 the differential dN (D) is found as (e.g. [START_REF] Mcdowell | The fractal crushing of granular materials[END_REF]:

dN (D) = CαD -α-1 dD (4.44)
which converts 4.43 in:

Comminution rule Fractal dimension α f (D) = D j D i sρD 3 CαD -α-1 dD D M Dm sρD 3 CαD -α-1 dD , (4.45)
as s, ρ, C and α are supposed constants for all the particles, the relative mass is finally given by:

D j D i D 2-α dD D M Dm D 2-α dD = D 3-α j -D 3-α i D 3-α M -D 3-α m (4.46)
Last equation represents the frequency density function of particles (in mass) of fractal type, having as maximum particle size D M . The comminution rule proposes that an initial uniform grain-size distribution of size D M becomes a fractal distribution after breakage. Then the comminution rule for a class M states that the distribution of broken particles with size greater than i and finer than j is given by:

f (D) = B ij,M = D (3-α) j -D (3-α) i D (3-α) M -D (3-α) m for D m ≤ D i < D j < D M (4.47)
Equation 4.47 should respect the mass conservation, then:

D M Dm B(D)dD = 1 (4.48)
which is verified.

Application of the model

In order to apply the model described in the previous section, the results of Nakata et al. (2001a,b) are used. These works have been selected because they present the evolution of the grain size distribution under mechanical tests, as well as breakage characteristics of individual grains. The model proposed here requires the estimation of different parameters as described in Table 4.2. It requires also the information concerning porosity and initial grain size distribution. However the required parameters are not directly obtained from the publication of Nakata et al. (2001a), for instance the fracture probability is given in terms of stresses and not in terms of energy as proposed in the present model. Therefore, a pre-processing of the data should be performed in order to obtain the required parameters.

The pre-processing procedure requires the assumption of several hypothesis due to the lack of crucial information.

Fracture probability estimation

In order to estimate the fracture probability of a class of particles two pdfs are needed. On one hand, the fracture probability distribution of individual particles, and, on the other hand the elastic energy distribution over the particles belonging to a class. In the works just cited, the authors estimated a fracture probability for individual particles of Silica sand under uniaxial compression tests. The probability is expressed in terms of tensile failure stress obtained by using the well known Jaeger's relationship between a compression force and a failure tensile stress, given by:

σ f = F f d 2 0 (4.49)
where σ f is the failure tensile stress, F f is the normal force at failure and d 0 is the initial diameter of the particle.

The values of σ f are plotted as a function of the particles diameter in Figure 4.12. The statistical analysis of these data leads to the fracture probability shown in Figure 4.13. Its shape suggests a form similar to the Weibull probability distribution. Moreover from Figure 4.12 the size dependency of fracture stress is evidenced. The results of Nakata et al. (2001b) are reinterpreted here in terms of energy. As for each particle under uniaxial compression a force-displacement (or stress -strain) curve is measured. The elastic energy is obtained by integrating this curve until fracture. However as the complete strain-stress curves are not available some hypothesis are made in order to estimate the required probability. In order to do that the following assumptions are made:

• The particles tested are supposed to have a spherical shape

• The mechanical behavior of particles is linear elastic until failure

• The contact between the piston and the particle is supposed to be an Hertzian contact, therefore there is no friction between the piston and the particle, then shear forces are neglected. According to Hertzian contact [START_REF] Gras | Tribologie, Principes et solutions industrielles[END_REF] the elastic energy in a spherical particle due to a normal force F N is given by:

U e p = 2 5 • K • F 5 3 N (4.50)
where

K = 3 4 (k 1 + k 2 ) 2 3 • d 0 2 -1 3 (4.51)
and

k i = 1 -ν 2 i E i ; i = 1, 2 (4.52)
in terms of displacement

U e p = 2 5 • K -3 2 • δ 5 2 N (4.53)
The index i represents the parameters for the piston and the silica particle, 1 and 2 respectively. Here the following values obtained from the literature are used: ν 1 = ν steel = 0.28, ν 2 = ν silice = 0.17, E 1 = E steel = 200GP a, for the Young's modulus of silica sand a fitting procedure was done with reference to the force-displacement curve published in Nakata et al. (2001b) and showed in Figure 4.14. The Young's modulus for silica was calculated around E 2 = E silice = 3GP a. Equations 4.50 and 4.49 are applied to data in Figure 4.12 in order to obtain the elastic energy stored until breakage. The statistical analysis of the elastic energy is given for three classes of particles in Figure 4.15. The first class from 0.25 to 0.3mm with a mean diameter d = 0.278mm, the second from 0.6 to 0.7mm with a mean diameter d = 0.649mm and the third 1.4 to 1.7mm with a mean diameter d = 1.47mm. From these data the breakage probability distribution is fitted according to the expression 4.25. The best fitted curve is shown also in Figure 4.15.

Elastic energy distribution estimation

From the works of [START_REF] Kruyt | Statistics of the elastic behaviour of granular materials[END_REF] the relative normal displacement between particles ∆ n is distributed according to a normal distribution. From the works of Hertz (see for instance Gras, Breakage Probability-P ) the elastic energy stored in a spherical particle is proportional to a power of the contact deformation. Combining these two facts the pdf for the elastic energy tends toward a lognormal distribution. Therefore, as a first approximation the elastic energy is supposed to follow a lognormal distribution of the form.

f P f = 1 -exp - d d 0 γ U e p U 0 mu d 0 = 0.335mm U 0 = 2.11e + 007 γ = 1.4 m u = 1.06 d 1.47mm d 0.649mm d 0.278mm
f u (U e p |λ, ζ) = 1 U e p • ζ √ 2π exp -(ln(U e p ) -λ) 2 2ζ 2 (4.54)
where λ and ζ are respectively the mean and standard deviation of the associated normal distribution. Equation 4.54 corresponds to the elastic energy distribution over a given particle size or class.

The fracture probability for individual particles has been defined previously (equation 4.25), and their parameters have been estimated from experimental results (Figure 4.15). Concerning the elastic energy distribution function the expression 4.54 is retained.

In order to apply this expression the parameters λ and ζ should be estimated. These parameters could vary from one class to another. At the author's knowledge there is not study involving the distribution of elastic energy by class in a granular assembly and due to the limitation in time this research can not be developed in the context of the present thesis. Therefore in order to assess λ and ζ let us propose that the mean elastic energy over the particles of a class ( Ū e p ) is proportional to the mass of the class (m i ), therefore:

Ū e p = m i m t U * (4.55)
where U * is parameter to be determined. By combining equations 4.40 and 4.55 the following expression is found:

U e REV = (1 -n) Nc i=1 m i m t • m i m t U * (4.56)
The parameter U * is supposed constant for all the classes of particles then it can be computed from equation 4.57 The negative solution is rejected, which leads to:

U * = U e REV (1 -n) Nc
λ = -1 + 1 + 2a 2 ln Ū e p a 2 (4.58)
In resume, under the simplifications adopted, the procedure to calculate the parameters for the elastic energy distribution function is:

Elastic energy at REV

The tests chosen to validate the model are those published by Nakata et al. (2001b). The material tested was a silica sand in two different gradations (Figure 4.16). The initial void ratio (e) for each one of the tests varies slightly, and as attended the uniform distribution has a higher voids ratio than the non uniform gradation. The tests were performed under vertical stress varying from around 20KPa until almost 100MPa (Figure 4.17 ). As in classical soil mechanics the initial linear branch of the curve e-log(σ v ) is considered representative of the soil's elastic behaviour. In terms of strainstress curves (Figure 4.18) the slope of the elastic behaviour (k) is considered as the inverse of elastic modulus. Therefore for each vertical stress σ v we can compute its associated elastic deformation ε e v and therefore to calculate the elastic energy as:

U e REV = σ ij dε e ij (4.59)
Therefore for each point in Figure 4.18 the elastic energy was computed as the area below the curve of the elastic behavior. The elastic energy for the last vertical stress is represented in Figure 4.19 by the red area. Once the elastic energy at the REV scale is known, the computation of the parameters describing the elastic energy distribution function can be calculated and the model can be applied. However, the change of grain size distribution curve is not instantaneous, but rather progressive. Therefore to take into account this effect the elastic energy for a given stress state was divided in subintervals in order to simulate a progressive increase in energy and hence, a progressive variation of grain size distribution.

Simulation results

For a given calculus step t a fixed increment in volumetric strain (∆ε v ) is added to the volumetric strain from t -1. The total volumetric strain ε v (t) is used to estimate the vertical stress from the stress-strain curve (Fig. 4.18). Afterwards the energy at the REV scale is calculated using equation 4.59. Then a new grain size distribution curve is obtained. The parameters α and a are obtained from a fitting procedure. Values corresponding to α = 2.1 and a = 0.6 were found. Using these values the model was calculated for different stress levels. The model has been applied for two different initial grain size distributions, and the agreement between simulation and experimental results seems to be very acceptable taking into account the different assumptions made. The results are shown in Figures 4.20a and 4.20b .

It is observed that the model overestimate slightly the amount of breakage for the largest particles. But in general trends the gradation curve is well represented.

Sensitivity analysis

In order to apply the model several simplifications have been adopted as consequence of lack of some required information. Therefore a sensitivity analysis concerning mainly the internal variables a and α have been carried out. The objective is to recognize the influence over these parameters on the model's response. Variable a is a dimensionless parameter relating the mean value and the standard deviation of the energy distribution function on a class (i.e., the coefficient of variation). Figure 4.21a shows that at lowers values of a a higher breakage is obtained. It is so because a narrow standard deviation suggest that almost all the particles have nearly the same amount of elastic energy, and once the rupture threshold is reached, almost all the particles will break. In contrast for a high standard deviation (i.e., larger coefficient of variation) there is a wide difference of energy among the particles, then only a few will break once the rupture threshold is reached. Concerning the fractal dimension (α). Its effect is well known, for a higher dimension, more amount of finest particles are obtained. Therefore the results found in Figure 4.21b are coherent.

It is believed that the coefficient of variation can change from one size to another, but by now it is taken as constant for all the size of particles.

Conclusions

In the previous chapter the reliable estimation of grading evolution led to improve the response of a mechanical model. In the present chapter a model describing the grading evolution has been proposed. The model includes the fracture size dependency of grains. It is based on the initial grain size distribution and assumes that a fracture criteria can be established in terms of elastic energy. From a micro-mechanical point of view the randomness of the variables describing the micro-mechanical state of the granular assembly can be described using probability density functions. Based on this a fracture criteria is established from a probabilistic fracture analysis. Despite the hypotheses made due to the lack of some crucial information, the description of grading evolution seems to be satisfactory for the study cases analysed. It is believed that a deeper understanding on the statistics in "highly"1 poly-dispersed materials can enlighten some of the hypotheses made here. We are especially interested in the distribution of elastic energy by size of particles in a poly-dispersed assembly.

Energetic Aspects in Continuum and Granular Materials

One of the main conclusions in Chapter 3 concerned the need to understand the link between energy dissipation by friction and fracture phenomena. It is believed that a thermodynamic analysis of experimental tests can enlighten this relationship. A thermodynamic framework allows to analyze and to take into account different phenomena occurring in a material, including the microscopic phenomena. In this chapter the principles of thermodynamics applied to continuum media are briefly reminded in order to settle the framework for further analysis on dissipation mechanisms.

Thermodynamics was initially developed as a theory that allowed to understand the heat's capacity to generate motion. Furthermore, it evolved into a theory that describes transformation of states of matter in general, motion generated by heat was then a consequence of certain transformations. A concise evolution of thermodynamics can be found in [START_REF] Kondepudi | Modern Thermodynamics, From heat engines to dissipative structures[END_REF].

Thermodynamics is founded on essentially two fundamental laws, one concerning energy and the other concerning entropy. Every system is associated with an energy and an entropy. When matter undergoes transformation from one state to another, the total energy remains unchanged or is conserved; the total entropy, in contrast, can only increase or, in idealized cases, remain unchanged.

In this chapter a very brief introduction to some main concepts in continuum media thermodynamics and energy dissipation in granular materials is done. This will allow to establish a general framework of analysis. Using this framework the soil sample in a triaxial test is idealized as a closed system exchanging only mechanical energy with the surrounding. The relation between the energy received (input power) and given back (output power) to the surroundings is analysed in a plot called "power-in/power-out". In this plot the pure frictional dissipation mechanism (materials without breakage of particles) is characterized by a linear relationship between the input and the output power. The inclusion of breakage as a phenomena leading to energy dissipation is identified by the excess of energy dissipation in reference to the frictional component. The relation between the energy dissipation by breakage is linked to the variation of the grain-size distribution curve.

General notions

The resolution of a problem using the finite element method 1 requires a rheological constitutive model developed for the continuum media. In the case of soils mechanics, these rheological constitutive models can be very complex, due to the complexity of phenomena involved. For a granular cohesionless material, the continuum media behavior is the consequence of particles mechanical behavior and interactions, i.e., the microscopic phenomena. In the framework of micro-mechanics, different techniques have been proposed to pass from an assembly of discrete particles to an equivalent continuum media (cf., section 4.3). From a thermodynamic point of view, the equivalent continuum media should represent the same characteristics of conservation and dissipation of energy as the discrete particles assembly. In order to establish a general framework of analysis let us to make an introduction to the thermodynamics of continuum media.

Thermodynamic systems

A thermodynamic system S is a part of the material universe whose energy exchange with the exterior (complementary part) is nothing but an exchange of heat and work done by volume or surface forces acting upon S . A system S can be i) isolated when it does not exchange energy or matter with the exterior, or ii) closed when it exchanges only energy, or finally iii) Open when it exchanges both energy and matter with the exterior. The state of a system is described by macroscopic quantities or state variables as volume, pressure, temperature, etc.

Thermodynamic state variables

A thermodynamic state variable α is a macroscopic quantity, which is a characteristic of the system S . It can be a scalar, a tensor, a vector. e.g., the stress tensor or temperature. A state variable is said to be extensive if it is proportional to the mass of the system, otherwise it is said to be intensive. It is accepted that a specific (i.e., per unit mass) intensive variable can be associated to an extensive variable, which allows to work only with intensive variables.

Thermodynamic state

A thermodynamic state E(S, α) of a system is the set of values of the thermodynamic state variables that characterize a system S at a given moment. If the system does not evolve with time it is said to be in thermodynamic equilibrium. On the other hand, if the system suffers a transition that makes it to evolve from one thermodynamic state to another because of a stimuli, this transition is called a thermodynamic process. A thermodynamic process is reversible if the inverse evolution of the system in time implies the reversal of the actions of the external stimuli, otherwise it is irreversible.

As mentioned by [START_REF] Germain | Cours de Mécanique des milieux continus. Tome I, Théorie Génerale[END_REF] some people say that classical thermodynamics has to be called thermostatics because it deals only with equilibrium states. Thermostatics is the science that compares systems in thermodynamic equilibrium (cf., [START_REF] Maugin | The Thermomechanics of Plasticity and Fracture[END_REF], i.e., it describes the transition from a state of equilibrium E 1 (S, α 1 ) to another state of equilibrium E 2 (S, α 2 ). Thermodynamics, in its whole sense, is the study of phenomena outside a state of equilibrium, but actually not far outside this equilibrium.

Thermostatics

As mentioned before the thermodynamic equilibrium of a system is established if the independent internal variables (α i , i = 0, 1, ..., n) and the absolute temperature (ϑ) are invariants with time. It is also possible to define variables depending of the independent internal variables, these variables are called dependent state variables or state functions. Gathered together, dependent and independent state variables, form a set of state variables defining the state of the system S . If the state is altered by the infinitesimal increments dα k , dϑ, the elementary work done on the system is of the form:

dW = A k dα k (5.1)
The coefficients A k are called the forces corresponding to α k . In general, a change in the state of the system also implies that heat is supplied or extracted. When a transformation takes place near the equilibrium, the variation of the state variables is given by its differentials, and in consequence a finite amount of work is done. The work done in a finite change of state depends not only on the initial and the final states but also on the manner in which the transition takes place.

First fundamental law

The first fundamental law of thermodynamics concerns the conservation of energy. The mechanical view of nature holds that all energy is ultimately reducible to kinetic and potential energy of interacting particles. i.e., the total energy is composed of kinetic and internal energy. In a closed system (no matter exchange) the first fundamental law states that there exists a state function E(α k , ϑ), called internal energy, such that:

dE = dQ + dW (5.2)
where dW is the infinitesimal elementary work done in reason of the evolution of variables α k , and dQ is the elementary heat supply. If a process is such that dQ = 0, i.e., if no heat is exchanged with the surroundings, the process is called adiabatic. If, on the other hand, dW = 0 i.e., no variation of internal state variables, the process is pure heating (or cooling). The first law may be considered as a balance equation. It expresses that the variation per unit of time of the internal energy is equal to the sum of the heat and work exchanged with the outside.

In general, the quantity dW is a sum of all the different forms of "work", each term being a product of an intensive variable and a differential of an extensive variable.

Second fundamental law

The second fundamental law of thermodynamics states that there exists another state function S(α k , ϑ) depending of the interval variables α k an the temperature ϑ; called entropy, this state function has the following properties [START_REF] Prigogine | Introduction à la thermodynamique de Processus irréversibles[END_REF]):

1. The entropy of a system is an extensive property. If a system is composed of different elements, then the total entropy is the sum of each element's entropy. e.g., the entropy of a mixture of gases is equal to the sum of each gas entropy.

2. The entropy's variation can be divided into two components. One the entropy flux due to the exchanges with the exterior and (entropy supply ), two the contribution of entropy due to the system's internal modifications (entropy production).

The entropy is then defined as:

ϑdS ≥ dQ (5.3)
If the previous relation holds with the equality sign, the process is referred to as reversible, otherwise as irreversible. In that sense relation 5.3 can be rewritten as:

dS = d (r) S + d (i) S
(5.4) where d (r) S = dQ ϑ is the reversible increment of S, called the entropy supply from outside, whereas

d (i) S ≥ 0 (5.5)
is the irreversible increment, referred to as the entropy production inside a system. A reversible process is characterized by the absence of entropy production; in an irreversible process the entropy production within the system is always positive.

From a thermodynamic point of view, the behavior of a continuous media can be described essentially by two functions: i) a thermodynamic potential and a dissipation function [START_REF] Germain | Cours de Mécanique des milieux continus. Tome I, Théorie Génerale[END_REF].

Thermodynamic potentials

A thermodynamic potentials is a function of state variables allowing to calculate the essential thermodynamic functions. As E is a function of state variables it is called thermodynamic potential, because from it we can deduce the variables characterizing all the thermodynamic properties of the system. For instance:

ϑ = ∂E ∂S , A k = ∂E ∂α k (5.6)
Let us take another thermodynamic potential, e.g., Z = Z(α k ) = Z(α 0 , α 1 , ..., α n ), it is possible to define an associated thermodynamic potential g(µ k ) through a partial Legendre transformation as:

g(µ k ) = Z(x k ) -µ k x k (5.7)
which leads to :

α k = - ∂g ∂µ k , and µ k = ∂Z ∂α k (5.8)
For instance, as described by [START_REF] Maugin | The Thermomechanics of Plasticity and Fracture[END_REF] pages 268-269, the Helmholtz free energy, enthalpy and Gibbs potential are associated thermodynamics potential of the internal energy.

The local accompanying state model

The approach of the local accompanying state model is based on the notion of local equilibrium (cf., [START_REF] Germain | Cours de Mécanique des milieux continus. Tome I, Théorie Génerale[END_REF][START_REF] Kondepudi | Modern Thermodynamics, From heat engines to dissipative structures[END_REF]. For a very large class of systems that are not in thermodynamic equilibrium, thermodynamic quantities such as temperature, concentration, pressure and internal energy remain well-defined concepts locally, i.e., one could meaningfully formulate a thermodynamic description of a system in which intensive variables such as temperature and pressure are well defined in each elemental volume, and extensive variables such as entropy and internal energy are replaced by their corresponding densities. A thermodynamic process then is conceived as a series of thermostatic equilibria. The response times which allow the thermostatic system to recover a new state of thermostatic equilibrium must be short. Through the use of material derivatives the first law can be written in terms of densities in the form:

ė = ϑ η + , = n β=1 τ β χβ (5.9)
where e is the specific internal energy, η is the specific entropy and τ β are state functions associated to the state variables χ β . It can be also introduced the specific Helmholtz free energy by: ψ = eϑη, then ψ = -η θ + (5.10)

Using the previous expressions together with the second law and the equations of motion a dissipation per unit volume (ϕ) can be obtained for the local form by:

ϕ = σ : D -ρ + ϑq.∇ 1 ϑ ≥ 0 (5.11)
which is known as the Clausius-Duhem inequality. In equation 5.11 σ is the Cauchy stress tensor, D is the deformation gradient, which under the hypothesis of small perturbations can be supposed as D ≈ ε, q is the heat flux vector, ∇ is the nabla operator and ρ is the mass density.

Laws of thermodynamics for continuous media

The first law will express that the total time rate of change in energy (kinetic and internal) balances the total supply of energy through external forces and heat. Which in a local form is expressed as:

ρ ė = σ : D -∇ • q + ρh (5.12)
where h is a term representing the source of heat per unit mass. The second law will state that, the total time rate of change is entropy is never less than the supply of entropy through heat. The application of this law in the form of equation 5.11 will lead to the proposition of dissipation expression like [START_REF] Maugin | The Thermomechanics of Plasticity and Fracture[END_REF]:

ϕ v = σ v : εe viscous dissipation (5.13) ϕ p = σ : εp + A • α plastic dissipation (5.14) ϕ q = ϑq • ∇ 1 ϑ thermal dissipation (5.15)
In the following of this thesis we focus on time independent, isothermal phenomena, therefore the viscous and thermal dissipations are neglected.

A trend in the study of granular materials, from a thermodynamic point of view, uses statistical mechanics and analogies to molecular thermodynamics to define thermodynamic state variables. In this field, different propositions have been made, e.g., in [START_REF] Edwards | Theory of powders[END_REF] it is proposed a "compactivity" which is a state variable analogue to the the temperature in classic thermodynamics and that is related to the density of the material. In [START_REF] Herrmann | On the thermodynamics of granular media[END_REF] the "contactropy" is defined in analogy to the entropy. The study of this trend is not the scope of this thesis. Therefore, a different approach is taken here, this approach is a more elementary one which relates the dissipation phenomena to the solid friction at the particles contacts. It is assumed that the system (REV sample) is composed of subsystems (particles) in interaction. Then the total dissipation of the system is the sum of all the individual contributions of the subsystems. Let us to review some theories having this approach for the description of mechanical behavior of granular materials.

Energy dissipation

Granular materials are mostly dissipative materials. The range of pure elasticity is restricted to a very low range of deformations. The mechanisms of dissipation in granular materials are divers. At the microscopic point of view we can point out: inter-granular friction, breakage of particles, chemical reactions, ageing. From a mesoscopic point of view (the REV scale) the dissipative mechanisms are represented by different classes of behavior, for instance plastic, thermoplastic or viscoplastic behaviors.

In this section some works on the friction dissipative mechanisms starting from micro-mechanics analysis are presented.

Rowe's dilatancy law

The works of [START_REF] Rowe | The stress-dilatancy relation for static equilibrium of an assembly of particles in contact[END_REF] established for different granular materials the strength-dilatancy characteristics using an energy approach. These works compared the analytic and experimental results for different granular assemblies. The cases under analysis concerned steel rods (Figure 5.1), uniform steel spheres in different assembly patterns, glass ballotini and sand. The analysis of the assemblies using the geometric characteristics and a given kinematic leaded Rowe to establish for regular packing (Figure 5.1) that: 1. Whatever the geometrical arrangement of the solids, the stress ratio at the peak strength and during subsequent states of deformation follows the law

σ 1 σ 3 = tan α tan(ψ + β) (5.16)
where the angle α corresponds to a geometric property of the packing, 2. The energy ratio for a fixed orientation of particle movement is given by the expression

Ė = tan(ψ + β) tan β = σ 1 σ 3 (1 + dV V ε1 )
(5.17)

where V is the initial volume of the sample and dV its variation2 .

3. Slip occurs well past the peak when the stress ratio reaches unity (σ 1 /σ 3 = 1)and the constant ambient pressure cannot be relieved.

In equations 5.16 and 5.17 β corresponds to the angle formed by the plane containing the contact points and the direction of loading, and ψ the friction angle between particles. [START_REF] Rowe | The stress-dilatancy relation for static equilibrium of an assembly of particles in contact[END_REF] proposes that by applying the principle of minimum work to a random mass of irregular particles the relocation of particles in a granular assembly would be achieved in the most economical way for a given particle arrangement, this means that the individual values of β are such that the rate of internal work done is a minimum. This may be expressed by the condition d Ė dβ = 0, which leads to find β = (45 -1 2 ψ) or :

σ 1 σ 3 (1 + dV V •dε 1 ) = tan 2 (45 + 1 2 ψ) (5.18)
During the experimental verification of equation 5.18 Rowe found that some materials were well described by this formulation while other not. In the light of experimental results [START_REF] Rowe | The stress-dilatancy relation for static equilibrium of an assembly of particles in contact[END_REF] proposes a modification to relation 5.18 as:

σ 1 σ 3 (1 + dV V •dε 1 ) = tan 2 (45 + 1 2 φ f ) (5.19)
where φ f is a friction angle which is not constant as seen in Figure 5.2. In fact, its value varies from ψ (φ µ in Rowe's notation) at minimum porosity to φ cv for loose packing, where φ cv corresponds to the friction angle at constant volume (nowadays referred as critical state). The work of [START_REF] Rowe | The stress-dilatancy relation for static equilibrium of an assembly of particles in contact[END_REF] introduced the principle of minimum energy ratio to get a stress-dilatancy relation for granular materials under axisymmetric stress conditions. According to the principle, the sliding contacts are restricted to those with the preferred angle β = (45 -1 2 ψ). This finding has been reviewed and verified by [START_REF] Horne | The behaviour of an assembly of rotund, rigid, cohesionless particles i and ii[END_REF], who extended this fact to non regular assemblies of cohesionless particles. However different criticisms have been expressed to this theory, e.g., De Josselin de Jong (1976) and [START_REF] Tokue | A consideration abou Rowe's minimum energy ratio principla and new concept of shear mechanism[END_REF] concerning the existence of the minimum energy ratio principle proposed by Rowe. At the author's belief, one point that should be pointed out concerns the existence of the β angle. More specifically, the physical possibility that the β = (45 -1 2 ψ) could stand for any irregular assembly of poly-dispersed particles. In other words, even if mathematically the solution for the minimum energy ratio was obtained for β = (45 -1 2 ψ) the conditions of the assembly could lead to a different solution. [START_REF] Frossard | L'approche energetique en mecanique des materiaux granulaires[END_REF] gives a perspective of the frictional dissipation by the introduction of a tensorial form for the analysis of the contacts. This quite general approach which includes Rowe's dilatancy formulas for particular conditions of axisymmetry stress and plane strain [START_REF] Frossard | L'approche energetique en mecanique des materiaux granulaires[END_REF][START_REF] Frossard | A micromechanical constitutive law resulting from energy dissipation by friction[END_REF]) is illustrated hereafter.

tan 2 45 + φmax 2 = σ 1 σ 3 max tan 2 45 + φr 2 = σ 1 σ 3 max -σ 3 εv ε1 Figure 5
.2: Variation of friction angles with porosity, after [START_REF] Rowe | The stress-dilatancy relation for static equilibrium of an assembly of particles in contact[END_REF]. φ max and φ r are defined in figure and φ f is given by equation 5.19.

The internal actions approach

The approach shown here concerns a summary of the works by [START_REF] Frossard | Dilatance, dissipation d'energie et critere de rupture tridimensionnel sous grandes deformations dans les materiaux granulaires[END_REF][START_REF] Frossard | L'approche energetique en mecanique des materiaux granulaires[END_REF]Frossard ( , 2005a[START_REF] Frossard | A new energy approach in granular media mechanics applications to rockfill dams[END_REF]. It starts by the contact of a couple of grains (see Figure 5.3), at this elementary contact the contact forces are represented by the vector f (a/b) acting at an angle ψ with respect to the normal vector of the contact. The relative slip velocity vector v(a/b) takes the direction of the tangential component of f (a/b), the work done is then given by the scalar product of these two vectors. Due to the geometric disposition of f (a/b) and v(a/b) the following equality (resulting from friction laws) is verified:

) ( ) f a b (a) Elementary Contact )
) )

f (a/b) • v(a/b) = sin ψ • f (a/b) • v(a/b) (5.20)
where • denotes the norm of a vector and ψ is the friction angle at the particles contact. It depends only of the mineral characteristics of particles. The last equation leads to define a symmetric work tensor given by:

p(m/n) = 1 2 { f (m/n) ⊗ v(m/n) + v(m/n) ⊗ f (m/n)} (5.21)
It can be shown that the eigenvalues of this tensor, called "tensor of elementary contact actions", verify the following properties;

p 1 + p 2 + p 3 = f (m/n) • v(m/n) (5.22a) p 2 = 0 (5.22b) |p 1 | + |p 2 | + |p 3 | = f (m/n) • v(m/n) (5.22c)
Thus the frictional dissipation equation (5.20) is expressed by:

p 1 + p 2 + p 3 = sin(ψ) • {|p 1 | + |p 2 | + |p 3 |} (5.23)
where | • | represents the absolute value. The tensor p can be viewed as a kind of "power tensor" of the equivalent system shown in Figure 5.4. Thus based on this power tensor it is possible to define the input power and the output power components as:

σ σ ε ε 1 1 σ σ
p + = 1 2   i p i + j | p j |   (5.24a) p -= 1 2   i p i - j | p j |   (5.24b)
It is worth to mention that:

p + + p -= i p i (5.25a) p + -p -= j | p j | (5.25b)
and in consequence, from equation 5.23:

p + p -= - 1 + sin ψ 1 -sin ψ = Cte (5.26)
The sum of the absolute values of the eigenvalues of the tensor p defines a norm [START_REF] Frossard | L'approche energetique en mecanique des materiaux granulaires[END_REF] that will be noted as

| • | N . Therefore, |p| N = |p 1 | + |p 2 | + |p 3 |.
In the light of this norm, equations 5.25 are rewritten as:

p + + p -= tr{p} (5.27) p + -p -= |p| N (5.28) tr{p} = sin(ψ) • |p| N (5.29)
In the last equation, the frictional dissipation at elementary contact is now expressed only as a simple condition on tensor p, instead of a mixed condition on two vectors in equation 5.20. For a discontinuous mass composed of set of particles, the total energy dissipation is the sum of the components at the elementary contacts, therefore:

m<n<N f (m/n) • v(m/n) = sin ψ • m<n<N f (m/n) • v(m/n) (5.30)
which leads to consider:

P (A) = m<n<N p(m/n) (5.31)
whence tensor P (A), called the tensor of internal actions in the granular mass has several properties as described in [START_REF] Frossard | L'approche energetique en mecanique des materiaux granulaires[END_REF], among them it has a set of eigenvectors. The set of eigenvalues of P (A) allow also to define the input and output power for the granular mass under deformation. A population effect appears naturally in the properties of this tensor P (A), due to a special type of interaction between adjacent moving contacts. This population effect has been shown to be globally defined by an "internal feedback rate" R(A). R(A) is a function of the distribution of elementary contact actions tensors:

   R(A) = 1 1-sin ψ • 1 - | A p| N A |p| N 0 ≤ R(A) ≤ 1
This internal feedback rate has also been related to the grade of disorder within the distribution of moving contacts orientations. A particular set of ordered distributions obtained for R(A) ≈ 0 corresponds to a minimal dissipation movement [START_REF] Frossard | L'approche energetique en mecanique des materiaux granulaires[END_REF][START_REF] Frossard | Localisation des déformations et chainons de contraintes en mécanique des milieux granulaires[END_REF]. With these definitions, the equation for energy dissipation resulting from the frictional contact of particles is given by:

tr{P (A)} = sin(ψ) 1 -R(A) • (1 -sin(ψ)) • |P (A)| N (5.32)
In the analysis of the passage between the discontinuous granular mass and its equivalent continuum, it has been shown that:

Equivalence énergétique avec le discontinu

( / ) u r π π ( / ) ur σ σ ( ) σ ε ε σ σ ε ε σ π π
Elementary Contact Discontinuous granular mass Equivalent continuum -REV Figure 5.5: Internal actions tensors from micro to meso-scale, after [START_REF] Frossard | L'approche energetique en mecanique des materiaux granulaires[END_REF].

• to the tensor P (A) corresponds an equivalent tensor of internal actions in the equivalent continuum π given by:

π = 1 2 σ ε + ε σ (5.33) verifying tr{π} = ij σ ij εij (5.34)
• the eigenvalues of this tensor π allow also the definition of the input and output power on the continuum equivalent media;

• when the discontinuous granular mass A is characterized by a continuum equivalent volume V (A). Then P (A) and π are related by: .35) this relation applies also to any representative elementary volume, provided that the granular material is sufficiently polydisperse with randomly irregular shaped particles as are the granular materials of general use in civil works

π(V (A)) = 1 V (A) V (A) πdV = 1 V (A) P (A) ( 5
The passage from particles contact to equivalent continuum is illustrated in Figure 5.5. The works developed hereafter focus on the verification of expressions like 5.18 and 5.32. In order to do so, it is required that other dissipation mechanisms like breakage be not present. It is though that such mechanisms are not of influence under the case of triaxial test under very low pressures.

Power in -Power out plot in triaxial test

In the subsequent of this chapter, a plot named Power in-Power out plot is introduced. This plot relates the input mechanical work with the output mechanical work, revealing the tendency of the media to dissipate or store energy. This representation is useful to identify the mechanical behavior of granular materials. To introduce this plot, let us to adopt cartesian coordinates. Then a triaxial test is seen like a cube where the stresses are applied over each face (Figure 5.6a). The shear stage of the triaxial test initiates from an isotropic state (σ 0 ). Any further increment of stress in the principal direction 1 will generate an increase of compressive strain in the same direction (both stress and strain tensors are coaxial), and in consequence an increase of the input power dP in = σ 1 dε 1 . Due to this increase of the input power the media is going to suffer extension strains in the others directions. This deformation will represent the power given back to the surroundings dP out = σ 2 dε 2 + σ 3 dε 3 . With σ 2 = σ 3 then dP out = σ 3 (dε 2 + dε 3 ). The sum of the increment of input power and the output power will represent the increment of internal work done by the REV (dW ). Hence an energy balance is represented by:

σ 1 = σ 0 σ 3 = σ 0 σ 2 = σ 0 (a) Initial condition σ 1 > σ 0 > 0 σ 3 = σ 0 σ 2 = σ 0 ε 2 2 ≤ 0 ε 3 2 ≤ 0 ε 1 > 0 (b) After deformation
σ 1 dε 1 + σ 3 (dε 2 + dε 3 ) = dW (5.36)
During the shear stage of a triaxial test the direction of dε 3 and dε 2 is opposed to σ 3 leading to a negative increment of power output. It is pointed out that hereafter the increase of output power should be read as the negative increment of the output power. In the case that the input power is equal to the output power, the increment of internal work is zero (dW = 0). It represents that the material neither stores nor dissipates energy. Else, if dW > 0 an amount of energy is stored and/or dissipated. The case when dW < 0 seem not to be feasible under the case of triaxial test. Therefore for a positive increment of internal work it is possible to establish:

σ 1 dε 1 + σ 3 (dε 2 + dε 3 ) = dW = dΨ + dΦ (5.37)
where, dΨ and dΦ are respectively the increment of free energy and dissipation. By integrating the last equation it is possible to express:

σ 1 dε 1 + σ 3 (dε 2 + dε 3 ) = dW = dΨ + dΦ (5.38)
where σ 1 dε 1 will be referred as the total input power (P in ) and σ 3 (dε 2 + dε 3 ) as the total output power (P out ). We refer then to the P in -P out plot the sequence of points having as abscissa the value of P out and as ordinate the value of P in . To establish a relation between the input and output power, let us to define the function S as:

σ 1 dε 1 -σ 3 (dε 2 + dε 3 ) = S(x i ) (5.39)
where x i represents different physical or internal parameters. According the relation between the input and output power the following cases are found (illustrated in Figure 5.7):

Case 1 dW = 0 → dW = 0. Therefore, there is no increment of internal work. The relation between the input power and the output power is a straight line crossing by the origin with a slope equal to one. Hence, S = 1.

Case 2 dW > 0 → dW > 0 ∴ dW = dΨ → dΦ = 0. The increment of internal work is stored in form of free energy, i.e., elastic energy. The relation between the input power and the output power is a monotonically increasing function. The development for the case of linear elasticity is shown in Appendix D and corresponds to:

S(E, ν, σ 0 , ε 1 ) = E -4νσ 0 ε 1 (5.40)
where E, ν, σ 0 and ε 1 are the Young's modulus, Poisson's coefficient, confining stress, and the deformation along the principal direction.

Case 3 dW > 0 → dW > 0 ∴ dW = dΦ → dΨ = 0. The increment of internal work is due only to dissipative mechanisms. Under a few simplifications it can be shown that, in the case of pure frictional dissipation (see Appendix D), the relation between the input power and the output power is a linear function crossing by the origin. Hence :

S(β, φ) = constante > 1 (5.41)
where β is a geometric argument, and φ a friction angle.

Case 4 dW > 0 → dW > 0 ∴ dW = dΨ + dΦ. The increment of internal work has two components, one due to the increase of free energy and the other due to dissipative mechanisms. Then the function S depends on the contribution of each mechanism.

dW = 0 dW = dΦ dW = dΨ σ 1 dε 1 -σ 3 (dε 2 + dε 3 )
Figure 5.7: P in -P out for cases 1 to 3 in a conventional triaxial test

P in -P out from experimental triaxial tests

The triaxial tests performed by [START_REF] Fukushima | Strenght and deformation characteristics of saturated sand at extremely low pressures[END_REF] are taken as reference in this section. Their research consisted in triaxial tests on Toyura sand with confining pressure from 2kP a to 392kP a.

These tests have been chosen because by selecting low confining pressures the effect of particles crushing phenomena is reduced. The stress-strain and volumetric response for two different confining pressures and two different voids ratio are shown in Figure 5.8. The deformation is taken negative in extension.

The material shows a behavior mainly dilatant which increases while the sample is denser. The dilatancy plot for the tests under analysis is shown in Figure 5.9. The path of the test seems to follow a "minimum energy line" [START_REF] Rowe | The stress-dilatancy relation for static equilibrium of an assembly of particles in contact[END_REF]. However the complete linearity is not evident.

The values of the Mohr-Coulomb envelope for the stress state at the end of the test are shown in Figure 5.10. Their mean value is φend = 41 o with a standard deviation of std(φ end ) = 3 o .

It should be pointed out that Mohr-Coulomb envelope's angles are not strictly those corresponding to the critical state. Because as seen in Figure 5.8 at the end of the test neither the deviatoric stress is constant nor the increment of volumetric strain is null [START_REF] Schofield | Critical State Soil Mechanics[END_REF]. Therefore the selection of these friction angles as those corresponding to the critical state could induce to an error (This is remembered because most of the constitutive models in soil mechanics are based on the concept of critical state, Hujeux (e.g., 1985); Manzari and Dafalias (e.g., 1997)). Let us now analyze these tests in the plane P in -P out . In a first time, the ratio between the increment of input and output is presented in Figure 5.11a. Whence, despite the fluctuations, it is observed how this relation lies around a mean value remaining almost constant. Regarding the total input-output power ratio in Figure 5.11b, it reveals a very good linear relationship between the total input-output power ratio. However, Figure 5.11a is not totally adequate to represent both confining pressures. Then let us divide the integral form of Figure 5.11b into Figures 5.12a and 5.12b. In Figures 5.12 [START_REF] Fukushima | Strenght and deformation characteristics of saturated sand at extremely low pressures[END_REF] expressions relating the dissipative frictional mechanism have been proposed :

S(β, ψ) = tan(β + ψ) tan(β) Rowe (5.42) S(R, ψ) = 1 + sin(ψ) -R(1 -sin(ψ)) (1 -sin(ψ))(1 -R)
Internal Actions approach (5.43) [START_REF] Fukushima | Strenght and deformation characteristics of saturated sand at extremely low pressures[END_REF] as described before the minimization of equation 5.42 leads to:

S(ψ) = tan 2 45 + ψ 2 (5.44)
Using equation 5.44 a first approximation to the inter-particles angle of friction (ψ) can be done. This angle correspond to 31.9 o and 32.1 for Figure 5.12a and 32.8 o and 30.8 for Figure 5.12b. These values lead to a mean value of ψ = 32 o with a standard deviation std(ψ) = 1 o . Which in practical terms represents the same value of ψ for all the tests. A value of ψ = 32 o seems to be slightly high for inter-particles friction. For instance [START_REF] Rowe | The stress-dilatancy relation for static equilibrium of an assembly of particles in contact[END_REF] reported a measured inter-particles friction angle of 23 o for coarse sand and 24 o for Leighton Buzzard sand. Under the assumption that equation 5.18 is valid we can examine the mobilized friction angle as a function of the dilatancy ( dεv dε 1 ), as shown in Figure 5.13. where the friction angle at peak (φ p ) and final (φ f ) states are plotted for the tests selected from [START_REF] Fukushima | Strenght and deformation characteristics of saturated sand at extremely low pressures[END_REF]. The friction angle φ ψ corresponds to the angle obtained from the expression: [START_REF] Fukushima | Strenght and deformation characteristics of saturated sand at extremely low pressures[END_REF], and friction angle from equation 5.45 with ψ = 32 o the evolution of the Mohr-Coulomb friction angle during all the tests (Figure 5.14). It is observed that the expression given by Rowe is near to the experimental results during a first stage. However, once the minimum ratio of dilatancy is reached, the values of friction angle get further from the theory of Rowe. In the light of this Rowe proposed to modify the expression as mentioned before (equation 5.19). Here let us remark that:

φ ψ = sin -1   (1 -dεv dε 1 ) tan 2 45 + ψ 2 -1 (1 -dεv dε 1 ) tan 2 45 + ψ 2 + 1   ( 5
• The expression in equation 5.18 is the consequence of minimization of equation 5.17 with respect to the angle β .

• β is an angle representing the sliding plane in a granular assembly. Then it is related to a probable kinematics of the assembly.

• All the assemblies of spheres that Rowe tested were under maximum density conditions and followed a given kinematical condition.

• The stress ratio -dilatancy plots [START_REF] Rowe | The stress-dilatancy relation for static equilibrium of an assembly of particles in contact[END_REF] and the power in-power out diagrams include not only the frictional component but an elastic component also. [START_REF] Fukushima | Strenght and deformation characteristics of saturated sand at extremely low pressures[END_REF] It is believed that the difference between Rowe's theory and the experimental results are due to an evolutive frictional dissipative mechanism. In the case of Rowe's theory it would mean that β angle changes during the test, mainly due to the evolution of the granular mass configuration. In the case of internal actions approach the internal feedback can change due to the same fact. However, let us admit that the linearity showed in the input-output power plot is representative of the frictional dissipation mechanism. This linearity is characterized by a slope S. Therefore any deviation from this tendency will represent a different mechanism of dissipation.

P in -P out in biaxial tests

In the case of biaxial test, again at low confining pressures with Toyura sand, the works of [START_REF] Tatsuoka | Strength and deformation characteristics of sand in plane strain compression at extremely low pressures[END_REF] have been analysed. The stress strain curve is shown in Figure 5.15, for a confining pressure of 0.05 kg cm 2 (5kP a). The input-output power plot is shown in Figure 5.16. An extraordinary linearity of the plot is remarked . However there seems to exist an inflection point dividing the curve into a bilinear relation. Until now there is no certainty about the signification of this inflection point, it could be due to a localization phenomena. This inflexion point is marked on the stress-strain curves also. The slope before the inflexion point corresponds to S = 3.0 and after the inflexion point to S = 3.3. Associated to these slopes a mean value of ψ = 31.6 o is obtained which is not much different from the value obtained from the triaxial test.

In the work of [START_REF] Tatsuoka | Strength and deformation characteristics of sand in plane strain compression at extremely low pressures[END_REF] a relationship between the peak stress ratio and the dilatancy shown in Figure 5.17 is found . The line corresponding to the dissipation equation in biaxial conditions given by:

σ 1 σ 3 = - dε 3 dε 1 tan 2 45 + ψ 2 (5.46) 0 2 4 6 8 10 1 2 3 4 5 6 7 ε 1 [%] σ 1 σ 3
Stress-strain curve propose a line with ψ = 33.7 o . From here, it is adopted that the frictional dissipation associated to the inter-particles friction reflects a linear relationship in the power in-power out plot. Therefore, any deviation from it will represent a different dissipative mechanism. This hypothesis underlies the fact that the friction angle ψ remains constant.

σ 1 σ 3 max =Peak Inflection Point

Energy dissipation including particles breakage

Energy balance during breakage in conventional triaxial tests

To understand the effects of grain crushing on the mechanical behaviour of granular materials let us associate the process depicted in Figure 5.18. This process underlies that particles are characterized by an elastic brittle behavior. Initially the external solicitations are transmitted into contact forces at the micro-mechanical level (Figure 5.18a). Then an arbitrary increment of stresses or strains will be transmitted as an increment of these contact forces. For the particles the change in contact forces leads to suffer a deformation (5.18b). In consequence, considering that particles behave elastically, to an increase in elastic energy. This is accompanied with an increase in dissipation, due to the frictional component of the contacts. The increment of stress or strain can continue in this way until that forces around a particle are near to its limit of strength. If forces are higher than this threshold, the particle breaks (5.18c). During breakage the elastic energy is converted into surface energy and kinetic energy.

However the system should reach a new state of static equilibrium. Therefore the kinetic energy is transformed again, a part is recovered as elastic energy while another is dissipated by friction (5.18d).

Finally a new pseudo-static state is reached, where the total surface energy has been increased as well as the total dissipation. The detailed energy balance of the processus just described is partitioned in Table 5.1. This approach is inspired from that developed for oedometric path in Frossard (2005a). The final energy balance suggests that : dW t = dU t + dP eq + + dP l + + dP l -+ dA t (5.47) where dW t is the total increment of mechanical work, dU t is the total increment of free energy associated to the elastic behaviour. dP l + and dP l -are respectively the power input and power output related to the frictional dissipation mechanism. dP eq + is the mechanical dissipation due to breakage and dA t is the total gain of surface energy.

This expression can be expressed in terms of stresses and strains for a conventional triaxial test. The component dP eq + corresponds to the increment of frictional dissipation required to equilibrate the system, initially in disequilibrium due to the creation of the new surfaces. Following Frossard (2005a), let us suppose that the strain component due to the breakage acts only on the first principal direction. This assumption is in accordance to several experimental observations were the particles Elasticity Frictional dissipation Kinetic Energy Surface Energy REV dU dΦ dK dS

Loading dU l dP l + + dP l - ≈ 0 - dW l Crushing -dU c - dK c dA c dW c Re-equilibrium
dU eq dP eq + -dK eq -dW eq Balance dU t dP eq + + dP l + + dP l - ≈ 0 dA t dW t

Table 5.1: Energy decomposition cracks were observed to be aligned mainly in the direction of the first principal stress. Under these assumptions equation 5.47 is rewritten as:

σ 1 dε 1 + σ 3 (dε 2 + dε 3 ) = σ 1 dε e 1 + σ 3 (dε e 2 + dε e 3 ) + σ 1 dε eq 1 + σ 1 dε l 1 + σ 3 (dε l 2 + dε l 3 ) + dA t (5.48)
As supposed before the component of the mechanical work along the third direction is not affected by particles' breakage. i.e., there is no deformation along the second and third direction related to particles breakage. Therefore σ 3 (dε eq 2 + dε eq 3 ) = 0. From the previous definitions and assumptions we obtain:

σ 1 dε 1 = σ 1 dε e 1 + σ 1 dε eq 1 + σ 1 dε l 1 + dA t
(5.49) using equations 5.39 and 5.40 for frictional dissipation the previous expression becomes:

σ 1 dε 1 -σ 1 dε e 1 -σ 3 (dε l 2 + dε l 3 )S(ψ) = σ 1 dε eq 1 + dA t (5.50)
Let us in a first approximation neglect the elastic component of the last equation and just to focus on the dissipative part of the increment of total work. Then equation 5.50 is rewritten as:

σ 1 dε 1 -σ 3 (dε l 2 + dε l 3 )S(ψ) = σ 1 dε eq 1 + dA t (5.51)
or integrating

σ 1 dε 1 -S(ψ) σ 3 (dε l 2 + dε l 3 ) = σ 1 dε eq 1 + dA t (5.52)
Finally, let us name the left part of the last equation the total breakage energy B E .

Then

dB E = σ 1 dε eq 1 + dA t (5.53)
It has been assumed that the components dε eq 1 and σ 1 are co-linear, therefore the increment of frictional dissipation due to particles breakage σ 1 dε eq 1 ≥ 0 is positive. Similarly, dA t ≥ 0. Under the framework established here (Table 5.1), it is coherent to think that a relation should exist between dA t and σ 1 dε eq 1 . In the following the analysis of the component σ 1 dε eq 1 is done in the input-output power plot. Afterwards, a proposition to link the dissipation due to breakage and the variation of the grain-size distribution curve is presented.

Energy dissipation due to breakage

Equation 5.52 is applied to the triaxial tests performed by [START_REF] Marachi | Strenght and deformation characteristics of rocfill materials[END_REF]. This work was a remarkable testing campaign performed at the end of the 60s during the development of a large dam construction program in California (USA). During this campaign different materials were tested (Oroville, Pyramide and Crushed Basalt) under different confining stresses (30, 140, 420 and 650 psi ) and for three different size of samples: small (S, sample diameter of 2.8 in), medium (M, sample diameter of 12 in) and large (L, sample diameter of 36 in). The analysis shown here concerns the Pyramid Dam material in size L, for which the stress-strain curves are shown in Figure 5.19. The dilatancy plot for the tests under analysis is shown in Figure 5.20. The path of the "minimum energy line" [START_REF] Rowe | The stress-dilatancy relation for static equilibrium of an assembly of particles in contact[END_REF] is not totally linear. The input-output power plot normalized by the confining pressure for Pyramid dam material is shown in Figure 5.21. The following observations are highlighted: • during the initial part of the curve a monotonically increasing function is observed (convex).

• after a given point the curve changes the tendency and passes from convex to a rather concave form.

• near the end of the curve the relation is rather linear, representing only frictional dissipation.

The previous observations can be seen clearer in Figure 5.22, where only the case for the highest confining stress is shown (blue curve). It is recalled that for materials without breakage the inputoutput power plot is a straight line crossing the origin. Taking the last portion of the experimental (blue) curve as representative of the frictional dissipation, a slope of S = 3.44 is found. Therefore, if the material was only frictional dissipative, the input-output power would be a line with slope S = 3.44, represented by the black line in The dissipation energy due to breakage ( σ 1 dε eq 1 ) can be plotted as a function of the deformation in order to observe its evolution all along the test. This evolution is shown in Figure 5.23, which presents at the same time the dissipation energy due to breakage and the deviatoric stress both normalized by their respective maximum values. It is remarked that the total breakage energy increases all along the test. For this particular test almost 90% of the total breakage energy is dissipated before the peak. From the resemblance of the curves presented in Figure 5.23, we can suggest that energy dissipation due to breakage is a positive function with decreasing slope, like an hyperbolic function. This observation agrees with the experimental observations of breakage evolution, in terms of Marsal's factor, by [START_REF] Indraratna | Modelling of particle breakage of coarse aggregates incorporating strength and dilatancy[END_REF].

max ; σ 1 dε eq 1 max( σ 1 dε eq 1 ) ε 1 [%]
σ 1 dε eq 1 max( σ 1 dε eq 1 ) q q max Figure 5.23: Evolution of breakage energy for Pyramid dam material size L. Triaxial test at a confining pressure of 650psi from [START_REF] Marachi | Strenght and deformation characteristics of rocfill materials[END_REF] The increment of dissipation energy due to breakage is plotted in Figure 5.25. It suggests that once breakage begins a peak is reached, and afterwards the increment of breakage energy has a tendency to reduce progressively. An explication to this fact can be related to the effect of particles breakage itself. The frictional dissipation component depends on the inter-particles friction angle and the arrangement of particles. It is the arrangement of particles that will be affected during particles breakage. Therefore, the increase in breakage reduces the component of the particles arrangement. In consequence it facilitates the dissipation by friction. As a consequence, once the arrangement of particles (microstructure) is changed it is easier to dissipate by friction. Therefore the increment of dissipation by breakage reduces.

The behavior of granular materials is affected by different mechanisms including elasticity, friction and breakage. It has been proposed that frictional dissipation in a granular materials depends on two principal factors. One concerning the inter-particles friction angle and the other concerning the particles packing characteristics (density, grain size distribution). It is thought that packing conditions max ; is affected highly by particles breakage. If it is so, the frictional dissipation is directly affected. i.e., the slope S changes with particles breakage. This fact can be observed, by tracing the dissipation factor S as a function of the breakage factor B defined by [START_REF] Marsal | Large scale testing of rockfill materials[END_REF].

σ 1 dε eq 1 max( σ 1 dε eq 1 ) ε 1 [%] σ 1 dε eq 1 max( σ 1 dε eq 1 ) q q max

Breakage energy and grain size distribution

It is possible to link the breakage energy with the variation of the grain size distribution curve. Hereafter, a simplified approach is shown and it is applied to the materials tested before. This link is developed under different hypotheses, mainly concerning the same particles' shape for all sizes. According to Griffiths's fracture criteria, the generation of a crack requires an amount of energy proportional to its surface to be present in the media. It is supposed that once the crack begins to propagate the entire particle breakes. Therefore the intention here it to relate the total particles surface of a media to the variation of the grain size distribution. For any initial granular assembly a finite surface is related. This surface is nothing more that the sum of the surfaces from each particle composing the assembly. Let us propose that the volume of each particle is a function of its characteristic diameter as:

v(D) = ςD 3 (5.54)
where v(D) is the volume of a particle with diameter D and ς is a constant related to the shape of the particle. In this case the surface (s(D)) of this particle is obtained by:

∂v(D) ∂D = s(D) = 3ςD 2 (5.55)
Therefore the total surface for a given group of particles of size D i is given by:

A(D i ) = N i s(D i ) = 3ςN i D 2 i (5.56)
where N i is the number of particles in the group. As mentioned in chapter 2 the retained grain size distribution (f (D)) is a statistical description of the particles sizes inside a sample. For a given diameter D i it is equal to the relative mass of the particles having this diameter in reference to the total mass of the sample, or:

f (D i ) = m i M t (5.57)
where m i is the mass of the particles having a diameter D i and M t represents the total mass of particles in the sample. The mass m i can be expressed as:

m i = N i ρv(D i ) = N i ρςD 3 i (5.58)
where ρ is the density of particles, taken as constant for all the particles of the sample. This can be done because particles are supposed to belong to the same mineral stock. Therefore a change in the retention grain size distribution is expressed by:

∆f (D i ) = ∆m i M t = ∆N i ρςD 3 i M t (5.59)
In a similar way the variation of the surface for a group of particles of size D i due the the variation in number is given by: ∆A

(D i ) = 3ς • ∆N i • D 2 i (5.60)
Combining equations 5.59 and 5.60 the variation of the retained relative mass for a diameter D i can be expressed by:

∆f (D i ) = ρD i ∆A(D i ) 3M t (5.61)
Therefore the increase in surface is given by:

∆A(D i ) = 3M t ∆f (D i ) ρD i (5.62)
and the total variation of particles surface is related to the variation of the entire grain size distribution by:

c i ∆A(D i ) = 3M t ρ c i ∆f (D i ) D i (5.63)
From the Griffith's theory the variation of fracture energy dΓ is proportional to the area created dA, therefore :

dΓ = γdA (5.64)
where γ is the specific fracture surface energy, often considered to be constant [START_REF] Ling | Fracture Mechanics,With an Introduction to Micromechanics[END_REF]. Then combining the last two equations one obtains:

Γ = γ c i ∆A(D i ) = 3M t γ ρ c i ∆f (D i ) D i (5.65)
considering that Mt ρ is the total volume of solids V s , then the factor 3γ

c i ∆f (D i ) D i
represents the total surface energy by unit volume of solids. It is expected that the breakage energy Γ be related to the total gain of surface energy dA t in equation 5.53. At the end of the previous section, the relation between dA t and the dissipation of energy due to breakage ( σ 1 dε eq 1 ) has been suggested. Therefore, it would be expected to find a relation between σ 1 dε eq 1 and 3γ c i ∆f (D i ) D i . This means that we could associate the variation of the grain size distribution curve due to crushing of particles, to the mechanical behavior of a crushable material. Therefore, let us compare the dissipation energy due to breakage obtained from the input-output plot and the factor 3

c i ∆f (D i ) D i
corresponding to the increase in surface energy calculated from the variation of the grain-size distribution curve. This comparison is shown in Figure 5.27. The parameter γ is not included in the aforementioned factor because it is not known. Values of γ between 72 and 178 J/m 2 have been reported by [START_REF] Ghaemmaghami | Large-scale testing on specific fracture energy determination of dam concrete[END_REF] for concrete materials, and between 9 and 88 J/m 2 for limestones and sandstones rocks by [START_REF] Friedman | Fracture-surface energy of rocks[END_REF].

In Figure 5.27 a linear relation has been plotted. This linear relation correspond to the best linear function that fits the data in the figure. As the value of γ is not known for these materials, we can not conclude about the equivalence between Γ and σ 1 dε eq 1 . However, as an exercise, let us suppose that Γ = σ 1 dε eq 1 . In this case a value of γ = 190J/m 2 is obtained from the best linear fit. From the values described before, it seems that this value is not nonsense for solid materials.

3 c i ∆f(D i ) D i [m -1 ] σ 1 dε eq 1 [kJ/m 3

From broken to non-broken material in conventional triaxial test

Now, an application looking for reconstruction of the stress-strain curve of non-broken material based on the materials experiencing breakage is done. The main objective is to evaluate if the approach developed before is consistent with the experimental observations. In chapter 2 several consequences were mentioned concerning the effect of particles breakage on the mechanical behavior of granular material. Mainly concerning the increase in deformability and the reduction of the maximum shear resistance. The estimation of the stress-strain curve without breakage effects is performed through the computation of a correction term, aimed at valuating the breakage process.

We start from normal relation without breakage:

σ 1 σ 3 = -(dε 2 + dε 3 ) dε 1 S(ψ) (5.66)
Then based on the assumption that breakage affects only the deformations along the first principal direction, the following relations are established:

dε 1 = dε e 1 + dε f 1 + dε b 1 (5.67) dε 2 = dε 3 (5.68) dε 3 = dε e 3 + dε f 3 (5.69)
where dε e , dε f and dε b are the strain increments due to elastic, frictional and breakage components. Therefore, by taking dε 3 and σ 1 /σ 3 from the measured values of the test and using equation 5.66 the reconstruction of a theoretical stress-strain curve without breakage is done. The result of the application of the previous procedure is illustrated in Figure 5.28. It is observed that the initial slope of the stress-strain curve increases representing a higher stiffness. Also that the volumetric strain is reduced.

Conclusions

During this chapter a thermodynamic framework has been adopted to analyze the characteristics of granular materials. Initially a plane plot named P in -P out has been proposed in order to characterize the behavior of granular materials. This plane relates the input power to the output power during a triaxial test revealing the tendency of the media to store or dissipate the energy. The form of elementary behaviors as elastic or pure frictional has been illustrated on this plane to be taken as reference for the real materials.

Afterwards, real granular materials which we consider as essentially dealing with solid friction dissipation have been analysed. They describe a completely linear relation in the input-output power plot. This agreement between the experimental analysis and the theoretical descriptions allows to define the normal trend without breakage.

In the latter part, the study of granular material experiencing breakage has been done using the P in -P out plane. The difference between the ideal frictional behavior and the real behavior is supposed to be attributed to breakage. As breakage is a non reversible phenomena, then one can determine that breakage is a dissipative mechanism. As mentioned during the last section, particles breakage affects the mechanical behavior of granular materials. By assuming that frictional dissipation is a function of the microstructure of the media, it is found that breakage affects this factor. From a mechanical point of view, breakage increases the deformability of the media and thus reduces the maximum shear resistance. By identifying the breakage component a reconstruction of a non-broken behavior can be done. This pursuit led us to verify that by canceling the component of breakage, a behavior with less deformation was obtained. However it did not allow to observe the reduction in maximum shear strength. It is believed that in order to obtain a reduction on the peak shear strength, a variation of the frictional factor (S) should happen during the test. This is actually true considering that breakage affects the microstructure and reduces the component of particles arrangement in the dissipation factor. However, by now, no technique can be proposed to identify this evolution.

From the analyses of the tests, it is thought that there exists a relationship between the three main phenomena concerned. i.e., elasticity, friction and breakage are related to the microstructure of the media. Even if the P in -P out plane can reveal the behavior tendency of the material, it does not allow to identify clearly each component or their interactions. Therefore more efforts should be made to research in this direction. In this sense, the identification of the previous components can be done experimentally by mean of mechanical tests including:

• unload path at different stages of the test, to identify the elastic component evolution.

• grain size distribution test before and after test, to identify if breakage is important.

• tests under different stress paths (triaxial, biaxial, simple shear), and different confining stresses.

• tests under different initial void ratios, to identify the importance of the initial state.

A detailed description of the behavior at low strain range would be also valuable to identify the elastic domain.

Chapter 6

Size-scale Effects Introduction

The mechanical characterization of materials is performed by testing a representative elementary volume (REV). However through experimental results it has been recognized that samples with different sizes, although geometrically similar, do not reflect the same mechanical characteristics, e.g., maximal strength. This suggests a dependency of the mechanical response on sample's size. This effect is called the size-scale effect. This dependency is rather embarrassing for mechanical analysis purposes, because mechanical analyses are based on the generalized behavior of a continuum media deduced from the behavior of a given REV (normally tested under well controlled laboratory conditions). But, since there exist mechanical variables depending on sample's size the definition of the REV is rather relative. The size-scale effect has been of main importance in structural engineering. The size-scale effect is rigorously defined through a comparison of geometrically similar structures of different sizes. According to classical failure theories, e.g., elastic analysis with limit stress or plastic analysis, the strength limit is independent of the structure size, for any geometry. In contrast failures governed by fracture mechanics exhibits a rather strong size effect as shown in Figure 6.1. Figure 6.1: Size-scale effects: (a) on the curves of nominal stress vs. relative deflection, and (b) on the strength in a bilogarithmic plot (adapted from ACI committee 446 1992), after [START_REF] Bazant | Fracture and Size effect in concrete and other quasibrittle materials[END_REF].

It is seen in Figure 6.1b that in the absence of size-scale effect the strength is rather constant (horizontal dashed line) for any size under same test conditions. In contrast, another trend is given by the linear elastic fracture mechanics (LEFM) which suggests a linearly dependency of the strength to the sample's size (in a log-log plain). The size-scale effect is present not only in limit strength but also on the stress-strain behaviour.

This chapter is devoted to the study of size scale effect in granular materials. As triaxial test under real scales can not be done the importance of this issue is capital for the granular materials and the rockfill mechanics . Then, reduced samples are used to estimate the mechanical parameters of the media. Some works by [START_REF] Fumagalli | Test on cohesionless materials for rockfill dams[END_REF] and [START_REF] Varadarajan | Testing and modelling two rockfill materials[END_REF] have resumed the different techniques to test reduced samples of rockfill-like materials. However, there exists still the need to be able to extrapolate the results obtained from small scale test to the real behavior of the material.

This chapter initiates by the descriptions of the causes of size-scale effect. Afterwards, it will be focused on the main causes applicable to granular materials. After, a consolidation of the works by Frossard (2005a[START_REF] Frossard | A new energy approach in granular media mechanics applications to rockfill dams[END_REF][START_REF] Frossard | Safety of large rockfill dams : Scale effect in rockfill shear strenght and in dams safety factors[END_REF]Frossard ( , 2010a,b) ,b) is presented as the general framework to take into account the size-scale effect in granular materials. Finally it is applied to three different materials. The materials concerned here correspond to those tested by [START_REF] Marachi | Strenght and deformation characteristics of rocfill materials[END_REF], which have already been described in section 2.5.6.

Sources of size-scale effect

Bazant and Planas (1998) number different size-scale effect sources, mainly from a structural point of view (applied to reinforced concrete structure), as:

1. Boundary layer effect, this effect is due to the fact that boundary layers of composite solids are dependent on a given geometrical characteristics, for instance a concrete layer adjacent to the walls has inevitably a smaller relative content of large aggregates, and a larger relative content of cement and mortar. Therefore the thickness of this surface layer is rather equal for different samples sizes, and it will affect differently a small or large sample because of its relative size inside the sample.

2. Diffusion phenomena, such as heat conduction of pore water transfer or hydration heat. Their size effect is due to the fact that the diffusion half-times are proportional to the square of the size of the structure. At the same time, the diffusion process changes the material properties and produces residual stresses which in turn produce inelastic strains and cracking. The extend and density of cracking may be rather different in small and large members, thus engendering a different response.

3. Statistical size effect, which is caused for the randomness of material strength and has traditionally been believed to explain most size effects in concrete structures. The theory of this size effect, originated by [START_REF] Weibull | A statistical distribution of wide applicability[END_REF], is based on the model of a chain. The failure load of a chain is determined by the minimum value of the strength of the links composing the chain, and the statistical size effects is due to the fact that the longer the chain, the smaller is the strength value that is likely to be encountered in the chain. This explanation is rather applicable to materials that fail at the initiation of a macroscopic crack, but not to materials failing after a large stable growth of cracking zones, as reinforced concrete.

4. Fracture mechanics size effect, due to the release of stored energy of the structure into the fracture front. Because the increase of crack's size (new surface) is proportional to the energy available in the media. According to [START_REF] Bazant | Fracture and Size effect in concrete and other quasibrittle materials[END_REF] it is the most important source of size-scale effect for concrete structures.

Therefore, despite the heterogeneity factor represented by the boundary layer effect, the rest of causes of size-scale effect are related to the creation and spreading of internal cracks. These phenomena are then determinant factors causing size-scale effects. As granular materials at high stress levels can experience particles breakage, it is very likely that such an effect be present. In cohesionless granular materials the cracks are developed and spread inside the particles, then the maximum crack length is restricted by the particle's size where it develops, leading clearly to particles splitting.

Size-scale effect in granular materials

It seems that the size-scale effect has been of much more interest in concrete technology, however recent problems experienced by high CFRDs have arisen its interest on the domain of granular materials. Size-scale effect have been studied in granular materials by Frossard (2005a[START_REF] Frossard | Du sol à l'ouvrage: une vision multi-échelles de la géomécanique[END_REF]. The importance of size-scale effect concerns two principal consequences, on one hand the reduction of strength and on the other hand the increase in deformability.

A granular material can be composed of many particles of different sizes. If the particles are all issued from the same mineral deposit, they should have similar physical characteristics (mineral composition, density, etc.) including cracks density. Supposing that the original material has a constant cracks density, then the total number of cracks inside a particle is proportional to its volume, and the Statistical size effect as described before can be applicable.

Size-scale effect in individual particles

It is possible to use the statistical size effect enounced in the previous section to perform an analysis of the size-scale effect. Let us suppose that two different populations of particles issued from the same mineral stock are tested. The difference between the two populations concerns only the particles diameter. Let us name the populations p a and p b , with corresponding characteristic diameters d a and d b . Therefore the survival probability given by Weibull probability for population a is given by:

P (d a ) = exp - d a d 0(a) 3 σ a σ 0(a) m (a) (6.1)
where σ 0(a) , d 0(a) and m (a) are parameters of the Weibull probability law for the population with characteristic diameter d a and σ a is the applied stress, consequently for population b the survival probability is given by:

P (d b ) = exp - d b d 0(b) 3 σ b σ 0(b) m (b) (6.2)
However, as both populations are issued from the same mineral stock they contain the same density of internal cracks and the parameters of Weibull distribution should verify:

d 0(a) = d 0(b) = d 0 ; σ 0(a) = σ 0(b) = σ 0 ; m (a) = m (b) = m (6.3)
Let us suppose that we would like to know the relation between the stress state σ a and σ b for a given but equal survival probability, i.e., P (d a ) = P (d b ), then by using equations 6.1, 6.2 and 6.3 we obtain.

d 3 a σ m a = d 3 b σ m b (6.4)
Therefore given a stress state σ b corresponding to a diameter d b acting in a population b and causing a given fracture probability. The same fracture probability will be obtained in a population a with diameter d a under a stress state such that:

σ a = σ b d b d a 3 m (6.5)
This shows the existence of a scale effect due to the particles sizes. From an experimental point of view, in section 4.6.1 it was shown how a Weibull distribution can represent the fracture probability for a sand tested by [START_REF] Nakata | A probabilistic approach to sand particle crushing in the triaxial test[END_REF]. In rockfill-like materials the works of [START_REF] Marsal | Embankment-Dam Engineering[END_REF] concluded that the force required to break the particles was proportional to a power of the diameter. The expression proposed in the last reference was of the type: fbr = η dλ (6.6)

where fbr and d are the mean crushing force and diameter respectively and η and λ are parameters of the material.

Size-scale effect in granular assemblies

The previous section showed the applicability of the Weibull's fracture probability for the assessment of the scale effect for individual grains. Hereafter the analytical extension is carried out to the applicability at the REV's scale. This procedure underlies several hypothesis as described by Frossard (2010a), and concerning mainly:

• The mineral stock of the particles is homogeneous.

• The physical friction angle at the contact of particles does not depend on particles sizes.

• Mode I of particle's rupture is predominant.

• Failure probability of particles can be represented by a Weibull probability density function.

• There exists a geometric similarity of the media under comparison.

The geometric similarity of two samples implies that both samples have the same pattern of particles' distribution in space. The similarity condition requires, in practice, to obtain the same void ratio with parallel grain size distributions and similar shape of particles. For two similar media a and b (with d a = d b ) the particles inside them will experience the same relative 1 displacement. Therefore, they will deform in a similar way under equivalent stress conditions. However, particles' breakage is one important mechanism of deformation. Therefore, in order to obtain the same strain state in both samples, it is required that the relative amount of breakage be the same. As particle's limit strength is size-dependent, the same relative amount of breakage is obtained at different stress states for samples a and b. Therefore, the main purpose of the size-scale analysis is to find the relation between those equivalent stress states.

Let us start by considering two similar media with parallel grain size distributions as shown in Figure 6.2 and equal void ratios. Then if a strictly geometrical similarity is accomplished for the Parallel Gradations

log d % P assing d M (b) d M (a) d m (b) d m (a) (b) (a)
Figure 6.2: Parallel Gradations media, the similarity should be also verified at the particles level, as schematically shown in Figure 6.3. This similarity implies that the displacement path for the particles will be similar in both media. In a granular media, as seen in section 4.3.1, the macroscopic stress tensor can be related to the forces at the level of particles and the vectors connecting the centers of the particles. Therefore the stress tensor for both media a and b can be expressed by:

σ a = 1 V a f a ⊗ l a σ b = 1 V b f b ⊗ l b (6.7)
where σ x is the stress tensor of the REV for media x, V x its volume, and f x is the force developed in a contact of the media x. As mentioned before the crushing force is related with particle's diameter through equation 6.8. However this equation is dimensional cumbersome. Then let us propose a variation using normalized parameters as:

f br = η d d 0 λ (6.8)
where f br corresponds to a vector force causing breakage, η is also a vector force, related to a reference diameter d 0 and the factor d d 0 λ stands for a scale factor of scalar nature. The following relations can be easily established between media a and b:

V b = d b d a 3 V a
Volumes relation (6.9)

d b = d b d a d a linear dimensions relation (6.10)
Then combining the previous equations, the stress tensor for both media can be rewritten as:

σ a = 1 V a d a d 0 λ η ⊗ l a σ b = d a d b 2-λ 1 V a d a d 0 λ η ⊗ l a (6.11)
By comparing the previous two equations, it is seen that the relation between both stress tensors is given by:

σ b = σ a d a d b 2-λ (6.12)
If a given homogeneous function of order one depends on a given stress tensor, then:

f (σ b ) = f σ a d a d b 2-λ = d a d b 2-λ f (σ a ) (6.13)
By a dimensional analysis it is possible to show that the λ power is related to Weibull's power m, by (Frossard, 2010a):

λ = 2 - 3 m (6.14)
which converts equation 6.12 into:

σ b = σ a d a d b 3 m (6.15)
Therefore, in order to estimate the size-scale effect for a material the Weibull's parameters should be estimated. The importance of equation 6.15 remains in its applicability to rockfill materials, where the response of a real material can be closely approached by testing a sample of materials in a reduced scale. In practice, the testing stress of the reduced scale sample shall be often much larger than the stress expected in the field for the full scale material.

Verification of size-scale effect

As mentioned before, the stress-strain characteristics of rockfill-like materials based on the size-scale effect is of upmost importance. Its application concerns the deformability analysis of geotechnical structures where the size-scale effect is likely to appear. At the author's knowledge only two experimental campaigns have been performed looking to understand different phenomena affecting rockfill-like 130 Table 6.1: Description of triaxial test for Pyramid dam material, after [START_REF] Marachi | Strenght and deformation characteristics of rocfill materials[END_REF] Sample D

(1) max Ø Sample (2) σ

(3) c

Small (S) 0.45 2.8 30, 140, 420 & 650 Medium (M) 2.0 12 30, 140, 420 & 650 Large (L) 6.0 36 30, 140, 420 & 650 (1) Maximum particle's size in inches (2) Sample's Diameter in inches (3) Confining pressures in psi, 1 psi≈689kPa materials, including the effect of particles' size. The first experimental campaign has been reported in [START_REF] Marachi | Strenght and deformation characteristics of rocfill materials[END_REF], the second one is still under execution and corresponds to the research project "ECHO" sponsored by the French National Research Agency (ANR) and coordinated by Tractebel Engineering France. The results concerning the second experimental campaign are still not available for public release. Therefore, the present study is focused on the first experimental campaign, where three different natures of rockfill-like materials have been tested under cylindrical triaxial conditions. For each material three sets of grain-size distributions were tested, all three being parallel. The initial void ratio was nearly the same for the different samples. Now, two main applications are shown. The first concerning the application over the shear strength envelope and the second, on the strain-stress curves. The former is an essential factor for stability analysis, and the interest on the latter application is upmost regarding the parameters identification for numerical modeling. It is worth to review the test's conditions for one of the materials where this verification of the size-scale effect is going to be applied. The material shown here corresponds to the Pyramid dam material. This material has been tested in three different sizes (cf., Table 6.1) under four confining pressures.

The objective is to estimate the mechanical behaviour for a material of size L based on the tests performed in sizes S and/or M . Afterwards to compare the estimated curves against the real experimental data for material L. From the comparison between estimated and experimental curves of size L the effectiveness of the method can be judged. In order to do that, the first step concerns the estimation of the Weibull's modulus m for the Pyramid dam material.

Estimation of Weibull's parameters

The Weibull's parameters are evaluated from the statistical analysis of strength data, for instance testing a representative population of particles under plate test. However for the cases presented here this information is not available. Nevertheless, there exist data concerning the breakage characteristics for the entire sample, this information can be used to estimate the Weibull's parameters using the method developed in the last section (this procedure has been developpen in : Frossard, 2010a). It is supposed that two samples geometrically similar will experience the same amount of breakage for two confining stresses related by equation 6.4. For instance, for Pyramid Dam material the amount of breakage given by Marsal's factor (cf., section 2.2.1), is plotted as a function of the confining stress for three parallel gradations, as shown in Figure 6.4. This procedure is illustrated schematically by the red dotted line in Figure 6.4. Using equation 6.4 to find Weibull factor m different combinations can be obtained among the three sets of data available. A mean factor around 8.5 is found. In order to verify the application of the previous procedure a prediction of the breakage curve for specimen diameter of 36in (dashed red line in Figure 6.5) is obtained based on the curve of specimen diameter 2.8in (continuous blue line in Figure 6.5). The estimated breakage curve seems to approach very well the experimental data for specimen diameter of 36in. It seems that the Weibull factor obtained is representative of the material under study. In a similar way the Weibull factor m for crushed basalt and Oroville dam material was found to be m = 9 and m = 7 respectively.

Size-scale effect on Shear-strength envelope

The shear-strength envelope represents the limit condition of stability. It constitutes the first factor to verify in every geotechnical problem. Therefore, its reliable estimation is crucial for a safe design. In the case of shear-strength envelope, represented for instance in the pq plane, both components are modified by the size-scale effect. For the particular case of conventional triaxial tests, p and q are linear functions of the stresses. e.g., q = σ 1σ 3 and p = (σ 1 + 2σ 3 )/3. Hence, equation 6.13 can be applied. Therefore, having the shear strength envelope for material "a" described by the couples of points (p a , q a ) the curve corresponding to the shear stress envelope of material "b" is obtained by:

(p b , q b ) = p a • (d a /d b ) 3 m , q a • (d a /d b ) 3 m (6.16)
For the present analysis, the material b corresponds to the material tested for large particles L. Whilst, material a will be S or M material. The application over the shear strength envelope is shown in Figure 6.6 for both cases. The results show a very good agreement between the experimental and the estimated curves for both cases. This is also observed for the materials of crushed basalt and Oroville dam material presented in Appendix E.

Size-scale effect on stress-strain response

Once the safety of the geotechnical structure is verified, the conditions of service are studied in order to determine if the structure can accomplish its function under normal solicitations. For instance, a structure can be safe, but experience inadmissible deformations, rendering it not useful. In geotechnical engineering this conditions of service are analysed using numerical tools, as finite element method. This numerical tools requires the definition of different mechanical parameters of a constitutive model representing the mechanical behaviour of materials. Hence, the definition of a reliable constitutive model and a representative set of mechanical parameters are the keys for a good estimation of the conditions of service in a structure. By reliable constitutive model, we mean a constitutive model able to represent the different mechanical characteristics of a given material (hardening, softening, time effect, etc.) in a particular problem. This requires evidently the estimation of mechanical parameters which represent the stress-strain characteristics. These mechanical parameters are obtained from mechanical test on representative elementary volumes under laboratory (or field) conditions. In the case of analogue materials and experiments with size-scaled effects these parameters should be corrected to represent the real in-situ material. In the next paragraphs we will explain how to treat such effects in the stress-strain response of a material tested in triaxial conditions. The objective is to estimate the stress-strain response of the largest sample L (in consequence largest particles size) based on the stress-strain curves of a smaller sample (S or M ).

Explicative example

To make clear the procedure, let us introduce the complete application of the size-scale effect to the stress-strain curve through an explicative example. Let us suppose a pair of given materials a and b fulfilling the conditions in section 6.2.2, with particle's diameter d a < d b . From laboratory tests the stress-strain response of a is obtained for two different confining stresses (σ a(l) for low confining stress and σ a(h) for a high confining stress). Then, each pair of points (ε a , σ a ) can be converted into a pair of points representing the behavior of material b (ε b , σ b ). The similarity conditions described in section 6.2.2 leads to ε a = ε b , and the problem is reduced to find the equivalent stress state in b holding for this strain condition. This stress state in b is found by applying the size-scale factor to the stress state of material a. Therefore, the strain characteristics of a material a tested at a confining pressure σ Moreover, let us suppose that the strain-stress response of media b is required for a given confining stress σ c(x) . Then in practice two cases are found: • when σ c(x) > σ b(h) , then σ c(x) is estimated by extrapolation

• when σ b(l) ≤ σ c(x) ≤ σ b(h) , then σ c(x)

Pyramid dam material

The procedure described before is applied on the Pyramid dam material. For which the Weibull factor m was found to be m ≈ 8.5 (see section 6.3.1). The application to the data from sample S is shown in Figure 6.8. For each confining pressure three curves are presented. The first corresponding to the experimental data for the material size L (continuous blue line). The second corresponding to the material size S, for the same confining stress (dashed red line), and finally the third corresponding to

ε 1 q σ a(h) σ a(l) σ b(h) = σ a(h) da d b 3 m σ b(l) = σ a(l) da d b 3 m
Figure 6.7: Explicative example of the size-scale effect application the estimation of the curve for a material of size L obtained from the data of size S (black lines). i.e., the application of the size-scale factor. For this last case two particular cases are differentiated, one when the curve is obtained from interpolation (continuous black line) and the other when the data are obtained from extrapolation (dashed black line). It is observed in Figure 6.8 that the estimated curves agree very well with the experimental curves for the interpolated cases. For the extrapolation, the exact shapes do not agree completely, but the approximation is considered to be good enough. The verification is also performed using material M as initial data. The results are shown in Figure 6.9. In this case, the extrapolation curves are closer to the experimental results both in shape and in values. The volumetric strain is also deduced by interpolation or extrapolation, using the same rule that for stresses. This procedure leads to the Figures 6.10 and 6.11. The agreement between experimental data and the correction by the scale factor is very satisfactory.

The same procedure has been applied to the crushed basalt and Oroville dam materials. The results are presented in Appendix E, and are also very satisfactory.

Conclusions

Particles breakage induce size-scale effects in granular materials. A methodology based on the physics of grain breakage, to extrapolate the mechanical characteristics of reduced gradations to approximate that of real gradations has been shown. This work is useful in the two main geotechnical analyses. On one hand, the limit equilibrium to check the stability of the structure, and on the other hand the deformability characteristics to check the correct state of service. However, the application of this scale correction can be done under the conditions mentioned in section 6.2.2. It requires in practice to use at least two sizes of samples to measure the scaling parameters and to extrapolate the rheological laws to the third one. For instance, triaxial tests on two materials with maximum particles size around 1.0cm and 5.0cm, will lead to samples diameters around 6-8cm and 30-40cm respectively. These diameters are still in the range of triaxial tests nowadays. The complete laboratory testing plan shall include:

• gradation tests, before and after mechanical tests,

• triaxial test at least for three confining stress levels,

• uniaxial compression of particles from different sizes, From these tests, one can analyze the change of distribution and find the Weibull parameters m. Therefore, the Weibull modulus m can be obtained. Afterwards, the results from the small size sample can be corrected and compared against those of the largest size sample. Granular materials like those used in civil works can be composed of different particles with different mineralogic origins, and therefore different mechanical and breakage characteristics. By now the work developed here does not cover these materials. Therefore, the application of the size-scale effect correction to these materials with different mineral compositions (therefore different Weibull moduli) is a perspective to study.

It is remarked that the parameter m was obtained by a back-analysis of the mechanical tests. There is no certainty that this factor is not affected by a "population effect". i.e., that the m value could be different if it is calculated from the statistical analysis of uniaxial compression of individual particles or from a back-analysis. Therefore, it is another perspective to evaluate from an experimental point of view.

Introduction

The current practice of the concrete face rockfill dams is totally empirical. For instance it is found in the guidelines (ICOLD, 1989) the following: "The design has been mainly empirical and based on experience and judgement. It departs from usual engineering practice in the sense that no sophisticated structural designs are required, but time tested solutions in prototypes have been the basis for design". Nowadays, numerical analysis is an important judgment tool to design these structures, specially when particular conditions are encountered e.g., complex geometry. The use of numerical tools as finite elements analysis joined to a high level constitutive model constitutes a powerful mean to assess reliably the behavior of geotechnical structures. However it demands the knowledge of several parameters describing the constitutive model. From a practical point of view these parameters can be assessed from laboratory tests, where boundary conditions are even imposed or well known. But, for rockfills materials it is not that easy, because the particles forming the material can reach 1m in diameter. Therefore, to elaborate a "representative elementary volume (REV)" it would lead to samples of the order of 10m diameter and 20m height (following the traditional rules of triaxial tests). However, testing such a sample is just not feasible. The largest triaxial test samples are in the order of 1m diameter and 2m height (e.g., [START_REF] Marsal | Research on the behavior of granular materials and rockfill samples[END_REF], [START_REF] Marachi | Strenght and deformation characteristics of rocfill materials[END_REF], among a few others). On the basis of these arguments the parameters assessment should be done in a different way to the traditional laboratory test.

This chapter regroups two different and independent numerical analyses of CFRDs performed during the development of this thesis. The first case concerns a numerical model in plane strain condition in order to identify the level of particles crushing in a high CFRD. The second case concerns the three dimensional case study of a well documented dam. This case study was the subject of theme B of the 10th Benchmark Workshop on Numerical Analysis of Dams held in Paris in 2009 under the auspices of the ICOLD1 . The principal aim in this case study is to define a methodology to deduce the mechanical properties of materials based on construction information, i.e., an inverse analysis. Similar approaches have been used in dam's projects e.g., [START_REF] Yu | An intelligent displacement back-analysis method for earthrockfill dams[END_REF], [START_REF] Zhou | Settlement analysis of the shuibuya concrete-face rockfill dam[END_REF] among others, and also in other geotechnical structures [START_REF] Tang | Application of nonlinear optimization technique to back analyses of deep excavation[END_REF]. However, all the inverse analyses have been performed using 2D finite element models, without take into account the 3D effect. This comparison is important because the measures obtained to calibrate the inverse analysis are affected by the real 3D configuration. In contrast, the work developed here includes the comparison between the 2D and 3D model.

All the numerical analyses were performed with GEFDyn software developed at Ecole Centrale Paris. It is worth to mention that for the construction of the numerical models, especially in the second case, some developments were performed in order to use GiD software as a pre-processing tool to generate GEFDyn input files for bi and three-dimensional models. It is hoped that this development could be useful for the daily works of modeling at Coyne et Bellier.

Stress path and crushing

As mentioned in Chapter 1 several high concrete face rockfill dams have had "structural" problems during reservoir impounding, leading mainly to important cracking of the concrete face. These problems are related to the deformations over the upstream face. Afterward in Chapter 2, it was shown, based on literature review, that a deformation mechanism is concerned by breakage of particles. Finally in Chapter 5 the identification of the breakage deformation mechanism was done for a particular experimental case showing the influence of particles breakage. It is also known that crushing of particles in a granular assembly is reflected in a variation of the grain-size distribution curve. Therefore, it would be expected that a variation of the grain-size distribution curve indicates the development of strains due to crushing of particles. The objective of this section is to use the crushing model developed in Chapter 4 gathered to a numerical finite element model in order to illustrate the level of crushing in a 200m high CFRD. The amount of variation of the grain-size distribution curve would indicate qualitatively2 the level of strains due to crushing of particles. In order to achieve the objective of this section the numerical model is presented initially. From the numerical model different points are selected to estimate the elastic energy at each calculus step. The elastic energy recovered is an input data for the grain-size distribution model described in Chapter 4. From the application of the grain-size distribution model a new grading curve is obtained for the end of construction and reservoir impounding.

Numerical model

The numerical model is a plane strain model of an hypothetical 200m high CFRD. Two stages of modeling are concerned representing respectively the construction and reservoir impounding stages. The model includes two domains, one representing the rockfill dam, and the other representing the foundation. For the rockfill dam the Hujeux's elasto-plastic multi-mechanism constitutive model has been adopted (some details are given in the next section). For the 200m deep foundation a linear elastic constitutive model has been supposed. The construction stage was simulated by the step by step activation of layers of elements, each layer having 5m depth. The reservoir impounding stage is simulated by a surface load over the upstream part of the dam and the foundation. The mesh and boundary conditions of the model are shown in Figure 7.1.

Mechanical parameters

The foundation was modeled with a linear elastic model with young's modulus E = 10GP a and Poisson's coefficient ν = 0.2. As mentioned before the rockfill dam is modeled using the Hujeux's model. This latter is a constitutive non-associated elasto-plastic model within the framework of 7.2. The set of parameters adopted produce the stress-strain and volumetric strain curves shown in Figure 7.3 in the case of modeling a drained triaxial test.

a cyc 1.0e-3 c iso 1.0e-2 c cyc 1.0e-2
The assessment of the variation of the grain-size distribution curve using the model described in Chapter 4 requires the definition of the total elastic energy in the media. In the case of Hujeux's model the elastic behavior is non linear. Therefore the increment of elastic strains are given by [START_REF] Hujeux | Une loi de comportement pour le chargement cyclique des sols. Genie Parasismique[END_REF]:

εe v = 1 K(p) ṗ (7.1) εe = 1 3G(p) q (7.2)
where ṗ and q are the increments of the isotropic and deviatoric stresses, εe v and εe are the increments of the volumetric and deviatoric elastic strains, and K and G are material parameters defined in Table [START_REF] Marachi | Strenght and deformation characteristics of rocfill materials[END_REF] 7.1. The total elastic energy is then found by:

U = t 0 p εe v dt + t 0 q εe dt (7.3)

Stress paths

Seven points inside the dam have been selected to perform the stress path analysis. The location of these points is shown in Figure 7.4. The seven points correspond to three different heights inside the dam. For each point the stress tensor has been recovered at each calculus step, during construction and reservoir impounding. Afterwards the principal stresses have been calculated and then the isotropic and deviatoric stresses have been obtained by:

p = σ I + σ II + σ III 3 (7.4) q = 3J 2 (7.5)
where

J 2 = 1 6 (σ I -σ II ) 2 + (σ I -σ III ) 2 + (σ II -σ III ) 2 (7.6)
and σ I , σ II and σ III are the first, second and third principal stresses. The stress paths for points 1 to 3 are shown in Figure 7.5a. These points are located near the bottom of the dam, at a height of around 5m above the foundation level. The stress paths show an evolution at almost constant q-p ratio during construction and reservoir impounding. It is also observed that the increment of stress due to reservoir impounding is more important at point 1 (upstream) and decreases toward the point 3 (downstream). A similar trend is observed for points 4 to 6, located at a height of around 75m and shown in Figure 7.5b.

Grain-size evolution model

The way of assessing the variation of grain-size distribution was described in Chapter 4. Therefore, two probability density functions should be defined, one concerning the fracture probability of particles and the other defining the distribution function of the elastic energy. Due to the lack of experimental information concerning the breakage probability of rockfill-like materials, let us assume that the dam is constructed with a material that can be represented by the fracture probability found in section 4.6.1 representing a silica sand material. The fracture probability is then defined by:

F (D, U e p |D 0 , U 0 ) = 1 -exp -D D 0 γ U e p U 0 mu
, with:

• U 0 = 2.11e7 J m 3 • D 0 = 0.335mm • γ = 1.4 • m u = 1.06
Moreover the distribution function of the elastic energy is a lognormal function as supposed in Chapter 4, with a variation coefficient a = 1.5. Concerning the comminution rule a value of α = 2.1 is adopted.

As seen in Chapter 1 very often the grain-size distribution of rockfill-like materials follows a Fuller distribution. In this case the Fuller distribution is characterized by a power n = 0.7 and maximum particle size D max = 0.50m. This leads to a grain-size distribution with a C u = 13. The results of the application of the grain-size evolution model to the stress paths of the points 1 to 3 are shown in figures 7.6a to 7.6c. In these figures three grain-size distribution curves are presented and corresponding to: the initial (in blue), after construction (in black) and after reservoir impounding (in red). The variation of the grain size distribution for the points near the upstream can be more important during reservoir impounding than during construction. For instance for point 1 (Figure 7.6d), the breakage factor given by Marsal's factor (cf., section 2.2.1) corresponds to B = 13% at the end of construction and B = 33% at the end of reservoir impounding, i.e., a variation ∆B = 20%.

In the case of the points near to the center of the dam, e.g., point 2, the main amount of breakage occurs during construction. For the particular case of point 2, we obtain a breakage factor corresponding to B = 42% for construction, this value increases only by ∆B = 2% after reservoir impounding.

For the downstream part of the dam, e.g., point 3, the variation of the breakage factor during reservoir impounding is negligible, ∆B ≈ 0% and the total breakage obtained is only due to construction B = 19%.

The breakage experienced during construction has not significant effect over the concrete face. Because most of the breakage happens before the concrete face construction. In the most critical case, the practice of dam's construction leads to construct the concrete face when the embankment is still under construction. e.g., in Mohale Project [START_REF] Johannesson | Lessons learned from mohale[END_REF] the construction of the face slabs was initiated when the dam reached 2/3 of the total height. Even in this case the increase of breakage due to the construction of the 1/3 of the dam left is not significant. For instance for the case in analysis the breakage factor for the stress state corresponding to 2/3 of the construction of the dam is B = 13% for point 1, B = 39.% for point 2 and B = 17% for point 3. By comparing these values with those presented previously it is seen, in the worst scenario (point 2 or 3), that the variation of breakage factor is around 3%.

In order to analyze the effect of crushing on the deformation of the upstream face, let us make one more time the hypothesis of strains due to breakage made in Chapter 5. It was proposed that the strains due to breakage follow the direction of the principal stress. Concerning the upstream face, during construction the direction of the principal stress is almost parallel to the orientation of the upstream face (Fig. 7.7a), therefore no perpendicular strains to the concrete face occur (due to breakage). During the reservoir impounding stage the direction of the principal stress on the upstream face rotates by increasing the water level. At the end of reservoir impounding the new direction of the principal stress has an important component along the perpendicular direction to the upstream face (Fig. 7.7b). Therefore it is expected that strains due to breakage are developed in this direction.

One solution to reduce the strains due to breakage during reservoir impounding could be to reduce the maximum particle size. This means to reduce the probability of fracture using the size dependency resistance of particles. For instance if instead of a maximum particle size of D max = 0.50m, the maximum particle size is restricted to D max = 0.15m much less breakage is obtained in the model as outlined by the small variation of curves in Figure 7.8.

Proposition for a new rockfill zone

A high level of deformation during reservoir impounding can be related to the particles breakage. Through a simplified model of grain-size evolution it has been shown that an important variation of the grain-size distribution can be observed during the reservoir impounding stage, especially in a zone near the upstream face.

Moreover, It has been shown that in conditions of plane strain the stress path in a dam follow mainly an almost constant qp ratio. Therefore, under this condition we can use only p as an index of the stress-state of a point in a dam. It is evidently a simplification. The usefulness of this simplification is to estimate a zone where the increase in the stress state is relatively important during reservoir impounding.

Let us define the factor of variation of stress-state as p i pc , where p c is the mean stress at the end of construction in a point, and p i is the mean stress at the end of reservoir impounding. From the analyses made before, let us define a factor p i pc ≥ 1.5 as a criteria beyond which an important crushing is expected. Therefore by plotting this factor for the end of reservoir impounding in the 200m dam, the Figure 7.9 is obtained. In this figure a "well defined" zone with p i pc ≥ 1.5 is observed. For this zone the increase of stress during reservoir impounding will generate an important amount of particles We can then propose a possible solution to reduce the level of crushing in a CFRD and consequently the deformations due to particles breakage. Based on the analysis made before, we can suggest that the inclusion of a new zone inside the zone 3B can reduce the amount of particles breakage. The zone proposed is named zone 3B B and it is shown in Figure 7.10. It is a subdivision of the already defined zone 3B. The characteristics proposed for this zone are:

Zone 3B B : zone with special gradation to reduce the risk of deformations by crushing of particles. It is composed of sand, gravel and small rockfill-sized particles. The gradation curve is parallel to the zone 3B, but with smaller maximum particle's size (D max ). It is a triangle shaped zone with the base going from the toe of the dam to a length of l 0 = 1 4 l t . The maximum particle size is proposed as D max ≤ 20cm. However, these dimensions can be modified in function of the nature of particles and the height of the dam. As D max < 40cm the gradation of the zone 3A should be adapted. An overview of the inclusion of this zone into the standard design of CFRD is shown in Figure 7.10.

The proposition of zone 3B B will generate a stiffer zone at the toe of the dam. As the particles size is reduced the amount of crushing should decrease. Leading to lower deformations than in the classical design.

Mohale case study

The dam's construction current practice makes the Concrete Faced Rockfill Dams (CF RD) one of the most popular types. This popularity is due to their relative construction simplicity, low cost of raw materials, general stability and short construction time. The construction techniques have evolved from a dumped material put in place before the 60's to a set of construction guidelines including compaction nowadays. However the current practice is still empirical and this approach has shown its restrictions. Therefore the use of numerical method is intended to give some highlights about the complex behavior of this kind of structure and the interaction of its components. In this context the theme B of the 10th Benchmark Workshop on Numerical Analysis of Dams presented an opportunity to compare different approaches of analysis of a well documented case study. The main purpose for the participants of this theme was to reproduce on a numerical model the cracking pattern observed on the concrete face of Mohale Dam. However at the end of the workshop the initial aim seems, from the author's point of view, not totally accomplished. The contributions to the benchmark presented important differences. In general all the participants used the finite element method. However, different constitutive models were employed. The mechanical parameters used by the participants were also estimated differently. The author participated to the contribution from Coyne et Bellier to the benchmark, and in the light of the results it was considered worth to define a methodology to perform this kind of analysis.

The case study refers to Mohale Dam where cracks observations were reported during reservoir im- pounding. These problems have been experienced in other cases, e.g., Campos Novos and Barragrande in Brazil [START_REF] Pritchard | Taking the empirical approach[END_REF] and Tianshengqiao-1 in China [START_REF] Ma | Key technical problems of extra-high concrete faced rockfill dam[END_REF] . The consequences among others are the increase in leakage flow, inter slabs as well as perimetral joint opening and slabs overlapping and fracture. The causes are still a research subject, and far of being simple. On one hand, the complex behavior of granular materials undergoing non linear behavior and even crushing of particles, and on the other hand their combination with the soil structure interaction (embankment -concrete slabs and embankment -abutments/foundation). From a modeling standpoint the different components imply a three dimensional multi-scale problem due to the sizes of the structures involved, as well as their interactions.

The present case study is divided into two stages, the parameters assessment by inverse analysis, and the three dimensional modeling of the dam. During the parameters assessment stage an optimization procedure is applied to the inverse analysis. The optimization is applied to an objective function. This stage has two main objectives : i) to define an efficient objective function to be used in this kind of analysis, and ii) to determine if by using this procedure we can estimate the deformations during reservoir impounding based on the data recovered during construction. Once the parameters for the numerical model have been assessed, the three-dimensional model is constructed in order to determine: i) if the results are comparable to the 2D analysis, and ii) if the level of strains obtained at the upstream face can generate face slabs cracking.

Optimization principles

The optimization method used comprehends the implementation of Sequential Quadratic Programming (SQP), which uses a combination of Quasi-Newton and line search algorithms [START_REF] Bonnans | Numerical Optimization[END_REF] to the minimization of a constrained objective function. This method is called constrained nonlinear optimization. The problem consists of minimizing a function f (x) given the following constraints :

min x f (x) such that =              c(x) ≤ 0 ceq(x) = 0 A • x ≤ b Aeq • x = beq lb ≤ x ≤ ub (7.7)
where x represents a set of variables. The objective function f (x) is minimized taking into account the constraints over x given by the nonlinear functions c(x) (inequalities) and ceq(x) (equalities), by the linear inequalities A • x ≤ b and equalities Aeq • x = beq, and restricted to the interval lb ≤ x ≤ ub. The set of x corresponds, in the present case, to the mechanical parameters of the constitutive model. The conditions are used to give physical meaning to the selection of parameters, e.g., to impose that elastic bulk modulus should be greater than zero. It is not intended to develop new trends in optimization procedures but to focus on the application of an existent procedure to a particular problem. In that sense the verification of suitability for two objective functions is studied.

Objective functions

Concerning the definition of the objective function there exist many possibilities. At the present case the objective function corresponds to a function defined by the difference between the numerical model and the field information. In the present case the following two propositions are studied:

S 1 = n i=1 (y m i -y c i ) 2 n i=1 (y m i ) 2 (7.8a) S 2 = n i=1 y m i -y c i y m i 2 (7.8b)
where y m i corresponds to the measured settlement at a point i and y c i corresponds to the calculated settlement at the same point. Expression 7.8a allows to give more importance to measures of high intensity, instead that expression 7.8b assumes that all differences have the same weight. The function S 1 or S 2 defined in equation 7.8 are implicitly functions of the mechanical parameters of the model. The inverse analysis is performed using the information during construction stages. Afterwards, the quality of the results is evaluated by comparing simulation results and field information for reservoir impounding. The field settlement measurements have been reported for 17 settlements cells located at the middle section of the dam [START_REF] Johannesson | Lessons learned from mohale[END_REF] for two main construction stages. For each construction stage the objective function (S 1 or S 2 ) is evaluated. Hence, a total value of the objective function has been defined as the sum of the objective function evaluated at both construction stages (i.e., S total 1 = S at construction stage 1 1 + S at construction stage 2 1 ).

Project description

The 145m high Mohale Dam is located in central Lesotho and it is part of the Lesotho Highlands Water Project. It rises from elevation 1940 up to 2085.5 m.a.s.l. with a parapet wall from elevation 2078 to 2085.5m.a.s.l. and with maximum water level at 2075m.a.s.l.. The upstream and downstream slopes are 1.4h-1v and the crest length reaches 600 m. The dam was constructed of basalt rockfill on a basalt foundation. The construction stage took place between early 1998 [START_REF] Gratwick | Mohale dam, lesotho[END_REF] and was completed in 2000 [START_REF] Johannesson | Lessons learned from mohale[END_REF]. The materials zoning is shown in Figure 7.11, the general characteristics of the materials are: Zone 1: Composed mainly of stone powered from crushing plant of fly ash, covered by overburden. The zones 1A and 1C correspond to impervious fill, the gradation requirement imposes that al least 30% of the material passes sieve #200. Zone 2: Composed of doleritic basalt. For zone 2A the gradation requirement restricts the D 15 3 between 0.3 and 0.7mm. For zone 2B the maximum particle size is imposed to 76mm. Zone 3: Composed of crushed doleritic basalt. The material in zone 3A is selected small from quarry run rock, the maximum particle size is restricted to 300mm. Zones 3B and 3C correspond to quarry run rockfill for which the gradation requirements restrict that no more than 10% of the material passes the sieve #200, and maximum 50% can be smaller than 25mm. Zone 3D imposes that more than 60% should have a size over 0.6m. Material of zone 3E corresponds to quarry run doleritic basalt for which the gradation requirements restrict that no more than 5% of the material passes the sieve #200, and maximum 40% can be smaller than 25mm. More specific design details are found in references [START_REF] Gratwick | Mohale dam, lesotho[END_REF] and [START_REF] Marulanda | Theme b -analysis of a concrete face rockfill dam including concrete face loading and deformation[END_REF].

Numerical model

The modeling hypotheses for the embankment can be summarized as: • There is no consideration of relative movements at the interface between the embankment and the abutments or foundation. This leads to fix displacements at foundation boundary nodes.

• The material zoning is simplified, only zone 3 is modeled as shown in Figure 7.12

• Due to the complex geometry of construction steps, this stage is also simplified and reduced to 6 main stages, all characterized by a lift height of 5 meters. As seen in Figure 7.12 the construction is not symmetrical (delimited by gray dashed lines), the first construction sequence

Constitutive model and parameters for rockfill materials

The constitutive model retained to represent the mechanical behavior of rockfill material is the Hujeux's model. This model has been already presented and commented. More details can be found in section B.1. Due to lack of laboratory or field data concerning the mechanical characteristics of the rockfill materials, a numerical back analysis has been carried out in order to assess the parameters of the constitutive model. The procedure consists of minimizing an error function who takes into account the deviations of the numerical model from the measurement data at some given points (points of settlement measures). This minimization is carried out by an iterative optimization procedure. Initially a 3D finite element model was built with data issued from previous experience in Coyne et Bellier. However this 3D model was found to be fairly expensive in computational time to be used during the back analysis. Therefore a 2D model was constructed to perform the back-analysis. This 2D model represents the middle B section of the dam. The initial values corresponding to the "start point" of the inverse analysis are shown in Figure 7.4. Two slightly different initial values for the bulk and shear moduli have been assumed for the materials of the zones 3B and 3C. This choice represents the difference of stiffness between the materials of the zones 3B and 3C due to compaction. Different cases have been studied, they differ on the number of materials in the dam and on the number of parameters took into account to optimize for each material.

Cases of analysis

As mentioned before, different cases of analysis has been defined. The different cases analyzed are resumed in The space of data where the solution is searched can be observed clearly for the first optimization case, where only two variables are selected to be optimized. The Figure 7.14 compares the search space for both estimation factors (S 1 and S 2 ), the iso-contours shown on it are more qualitative than quantitative because their calculus depends on the distribution of points available for the interpolations. It should be remarked that the magnitude of both estimators can not be compared, because they are relative to two different quantities. In Figure 7.14 it is observed that by using the S 1 estimator the algorithm searches in a wider range, instead for S 2 it seems to find a narrower region corresponding to the minimum. The solution obtained, concerning the parameters values, for both cases are relatively close.

A fact to verify in this kind of optimization procedure concerns the number of iterations (or evaluations) needed to reach a minimum, or better, the number of steps beyond which no improvement in optimization is obtained. This evolution is shown in Figure 7.15, where a minimum point is observed around the step 45 for estimator S 1 and 22 for estimator S 2 .

Results comparison

It is highlighted that the objective function (any of equations in 7.8) is estimated only over the construction stages, and with the parameters obtained a forward analysis is done for the reservoir impounding stage. This is done in order to answer the question: can we reliably predict the deformations during reservoir impounding based only on construction information?. It is clear that a back-analysis can be done including also or only the reservoir impounding stage, but in that case the inverse analysis does not allow to answer the question just proposed.

The cases in Table 7.3 are compared in Figure 7.16 for the first and second construction stages and for reservoir impounding. In a global view no significant improvement is found between the first and fourth cases. However a look closer to the reservoir impounding stage in Figures 7.17 the similarity between the model and measures also increase . In contrast for the downstream zone no significant changes are observed. It is reminded that the upstream zone is in fact the zone of interest, because of the presence of the concrete face slabs.

Finally, the results from optimization case number five show significant similarities to the measurement data. This is so, because the constitutive model is able to represent a wide variety of soils behavior, under the condition of a careful assessment of the constitutive parameters. As is the case of the fifth case.

It is important to compare the numerical results obtained against well known mechanical characteristics of materials. For instance in the present case the material can be compared against typical behavior of a dense cohesionless granular materials. In order to do that numerical simulations have been performed to simulate drained triaxial test (Figure 7.19). The stress-strain curves show a softening behavior with an important dilatancy characteristic for low confining stresses, which is coherent with the behavior of dense granular materials. To higher stresses the peak of the stress-strain curve disappears and the volumetric behavior is mainly contracting. The parameters' assessment by inverse analysis of construction deformations has given, for the present case, a good estimation concerning the reservoir impounding deformations. An objective function of the type S 2 shows a faster convergence to a minimum than a function of the type S 1 . In the same way, function S 2 gives a better result concerning prediction of reservoir impounding deformations.

The reservoir impounding displacements are reproduced rather well in comparison with the measured data.

From the parameters assessment procedure it has been remarked that the elastic modulus for the same zone (e.g., 3B) is very different for the first and the second set, even if they are supposed to be akin materials (same mineral stock, same construction procedure, same grain size distribution, etc.). This represents a higher compressibility in the material with higher stresses, which is mainly associated to the breakage of particles. However, as attended the values for zone 3B are higher than for zone 3C. It seems that particles breakage was a main issue in this project, as signaled by [START_REF] Johannesson | Lessons learned from mohale[END_REF] who showed that above a given stress, particles breakage leads to a fall in stiffness modulus.

The mechanical response of the materials shows a dilatant behaviour, and associated to it a peak in stress, followed by a strain softening tendency at relative low stresses. At high stresses (1MPa) just one of the materials presents a low dilatant behaviour and all the stress-strain curves seems to reach a plateau without a significant peak stress. The parameters defining the size and shape of the yield function, P c0 and b, are of great importance. It has been observed by the results improvement of optimization case No 5 (Table 7.3). This observation is also supported by the stress paths shown in Figure 7.7, which follow mainly a k 0 path.

The potentiality of a complex constitutive model emerges by a careful estimation of all the parameters. And it is well known that this type of models have some parameters that can not be assessed

HS-B14 HS-B15 HS-B16 HS-B17 HS-B18

Water However, a critical analysis of the results should be done, to ensure that the results obtained form the inverse analysis have a physical meaning. For instance by comparing the parameters obtained with a database of parameters adjusted or laboratory tests on real rockfill. 5.0e-2 5.0e-2 5.0e-2 5.0e-2 5.0e-2 r mob 0.9 0.9 0.9 0.9 0.9 r ela iso 1.0e-4 1.0e-4 1.0e-4 1.0e-4 1.0e-4 The results of the bi-dimensional and three-dimensional models are compared in Figures 7.21 and 7.23. In general the settlements obtained from the 3D numerical model are lower than those calculated from the 2D model and also from those measured in situ. The settlements estimated from the 3D model for the second construction stage are underestimated in average by 28% for the B section. The 3D model seems to be affected by the influence of abutments due to the cross-valley arching effect. This influence is due to the slope of the abutments. As mentioned by [START_REF] Lefebvre | Three-dimensional finite element analyses of dams[END_REF] when the slope of the abutment is higher than 3:1 the 2D models show higher stresses than the 3D models at the bottom of the valley. In the present case study the abutments slopes are in the order of 1.6:1, which represents a stiff slope and an important 3D effect.

Three dimensional analysis

B1 B2 B4 B5 B6 B7 B8 B9 B10 B11 B12 B13 -1.4 -1.2 -1 -0.8 -0.6 -0.4 -0.2 0 Settlement [m] Settlement Cell Calculated 3D Calculated 2D Measured

Settlements during reservoir impounding

It is reminded that the principal objective during this section is to identify if the level of deformations and strains during the reservoir impounding stage can lead to cracking problems at the upstream face.

Initially the settlements during reservoir impounding are compared against the measured and the 2D analysis. This comparison is shown in Figure 7.25 for 5 settlement cells located at a height of 100m.

As for previous stage the settlements obtained from the 3D model are lower than those from the 2D analysis and also from those measured. The difference between the 3D model and the measurements increases from upstream to downstream. It is likely that the cross-valley arching effect affects this stage also which can explain the reduction of accuracy in the upstream-downstream direction.

Upstream face strains

The 3D model does not include elements modeling the concrete face. Therefore, in order to study the eventual effects over the concrete face the strains at the upstream face of the dam are studied. Under the hypothesis that there is no relative movements between the concrete face and the rockfill dam, the strains at the upstream face should be near to those of the concrete face. Therefore, the distribution of strains shown in the following figures should be regarded as the strains over a plane parallel to the upstream face, but very near to it. The first representation, Figure 7.26, concerns the distribution of the first principal strain, it is remarked that positive values represent tensile strains. This figure should be complemented with Figure 7.27 which shows the direction of this principal strain. It is remarked that a high concentration of tensile strains is located at the lower boundaries of the dam. However, it is thought that a part of these strains is relieved by the presence of the perimetral joint. This proposition is supported by the fact that: i) the direction of the strain is oriented in the direction of opening of the perimetral joint, and ii) because perimetral joints are designed to follow this type of displacements. Assuming a maximum allowable tensile strain of 0.5% it is remarked that the lower part of the concrete face can suffer cracking. The dark blue zone represents compression strains, and the maximum value is located at the center of the lower third of the upstream face. At this location compression cracks has been observed in the case of Campos Novos in Brazil (cf., Section 1.2). The second principal strain (Figures 7.28 and 7.29) over the upstream face is mainly in compression and affects especially the lower third part of the upstream face. The maximum value in compression is near to 0.3%, i.e., above the limit of reversible strain in unconfined compression of the concrete. The directions of this principal strain is nearly horizontal towards the vertical midsection of the upstream face, and they get aligned parallel to the abutments slope near the left and right flanks. The results obtained here do not allow to explain the failure in compression observed in Mohale dam at the midsection of the crest. Therefore it is possible that this phenomena can be a consequence of a rather complex rockfill-concrete face interaction, or the consequence of a progressive failure initiated in the lower third of the concrete face and developed upwards up to the crest.

Case study conclusions

During this section the case study of Mohale dam has been analysed. The unknown mechanical parameters have been estimated using an inverse analysis technique applied to a bi-dimensional finite element model. This has also allowed to define an optimal objective function suited to perform this Figure 7.28: Second principal strain on the upstream face (positive values represent tensional strains) Figure 7.29: Direction of the second principal strain on the upstream face type of analysis. The results for the 2D case were judged to be accurate enough in comparison to the measurements in situ. Afterwards a three-dimensional (3D) model has been built to reproduce the behavior of the dam. However, the comparison between the field information and the results from the 3D model in terms of settlements is less accurate than those obtained by the 2D model. This is due mainly to the cross-valley arching effect which is very important in this case. According to previous studies by [START_REF] Lefebvre | Three-dimensional finite element analyses of dams[END_REF] and the experience obtained during the present case study, it is though that:

• The inverse analysis is a good approach to estimate a set of mechanical parameters based on insitu measurements. However, representative simple finite element models are required to reduce the cost of calculus time, especially when the finite element models use non-linear elasto-plastic constitutive models.

• The advantage of this procedure remains on the fact that the parameters of the constitutive model will represent accurately the behaviour of the soil (rockfill) in-situ, taking into account the particular boundary conditions of the structure.

• The differences between the 2D and 3D model are explained by the cross-valley arching effect, which does not exists in the 2D case.

Therefore, it is suggested that:

• The in-situ measurements should be corrected (magnified) in order to use a 2D finite element model with an inverse analysis. The magnification factor applied to settlements is related to the valley shape. For instance in V shape valleys with abutments slope of 1:1 the factor can be around 1.35. For the present case (with slope 1.6:1 and horizontal flat bottom) it is thought that a factor around 1.25 could improve the 3D model results.

General Conclusions

This research has been motivated by the different problems encountered in recent high concrete face rockfill dams. The main problem concerns the development of cracks in the concrete face. The cracks over the concrete face are generated by important deformations occurring at the upstream face due to the water load. Therefore, initially this research was focused on the influence of breakage as a deformation mechanism in granular materials. The literature review has shown that the particles breakage has been identified quite a long time ago as a key feature in granular materials. However, its influence on the mechanical behavior has remained as qualitative descriptions.

Particles breakage is dependent of different factors as: grain size distribution, particles shape, dependency of particles strength on size, contact conditions, hydric conditions, mineral composition, and others. Different constitutive models have been modified in order to represent the influence of breakage on the mechanical behavior of granular materials. However, the parameters included to fulfill this aim are not related to mechanical parameters of the particles or the assembly. A special interest was presented to a relatively new theory called "continuum breakage mechanics (CBM) " because it takes into account the grain-size distribution of particles and it is developed into a thermodynamic framework. Different analyses performed using the CBM have shown some applicability limitations. However, these analyses allowed to identify two principal points. First, the need to model the variation of the grain-size distribution curve associated to the increase of breakage, and second, that it was necessary to identify the relationship among elasticity, friction and breakage.

Therefore, this research was oriented toward the identification of the breakage component in granular soils and focused on triaxial test conditions. In a thermodynamical framework it is assumed that the soil sample is a system that gets mechanical energy from the surroundings (power-in) but also that it gives back a part (power-out). The difference between the power-in and the power-out will be the net work made by the soil. This work is composed of a recoverable strain energy and a dissipation energy. Under isothermal conditions, the recoverable strain energy is related to the elastic behavior of the soil while that dissipation is related to the frictional and breakage components. The frictional component depends mainly on the physical friction coefficient between particles and the microstructure.

A representation called power-in/power-out plot has been proposed. This representation is obtained by plotting the total output power against the total input power. This plot allows to identify the tendency of the material to store or dissipate the input energy. If the material is mainly dissipative the form of the power-in power-out plot is a straight line, under condition of no significant change in microstructure. In contrast, if the material tends to store the energy the form of the power-in power-out relation is an increasing function with an increasing slope. The breakage component is then observed as a rather increasing function with a decreasing slope.

It has been shown that granular materials are mainly concerned by frictional dissipation. Therefore represented by a straight line in the power-in power-out plot. The inclusion of breakage as an additional dissipative mechanism will generate a different form in the power-in power-out plot. The deviation from the pure frictional straight line is then considered to be associated to the particles breakage. The formation of cracks leading to particles breakage represents physically the formation of new surfaces. From a thermodynamical point of view the creation of solid surfaces are related to a specific surface energy and also to a dissipative phenomena. The energy required to create the new surfaces should be present in the media. In this sense a model describing the variation of the grain-size distribution curve has been proposed. This model is based on a probabilistic fracture analysis based on the elastic energy available on the media.

A proposition to link the variation of the grain-size distribution curve to the deformations due to breakage has been done. However, this subject is left rather as a perspective for future works.

The size dependency of particles strength leads to a similar size dependency at the continuum scale. An already proposed size-scale law has been proven to be applicable to the strain-stress curves for three different rockfill materials. This fact is very useful specially concerning the assessment of mechanical parameters for rockfill materials based on reduced grain-size distributions and laboratory tests.

Another technique has been also proven to be useful for the estimation of mechanical parameters in rockfill materials. This technique concerns the inverse analysis using finite element analysis. The use of this technique with 2D finite element models led to find good agreements between model and in-situ measurements for a set of parameters. However, the application of these parameters to 3D models was less representative than for the 2D model. The reason of this difference is that 3D conditions modify the distribution of stresses and displacements with respect to a 2D condition. Therefore, to perform an inverse analysis using a 2D model, it is proposed that measured displacements should be multiplied by a factor issued from an initial comparison between the 3D and the 2D model.

From these research works, a variation to the classical design guidelines has been proposed. This variation might reduce the deformations on the upstream face due to particles breakage, and therefore to reduce the cracks developed on the upstream concrete face. > 0:extension, dεv dε 1 < 0:compression ; (2) : Compression strength of particles (3) :Initial Porosity (4) : Saturation condition (5) : Solids specific weight 1. [START_REF] Indraratna | Large-scale triaxial testing of grey wacke rockfill[END_REF] 2. [START_REF] Charles | The influence of confining pressure on the shear strength of compacted rockfill[END_REF] 3. [START_REF] Marachi | Strenght and deformation characteristics of rocfill materials[END_REF] 4. [START_REF] Marsal | Research on the behavior of granular materials and rockfill samples[END_REF] 5. [START_REF] Varadarajan | Testing and modelling two rockfill materials[END_REF] 6. [START_REF] Chávez | Estudio del Comportamiento triaxial de materiales granulares de tamaño medio con énfasis en la influencia de la succión[END_REF] 7. [START_REF] Marsal | Large scale testing of rockfill materials[END_REF] Adopting the soil mechanics sign convention (compression positive), the deviatoric primary yield surface of the k plane is given by:

f k (σ , ε p v , r k ) = q k -sin φ pp • p k • F k • r k (B.6)
with:

F k = 1 -b ln p k p c (B.7) p c = p co exp(β ε p v ) (B.8)
where, σ is the effective stress tensor and φ pp is the friction angle at the critical state. The parameter b controls the form of the yield surface in the (p , q) plane and varies from b = 0 to 1 passing from a Coulomb type surface to a Cam-Clay type one [START_REF] Schofield | Critical State Soil Mechanics[END_REF]) (figure B.1). β is the plasticity compression modulus and p co represents the critical state stress corresponding to the initial voids ratio. The internal variable r k , called degree of mobilized friction, is associated with the plastic deviatoric strain. This variable introduces the effect of shear hardening of the soil and permits the decomposition of the behaviour domain into pseudo-elastic, hysteretic and mobilized domains. It is given by: a 1 , a 2 and m are model parameters and r hys and r mob designate the extend of the domain where hysteresis degradation occurs. When the plastic strains grow dramatically in the soil, the function r k reaches its maximal value asymptotically: lim ε p →+∞ r k = 1 or lim

r k = r el k + |ε p -ε p o | a + |ε p -ε p o | (B.
ε p →+∞ q k = sin φ pp • p k • F k (B.12)
The isotropic yield surface is assumed to be : .13) with :

f iso = |p | -d p c r iso (B
r iso = r ela iso + t 0 |( εp v ) iso | dt c pc p ref + t 0 |( εp v ) iso | dt (B.14)
where d is a model parameter representing the distance between the isotropic consolidation line and the critical state line in the (ε vln p ) plane (i.e. d = |p iso |/p c ) and c controls the volumetric hardening.

In the model, an associated flow rule in the deviatoric plane (k) is assumed, and the Roscoe's dilatancy law [START_REF] Roscoe | On the yielding of soils[END_REF] is used to obtain the increment of the volumetric plastic strain of each deviatoric mechanism so that :

εp vk = λp k • α ψ • α k (r k ) sin ψ - q k p k (B.15)
where ψ is the characteristic angle [START_REF] Luong | Phénomènes cycliques dans les sols pulvérulents[END_REF] defining the limit between dilatancy and contractance of the material (figure B.2) and α ψ a constant parameter. The parameters of the elastoplastic model concern both the elastic and plastic behaviour of the soil and they are separated into two categories: those that can be directly measured from either in-situ or laboratory test results and those which, cannot be directly measured (Table B.1).

B.2 Non-linear Model of Duncan and Chang

The nonlinear behavior of soils can be explained by its dependence to different factors such as: density, water content, the fabric, drainage conditions, strain conditions, duration of the loadings, stress history, confinement pressure and the shear stress. For modelling such a behavior Kondner and his co-workers ψ a 1 , a 2 , α ψ , m, c Threshold domains r ela , r hys , r mob , r ela iso proposed a relation between the deviatoric stress and the axial deformation drawn from the behavior of the soils tested with the triaxial device. This relation known by its hyperbolic form was taken again by [START_REF] Duncan | Nonlinear analysis of stress and strain in soils[END_REF], in theory it is represented by:

σ 1 -σ 3 = ε 1 a + bε 1 (B.16)
An advantage of this approach is that parameters a and b represent physical properties, a being the inverse of the initial elastic modulus and b the inverse of the asymptotic value of stress of the hyperbole. Another advantage of the equation B.16 is the facility of linearization, which means that the parameters a and b are calculated easily starting from the triaxial compression tests data. However, it is found that the asymptotic value is greater than the real final value, for this reason a factor R f is introduced, this factor is the quotient between: the difference in principal stress in the failure and the difference at the end of the test. This quotient varies between 0.75 and 1 and is independent of the confining pressure [START_REF] Duncan | Nonlinear analysis of stress and strain in soils[END_REF]. Then the expression B.16 is also expressed in the form:

σ 1 -σ 3 = ε 1 1 E i + ε 1 R f (σ 1 -σ 3 ) f (B.17)
Moreover the dependence of the elasticity modulus with the confining pressure is given by an expression of the form:

E i = K • p a σ 3 p a n (B.18)
The expressions B.17 and B.18 gather with a failure criteria as the Mohr Coulomb (equation B.19) give a simple formulation for the assessment of the non linear behavior of soils.

(σ 1σ 3 ) f = 2c cos φ + 2σ 3 sin φ 1sin φ (B.19)

The solution of a problem given by the equations B.17 to B.19 can be assessed by an incremental procedure. In this case a tangent modulus derived from equation B.17 is used. But as we can not know the initial state of strains it is suggested to express this modulus as a function of the stress state. In order to do that the expression B.17 is rewritten as: [START_REF] Weibull | A statistical distribution of wide applicability[END_REF] which leads us to find a tangent modulus of the form:

ε 1 = σ 1 -σ 3 E i 1 - R f (σ 1 -σ 3 ) (σ 1 -σ 3 ) f (B.
E t = 1 - R f (1 -sin φ)(σ 1 -σ 3 ) 2c cos φ + 2σ 3 sin φ 2 • K • p a σ 3 p a n (B.21)
This approach has been used to estimate the behavior of several types of structures [START_REF] Lade | Elastoplastic stress-strain theory for cohesionless soil[END_REF]) using the generalized Hooke's law as constitutive equations, however it presents some disadvantages. First there is no relationship between the volumetric deformation and the shear stress, which does not allow a good assessment in the dilatant behavior, then it neglects the action of the second principal stress and the dependence of the stress path, we find also that the behavior of materials exhibiting a peak is not possible to model. Relations given by equations B.18 and B.21 have been used by [START_REF] Khalid | Nonlinear analysis of concrete face rockfill dam[END_REF] for the analysis of Cethana dam in Australia, in this case the Poisson's ratio was defined by the relationship given by Kulhawy as follows:

v t = G -F log σ 3 pa    1 - d(σ 1 -σ 3 ) 1- R f (1-sin φ)(σ 1 -σ 3 ) 2c cos φ+2σ 3 sin φ 2 •K•pa σ 3 pa n    2 (B.22)
Where G, F and d are dimensionless parameters. However, the equation B.22 should be limited because at low confining stresses the value of v t may exceeds 0.5, to avoid that [START_REF] Khalid | Nonlinear analysis of concrete face rockfill dam[END_REF] bounded the range with an upper and lower limit of 0.18 and 0485 respectively. Another trend is given by the models issued from the plasticity theory, established by a failure criterion, one (or more) yield surface, a flow rule and a plastic potential. In this context there are several models describing the behavior of geomaterials, some of them, the most used are described by [START_REF] Mestat | Lois de comportement des geomateriaux et modelisation par la methode des elements finis[END_REF]. This type of modelling has the advantage of representing the behavior of materials with peak and are based on the concept of the critical state. Thereafter we will briefly describe some models, obviously the subject of this report is not to make a statement of all existing models, this may be quite long, but to describe those who are used to represent the behavior of rockfills and finally conclude about the advantages and disadvantages.

B.3 One dimensional deformation of rockfills

In the state of current knowledge, it was identified that the fracture of aggregates generates plastic deformations in rockfills. Especially if rockfills are subject to high stress levels. This fracture is increased by increasing the water content in the pores of materials, according to the findings of Duncan and Nobari ( 1972 ) from oedometric and triaxial tests and [START_REF] Oldecop | A model for rockfill compressibility[END_REF] from oedometric test with moisture content control.

The volumetric strain evolution is proposed to be linked to the propagation of cracks. This mechanism provides, quantitatively, an explanation of physical phenomena, the long term and collapse deformations in rockfills and its dependence on the stress level and water content. The model of crack propagation is based on the approach of linear elastic fracture mechanics which provides a framework for analyzing the field of stress at the vicinity of a crack. In this context the mode I (shear normal to crack tip) has been identified as the predominant mechanism in the fracture of granular materials. This mode establishes the crack spreading by tensional stresses. This tensional stresses inside the particles are generated indirectly by the compressive forces at the particles surface, as in the Brazilian test. In general a crack in a particle spreads when tension stress exceeds a certain threshold. This threshold is reached when the stress intensity factor (K) reaches itself a hardness factor of the material (K c ). The intensity factor constraints K is given by:

K = β • σ √ π • a (B.23)
where λ d 0 is the maximum clastic compressibility index and α 0 = κw+χw σ 0 -σy and χ w is the slope in the space log wε collapse determined for a level of stress.

In the model for one-dimensional behavior presented by [START_REF] Oldecop | A model for rockfill compressibility[END_REF] the parameter χ w is supposed to remain constant. But there is a question still unanswered about the validity of this hypothesis with the evolution of the load. Specially because this parameter will determine the magnitude of the plastic strains and the hardening behavior.

B.4 Extension of Barcelona's Model to Triaxial Path

In the last section we have showed a little overview of the Barcelona's model, which is presented completely in [START_REF] Oldecop | A model for rockfill compressibility[END_REF]. Now we will mention some works based on that model but extended to triaxial solicitations. The complete reference is explained in [START_REF] Chávez | Estudio del Comportamiento triaxial de materiales granulares de tamaño medio con énfasis en la influencia de la succión[END_REF] and summarized in [START_REF] Chávez | A constitutive model for crushed granular aggregates which incluyes suction effects[END_REF]. The development of this model is based in some important works, specially [START_REF] Wan | A simple constitutive model for granular soils: Modified stress dilatancy approach[END_REF] and [START_REF] Wood | Strain softening and state parameter for sand modeling[END_REF] and in a extensive experimental work. The experimental work gathers current classification tests and more particular triaxial tests with control of water content. It means with controlled suction. The confinement stress is variable between 0.1 and 0.8 MPa. As a general conclusion they establish that the increase in suction contributes significatively to keep a higher void ratio at a given stress level, and as consequence the critical state depends strongly on the level of suction (or relative humidity). The formulation of the model is made in the classical elasto-plastic theory applied to unsaturated soils. The yield surfaces are showed schematically in figure B.5 and they are establish by the following equations f q = q p η s = 0 (B.35a) The plastic deformations are linked to the breakage of particles, and in that sense the authors propose that the breakage is a function of the effective plastic work, that by definition is the difference between the total and the elastic work. The hardening behavior is represented by the proposition originally of [START_REF] Wood | Strain softening and state parameter for sand modeling[END_REF] given by the following equation.

f
η s M p = ε p q b + ε p q (B.37)
Where ε p q is the deviatoric plastic strain, b is a parameter related to the initial strength of the material and M p is given by: With M c is the slope at the critical state in the space of invariants pq, k is a constants and ψ e is the state parameter The following paragraphs summaries the works done by Kohgo et al. (2007a), who propose an elastoplastic model for the application to rockfills, this model is based on oedometric and triaxial tests conducted on two types of material (volcanic tuff and Andesite ) tested both in dry and saturated state (Kohgo et al., 2007b). The formulation is a modification of a first model developed for unsaturated soils [START_REF] Kohgo | Verification of the generalized elastoplastic model for unsaturated soils[END_REF].

M p = M c (1 -k • ψ e ) (B.38)
During the testing stage the authors considered three states of the materials in order to observe the mechanical behavior. The first state is a saturated state, the second is a unsaturated state and the third is a state where the initial dry material is soaked at a point and then the response is observed.

The model takes into account the effect of water content via the water retention curve. In order to do that some state surfaces have been defined and showed in figure B.9. The application of this model have shown satisfactory results, representing well the experimental behavior observed over all the tests. One of the results is showed in figure B.11

B.7 Fractal dimension

A fractal is an object or quantity that displays self-similarity, in a somewhat technical sense, on all scales. The object need not exhibit exactly the same structure at all scales, but the same "type" of structures must appear on all scales. A plot of the quantity on a log-log graph versus scale then gives a straight line, whose slope is said to be the fractal dimension.

  Figure 1.2: Typical section of a CFRD.

  Figure 1.3: Plinth Schema.

  Figure 1.4: Cracks between slabs 17-18, after Johannesson and Tohlang (2007).

Figure 1

 1 Figure 1.6: Schema of deformation of a CFRD in 2D under the action of water load, after Giudici et al.(2000)

  Figure 2.1: Representations of particles' size distribution

Figure 2 . 2 :

 22 Figure 2.2: Calculus of Marsal's breakage index.

  Figure 2.3: Hardin's parameters

Figure 2

 2 Figure 2.6: Variation of breakage index B with gradation, according to data in Appendix A. Greater breakage is observed for uniform gradations for triaxial drained tests.

Figure 2 . 7 :

 27 Figure 2.7: Variation of breakage index with confining stress, according to data in Appendix A

  Figure2.10: Observations from triaxial test on sands at high stresses, after[START_REF] Vesić | Behavior of granular materials under high stresses[END_REF].

Figure 2 .

 2 Figure 2.14: Model of McDowell et al.(1996)

  Figure 2.15: Gradation curves for Pyramid dam material, after[START_REF] Marachi | Strenght and deformation characteristics of rocfill materials[END_REF]

  nombre plus grand que)

Figure 2 .

 2 Figure2.17: Fractal dimension evolution for three materials, after[START_REF] Marachi | Strenght and deformation characteristics of rocfill materials[END_REF]

Figure 2 .

 2 Figure2.18: Grain size distribution curves for dense material of "La Soledad", after[START_REF] Marsal | Research on the behavior of granular materials and rockfill samples[END_REF]

  variation of fractal dimension with σ3

Figure 3

 3 Figure 3.1: Modified definition of Hardin's breakage index afterEinav (2007a)

  Figure 3.2: Results of tests on petroleum Coke, after[START_REF] Bard | Comportement des Materiaux Granulaires Secs et a Liant Hydrocarboné[END_REF]

  Figure 3.3: Application flowchart for Einav's Model

Figure 3

 3 Figure3.4: Average curve and experimental curves for coke material[START_REF] Bard | Comportement des Materiaux Granulaires Secs et a Liant Hydrocarboné[END_REF]

  Figure 3.5: Gradation curves by[START_REF] Bard | Comportement des Materiaux Granulaires Secs et a Liant Hydrocarboné[END_REF] and hypothesis of final fractal distribution for α = 2.6

  Figure 3.6: Estimation of ω parameter keeping p y0 = 36kP a

Figure 3 . 7 :

 37 Figure 3.7: Comparison between experimental data and model results.

  variation of ω with pressure.

  . Model ω variable B.F. ω-5% B.F. ω+5% (b) Simulation of coke test with ω variable.

Figure 3

 3 Figure 3.8: Assessment of variable w.

  J 2 u [%] Error [%](b) Error E of the model vs variation in J2u.

Figure

  Figure 3.9: Sensibility to J 2u

Figure 3 .

 3 Figure 3.10: Simulation results including the experimental final gradation for coke material.

  Figure 3.11: Simulation results using experimental breakage and final grain-size distribution

Figure 3

 3 Figure3.12: Simulation results using variable ω with experimental breakage and final grain-size distribution. (experimental results from[START_REF] Bard | Comportement des Materiaux Granulaires Secs et a Liant Hydrocarboné[END_REF]

Figure 3 .

 3 Figure 3.13: Simulation results for test in Silica sand, after Nakata et al. (2001a).

Figure 4

 4 Figure 4.1: Schematic rupture modes, after[START_REF] Guyon | Du sac de billes au tas de sable[END_REF].

Figure 4

 4 Figure 4.2: Crack propagation modes in linear elastic fracture mechanics after Maugin (1992).

Figure 4

 4 Figure 4.3: Brazilian test principle

Figure 4

 4 Figure 4.4: Crushing of individual particles under uniaxial compression, after Nakata et al. (1999).

  Figure 4.5 shows three samples, two of discs and the other of a sphere, and Figure 4.6 shows schematically the way of crack patterns.

  (a) Glass discs failure at 18kN. (b) Glass discs failure at 46kN. (c) Glass spheres failure at 7.3kN.

Figure 4

 4 Figure 4.5: Crack evolution in particles, after Schonert (2004).

Figure 4

 4 Figure4.6: Schematic crack pattern in a compressed sphere, after[START_REF] Schonert | Breakage of spheres and circular discs[END_REF]

Figure 4

 4 Figure 4.8: Network of normal forces in sample D, after[START_REF] Radjai | Force distributions in dense two-dimensional granular systems[END_REF]

Figure 4 .

 4 Figure 4.10: Schematic elementary particles contact c.

Figure 4 .

 4 Figure 4.11: Probabilistic fracture analysis schema.

Figure 4 .

 4 Figure 4.12: Fracture tensile stress for silica sand after Nakata et al. (2001a)

Figure 4 .

 4 Figure 4.13: Fracture probability for two sands after Nakata et al. (2001b)

Figure 4

 4 Figure 4.14: Experimental and fitted curve for silica sand (Nakata et al. (2001b))

Figure 4 .

 4 Figure 4.15: Fitted curves and parameters of breakage probability

  55 and 4.57 relate the elastic energy to the mean elastic energy for the particles in a class c, i.e., to the mean of the distribution function, which leads to E(U e p ) = exp λ + ζ 2 2= Ū e p . In order to estimate the parameters λ and ζ for the distribution of energy a relation between these two parameters is established. It is proposed that standard deviation ζ is supposed to be a fraction of the mean ζ = a • λ. i.e., a is the coefficient of variation of the associated normal distribution. Then the following procedure can be developed: possible, one negative and one positive.

Figure 4

 4 Figure4.16: Initial Gradations, after Nakata et al. (2001a)

  Figure4.17: Oedometric Test, afterNakata et al. (2001a)

Figure 4 .

 4 Figure 4.18: Stress-strain curve, afterNakata et al. (2001a)

  Silica Sand 0.25-2.0mm

Figure 4 . 20 :

 420 Figure 4.20: Experimental and numerical results comparison based on data from Nakata et al. (2001b). Points represent experimental data and continuous lines the model. 90

Figure 5

 5 Figure 5.1: Rowe's Test configuration, after Rowe (1962).

  Figure 5.3: Elementary contact after[START_REF] Frossard | L'approche energetique en mecanique des materiaux granulaires[END_REF]

Figure 5 . 4 :

 54 Figure 5.4: Pseudo continuum volume.

Figure 5 . 6 :

 56 Figure 5.6: Representation of conventional triaxial test

  Figure 5.8: Triaxial test[START_REF] Fukushima | Strenght and deformation characteristics of saturated sand at extremely low pressures[END_REF]

  Figure5.9: Dilatancy plot for Triaxial test in[START_REF] Fukushima | Strenght and deformation characteristics of saturated sand at extremely low pressures[END_REF]

Figure 5

 5 Figure 5.11: P in -P out relation for Triaxial tests in[START_REF] Fukushima | Strenght and deformation characteristics of saturated sand at extremely low pressures[END_REF]

  Figure 5.13: Peak and final friction angles for Triaxial tests in Fukushima and Tatsuoka (1984), and friction angle from equation 5.45 with ψ = 32 o

Figure 5

 5 Figure 5.14: Mobilized Mohr-Coulomb friction angle for Triaxial tests in[START_REF] Fukushima | Strenght and deformation characteristics of saturated sand at extremely low pressures[END_REF]

Figure 5

 5 Figure5.15: Biaxial test[START_REF] Tatsuoka | Strength and deformation characteristics of sand in plane strain compression at extremely low pressures[END_REF]

Figure 5

 5 Figure5.17: Biaxial test on Toyura sand, after[START_REF] Tatsuoka | Strength and deformation characteristics of sand in plane strain compression at extremely low pressures[END_REF].

  Figure 5.18: Particles crushing pattern

Figure 5 Figure 5 . 20 :

 5520 Figure 5.19: Triaxial test Pyramid Dam Material size L[START_REF] Marachi | Strenght and deformation characteristics of rocfill materials[END_REF]

FigureFigure 5 Figure 5

 55 Figure 5.21: Power in-Power out plot for Pyramid Dam Material size L Triaxial tests, after[START_REF] Marachi | Strenght and deformation characteristics of rocfill materials[END_REF]

Figure 5 .

 5 Figure 5.24 gathers all the tests for sample L of Pyramid dam material showing again a continuous increase of the breakage energy.

Figure 5 .

 5 Figure5.24: Evolution of breakage energy for Pyramid dam material size L. Triaxial test from[START_REF] Marachi | Strenght and deformation characteristics of rocfill materials[END_REF]

Figure 5 .

 5 Figure 5.25: Evolution of increment of breakage energy for Pyramid dam material size L. Triaxial test from[START_REF] Marachi | Strenght and deformation characteristics of rocfill materials[END_REF]

  Figure 5.26: Variation of S with breakage index for Pyramid dam material size L.

]Figure 5

 5 Figure 5.27: Dissipation energy due to particles breakage vs surface energy by unit volume of solids

Figure 5

 5 Figure 5.28: Reconstitution of non-broken curves from Triaxial test Pyramid Dam Material size L σ c = 4.5M P a after[START_REF] Marachi | Strenght and deformation characteristics of rocfill materials[END_REF]

Figure 6

 6 Figure 6.3: Geometrical similarity at particles level.

Figure 6

 6 Figure6.4: Breakage factor evolution as a function of the confining stress for three sample's size of Pyramid dam material, after[START_REF] Marachi | Strenght and deformation characteristics of rocfill materials[END_REF], and Weibull coefficient estimation for size-scale effect analysis.

  Figure6.5: Estimation of breakage curve for Pyramid Dam Material tested by[START_REF] Marachi | Strenght and deformation characteristics of rocfill materials[END_REF].

  Experimental d = 0.45in Experimental d = 6in Model d 0.45 → d 6 (a) From maximum particle size d = 0.45in to d = 6in. Experimental d = 2in Experimental d = 6in Model d 2 → d 6 (b) From maximum particle size d = 2in to d = 6in.

Figure 6

 6 Figure 6.6: Size-scale effect application to Shear-strength envelope of Pyramid Dam Material

  c(a) are the same that those of a material b tested at σ c(b) = σ c(a) (d a /d b ) 3 m . This procedure is schematically shown in Figure 6.7. The new set of curves for media b are suppose to have the same strain characteristics as media a. Let us clarify by a numerical application. Let us take σ c(a) = 1000kP a, and (d a /d b )3 m = 0.3. Therefore, the strain characteristics of material a will be observed in a material b tested at a confining stress of σ c(b) = 1000kP ax0.3 = 300kP a.

  is estimated by interpolation between the curves σ b(l) and σ b(h)

  Figure 6.8: Size-scale effect application to Stress-strain curve for Pyramid Dam Material, S → L

  Figure 6.9: Size-scale effect application to Stress-strain curve for Pyramid Dam Material, M → L

  Figure 6.10: Size-scale effect application to volumetric strain curve for Pyramid Dam Material, S → L

  Figure 7.1: Numerical model of a 200m high CFRD isotropic hardening plasticity using effective stress and characteristic state concept. In this model, four mechanics are defined: three plane strain deviatoric mechanisms associated to three perpendicular planes, and a mechanism for the isotropic behavior. The four mechanisms are coupled through the volumetric strain. A description of the model parameters is shown in Table 7.1 and represented in Figure 7.2, more details are found in appendix B.1.

Figure 7

 7 Figure 7.2: Hujeux's Model

Figure 7

 7 Figure7.3: Stress-strain response for Pyramid dam material size M, according to experiments by[START_REF] Marachi | Strenght and deformation characteristics of rocfill materials[END_REF]

FigureFigure 7

 7 Figure 7.4: Selected points for stress path analysis

Figure 7 . 6 :

 76 Figure7.6: Evolution of grain-size distribution curve and breakage factor as defined by[START_REF] Marsal | Large scale testing of rockfill materials[END_REF]

  (a) at the end of construction (b) at the end of reservoir impounding

Figure 7 . 7 :

 77 Figure 7.7: Directions of the principal stress

Figure 7

 7 Figure 7.8: Evolution of grain-size distribution curve

Figure 7 . 9 :

 79 Figure 7.9: Factor p i pc . Darkest red zones represents p i pc ≥ 1.5

Figure 7 .

 7 Figure 7.11: Zoning of Mohale Dam after[START_REF] Marulanda | Theme b -analysis of a concrete face rockfill dam including concrete face loading and deformation[END_REF]

Figure 7 .

 7 Figure 7.12: Material Zoning and settlement cells distribution

  (CS01) lifts up to the height of 2010 m.a.s.l., afterwards the second construction stage (CS02) located at upstream reaches 1995 m.a.s.l. followed by the CS03 at downstream to the same height. Next, the CS04 reaches 2030 m.a.s.l. and finally the CS05 and CS06 lift up to 2040 and 2078 m.a.s.l. respectively.

Figure 7 .

 7 Figure 7.13: 2D Mesh for inverse analysis of Mohale Dam

Figure 7

 7 Figure 7.14: First case search space and iso-contours for a)estimator S 1 and b) estimator S 2

Figure 7

 7 Figure 7.17: Results comparison for indicator S 1 for reservoir impounding

Figure 7 .

 7 Figure 7.18: Results comparison for indicator S 2 for reservoir impounding

AFigure 7

 7 Figure 7.19: Stress strain curves for parameters obtained in case five

Figure 7 . 20 :

 720 Figure 7.20: Finite Element mesh for 3D model

Figure 7 .

 7 Figure 7.21: B section settlement for First Reading Construction stage

Figure 7 Figure 7

 77 Figure 7.22: A and C sections settlement for first reading construction stage

Figure 7 .

 7 Figure 7.26: First principal strain on the upstream face (positive values represent tensional strains)

  Figure B.1: Influence of parameter b on the yield surface shape.

  9)ε p is the plastic shear strain accumulated during the shearing (ε p = t to εp dt); ε p o is the plastic shear strain accumulated from the very beginning of the loading until the last unloading (ε p o = to 0 εp dt).Parameter a is defined as :a = a 1 + (a 2a 1 ) α k (r k )

  Figure B.2: Critical State and Chacteristic State Lines.

  p = p(Ψ)p * = 0 for p < p y or (B.35b)f p = p[λ i + λ d (Ψ)κ]p y λ d (Ψ)p * (λ iκ) for p > p y (B.35c)We see the relation between equations B.35 and equations B.31 and B.32. Then p y is equivalent to σ y and represents the clastic yield stress. By analogy λ d (Ψ) is given by:λ d (Ψ) = λ d 0 + α Ψ ln Ψ + p atm p atm (B.36)

Figure B. 5 :

 5 Figure B.5: Yield surfaces of the model proposed by Chávez and Alonso (2003).

  1 . As M c and b are functions of the breakage of particles they are also function of the effective plastic strain, moreover M c and ψ e are function of the suction as well. Then by combining equations B.37 and B.38 we obtain:η s M c (W P E , Ψ)(1k(vv cr )) vcr is the specific volume at the critical state and is a linear funtion of the mean stress and of the suction. v cr = y(Ψ)m cr p (B.40) where y(Ψ) = y 0 + α y ln Ψ + p atm p atm (B.41) M c (W P E , Ψ) = M cr (W P E ) 1 + k s ln Ψ + p atm p atm (B.42) and M cr (W P E ) = M cres + (M c 0 -M cres ) exp(-a.been made to M c and b in Chávez (2004). To represent the isotropic behavior Chávez and Alonso (2003) propose: dp * = dε p q λ iκ (B.45)

  Figure B.9 shows the possible shapes of state surfaces passing from a saturated state to a dry state (or inverse). The states surfaces may be expressed as: e = -λ * log p + Γ * , s * is the suction, λ * is the slope of elog p curves, Γ is void ratio of elog p curves at p =unit, λ and Γ are values of λ * and Γ * at saturation, respectively, s * m is suction where the state surface has the minimum or maximum slope values,λ * = λ * m at s * = s * m ,Γ * = Γ * m at s * = s * m, and e 0 01 , e 0 02 , λ * f 1 , λ * f 2 , a * 1 and a * 2 are material parameters.f (I 1 , J 2 , θ) = α * cs I I 1 , J 2 , θ) = b * 2 (I 1 -I 0 ) 2 + a 2 J 2 g(θ) 2a 2 b * 2 θsin θ sin φ ) (B.63)Where f is the critical state line, f 1 is the failure surface, f 2 is the cap surface, I 1 is the first stress invariant, J 2 is the second invariant of deviator stress, θ is Lode angle and I 0 , a, b * , K * , α * and α * cs are defined in figure B.10. The material parameters are summarized in table B.5The plastic potential are defined byψ 1 = b 2 (I 1 -I 0 ) 2 + a * 2 J 2a * 2 b 2 = 0 (I 0 < I 1 ≤ P 2 ) (B.64) ψ 2 = b 2 (I 1 -I 0 ) 2 + a 2 J 2a 2 b 2 = 0 (I c < I 1 ≤ I 0

  Figure D.1: P in -P out For linear elasticity

  Figure E.4: Size-scale effect application to stress-strain curve of crushed basalt, from maximum particle size d = 2in to d = 6in.

  Figure E.7: Size-scale effect application to stress-strain curve of Oroville dam material, from maximum particle size d = 0.45in to d = 6in

Table 1 .

 1 2: Parameters influence on particles' breakage.

				Increase in
	Particles strength	Friction Angle	Voids Ratio	Shape d i d j c	Water Content	Gradation	Breakage	Stress State
	leads to:

Table 2

 2

		.1: Factors affecting Particles' breakage
	Particles	Assembly	Condition
	Mineral nature of particles	Gradation	Water content conditions
	Particles' shape	Compacity	Stress level
	Particles' size		Stress path
			Time
	breakage increased with the amount of the calcareous sand all other conditions kept constant. The
	works of

Table 2

 2

	.2: Constitutive Equations for Barcelona's Model
	Yield Functions

variation with σ3 Figure 2.16: Fractal dimension analysis for material size S of Pyramid dam, based on data given by

		-2.3
		-2.4
		-2.5
		-2.6
	α
	α=2.47	-2.7
	α=2.57 α=2.71	-2.8
	α=2.87	-2.9
	α=2.98
		0 -3	1000	2000	3000	4000	5000
				σ 3 [kP a]
	(c) Fractal dimension α	(d) Fractal dimension (α)

Table 4 .

 4 1: Description of DEM Models in[START_REF] Radjai | Force distributions in dense two-dimensional granular systems[END_REF].

	Model Number of Particles Particles size [mm] Particles Distribution
	A	500	3.8,7.5	uniformly
	B	1200	3.8,7.5	uniformly
	C	4025	1.5,7.5	uniformly
	D	1024	0.65,1.05,1.6	512,320,192

  The work of[START_REF] Mcdowell | The fractal crushing of granular materials[END_REF] proposed a fractal pattern, indeed they proposed that at every time the distribution of particles size was characterized by a fractal dimension. This proposition was an emergent characteristics of a defined self similar pattern of breakage for individual grains.McDowell et al. mentioned that this was in accordance with previous observations by Turcotte, Epstein, Sammis et al. (cited by McDowell et al.,

Table 4

 4

	.2: Grain size evolution model parameters
	Fracture Criteria
	Fracture Probability Distribution (equation 4.25) d 0 , U 0 , γ, m u
	Energy Distribution (equation 4.40)	Ū e p(i)

  Table 7.1 and represented in Figure 7.2, more details are found in appendix B.1.

		Elasticity
	K ref Reference bulk modulus G ref Reference shear modulus	K(p) = K ref G(p) = G ref	p p ref p ref p	ne ne
	n e	Power of elastic behaviour
	p ref	Reference pressure
		Critical State and Plasticity
	φ pp	friction angle at critical state
	β	coefficient related to the slope of the critical state line (CSL) in isotropic loading
	d	coefficient related to the isotropic mechanism
	b	curvature coefficient for the yield function
	p co	Initial pressure of consolidation
		Flow Rule and Isotropic Hardening
	ψ( • ) characteristic angle
	α ψ	Dilatancy coefficient for deviatoric mechanisms
	a m	coefficient of the deviatoric mechanisms for monotonic loading
	a cyc	coefficient of the deviatoric mechanisms for cyclic loading
	c iso	coefficient of the isotropic mechanism for monotonic loading
	c cyc	coefficient of the isotropic mechanism for cyclic loading
		Table 7.1: Parameters of Hujeux's model
	In this case only one set of parameters has been defined, representing in this way that there is

Table 7

 7

	.2: Hujeux's parameters for hypothetical 200m high CFRD
		Elasticity
	K ref (M P a)	65.7	G ref (M P a)	49.3
	n e	0.3	p ref (M P a)	1.0
	Critical State and Plasticity
	φ pp ( • )	38	β	46.2
	d	2.15	b	0.1
	p co (M P a)	9.8
	Flow Rule and Isotropic Hardening
	ψ( • )	36	α ψ	1
	a m	0.005

Table 7

 7 For instance the optimization case No 1 is composed of two materials. One in zones 3B(1) and 3B(2) and another in zones 3C(1) and 3C(2), and only one parameter is optimized (G ref ), therefore the total number of variables are 2.

	Optimization No	Materials	Parameters Optimized	Total number of optimization variables
	1	3B(1)=3B(2)=3B 3C(1)=3C(2)=3C	G ref	2
	2	3B(1),3B(2) 3C(1),3C(2)	G ref	4
	3	3B(1),3B(2) 3C(1),3C(2)	G ref n e	8
	4	3B(1),3B(2) 3C(1),3C(2)	G ref , n e φ pp , β, d	20
	5	3B(1),3B(2) 3C(1),3C(2)	G ref , n e , φ pp , β, d P c 0 , b, ψ, a m , c	40
		Table 7.3: Optimization cases

.3. The cases differ on the number of materials (each material is associated to a zone of the dam, i.e., 3B(1), 3B(2), 3C(1), 3C(2)) and on the parameters to optimize for each material.

Table 7

 7 As for previous bi-dimensional model, the construction settlements were assessed for B section, but also for A and C. The results are shown in Figures 7.21 to 7.24.

	.4: Hujeux's parameters

Table A .

 A 1: Compilation of shear strength data for coarse granular materials

				Table A.1 -Continued
	No 1 No No No 3 4 5 2 Continued on Next Page. . . Material description Greywacke: quartz, feldspar, mica Material description Material description Material description Malpaso Material: Angular quarried material Purulia: Angular fragments of metamorphic rock with minerals including quartz, biotite and feldspar. Sandstone Amphibolite containing augite, biotite and quartz -rounded particles Slate (strong) Slate (weak) Basalt Pyramid dam material: argillite Crushed Basalt: fine grained olivine basalt, containing some phenocrystals of olivine and plagioclase feldspar Continued on Next Page. . . andesitic fragments El Infiernillo rockfill: Calas Material,angular particles of silicified conglomerate Ranjit Sagar: Rounded to subrounded alluvial, conglomerate,Sandstone, recent alluvium materials of older rocks and grits of chart and jasper, other quartzite, shale, claystone, 5 including sedimentary materials Rockfill material, sub-angular La soledad gravel: Angular dioritic fragments 4 Contreras gravel of alluvial origin, rounded particles, an-desitic and basaltic origin Gravelly Material From Pin-zadaran: Alluvial material 6 Crushed Cambric slate: angu-lar material containing quartz, muscovite, dolomite El Infiernillo rockfill:quarried by explosives, angular diorite fragments Basalt fragments, angular 7 Granitic gneiss particles (1) : dεv dε 1	dmax [mm] 25.4 dmax [mm] dmax [mm] dmax [mm] 38.1 51 152 200 25 50 38 11.4 80 38 11.4 51 152 11.4 200 80 50 25 200 51 152 70 200 200 40 200 200 200	Cu d 60 d 10 6 Cu d 60 d 10 Cu d 60 d 10 49 98 Cu d 60 d 10 38 6 7.1 7.1 82 42 33 58 77 35 28 12 37 49 -31 82 91 12 136 95 ---7.1 7.1 7.1 7.1 1.26 153 160 120 1.26 25 30 2.1 1.75 124 94 90 107 92 4.5 4.8 7.6 7.0 4.5 4.9 4.7 6.8 19 14	σc [M P a] 0.1 0.2 0.3 0.4 0.5 0.6 Table A.1 -Continued ( σ 1 σ 3 )max φ peak [ o ] 5.12 5.02 4.76 4.63 4.25 4.26 42.3 41.9 40.7 40.1 38.2 38.3 σc ( σ 1 σ 3 )max φ peak [M P a] [ o ] σc ( σ 1 σ 3 )max φ peak [M P a] [ o ] Table A.1 -Continued 0.039 0.088 7.75 7.22 50.5 49.2 σc ( σ 1 σ 3 )max φ peak [M P a] [ o ] 0.21 6.84 48.2 0.184 6.86 48.2 0.30 4.82 41.0 0.1 0.2 0.3 0.4 0.5 0.6 5.72 5.27 4.83 4.59 4.37 4.25 44.6 42.9 41.1 39.9 38.9 38.2 0.97 2.90 4.48 5.10 4.32 4.04 42.2 38.6 37.1 0.21 0.97 2.90 6.74 4.78 4.01 47.9 40.8 36.8 0.479 0.479 0.979 0.979 1.65 2.44 5.94 4.47 4.40 4.09 5.41 4.00 45.4 39.4 39.0 37.4 43.4 36.8 0.60 0.90 1.2 4.60 4.53 4.48 40.0 39.6 39.4 0.30 5.20 42.6 0.60 0.90 4.31 4.13 38.5 37.6 4.48 3.84 35.9 2.44 4.05 37.1 1.2 4.02 37.0 0.21 7.96 51 2.44 4.06 37.2 0.30 5.11 42.3 0.027 0.092 0.100 0.282 0.695 13.75 7.79 7.25 5.31 4.32 59.8 50.6 49.2 43.1 38.6 0.97 2.90 6.34 5.14 46.7 42.4 2.44 2.44 4.11 4.50 37.4 39.5 0.60 0.90 4.14 3.90 37.6 36.3 4.48 4.91 41.4 0.040 8.97 53.08 1.2 3.75 35.3 0.21 7.02 48.6 0.088 8.56 52.2 0.1 6.42 46.9 0.100 10.78 5.36 12.73 58.7 56.1 43.3 0.21 0.97 2.90 4.48 7.75 5.33 4.47 4.21 50.5 43.1 39.3 37.8 0.21 0.97 2.90 4.48 6.63 4.90 4.11 3.89 47.3 41.2 37.4 36.2 0.21 0.97 2.90 4.48 6.45 4.61 4.00 3.77 46.7 39.8 36.7 35.4 0.21 0.97 2.90 4.48 8.44 5.43 4.60 4.40 51.9 43.5 39.9 38.9 1.65 4.15 37.7 0.979 4.46 39.3 1.4 6.79 48.0 0.034 0.066 0.184 0.479 0.979 1.65 2.44 0.184 0.479 7.31 6.18 6.00 4.81 4.48 4.18 4.03 6.02 5.10 49.4 46.2 45.6 41.0 39.4 37.9 37.0 45.6 42.2 0.70 1.1 7.00 6.80 48.6 48.0 0.35 7.53 50.0 1.4 5.54 43.9 0.70 1.1 5.67 5.55 44.4 44.0 0.35 5.92 45.3 0.70 1.1 1.4 4.41 4.15 4.13 39.1 37.7 37.6 0.069 6.29 46.5 0.35 4.92 41.5 0.97 2.90 4.48 5.62 4.82 4.59 44.3 41.0 40.0 0.21 0.97 2.90 4.48 6.49 5.39 4.53 4.25 47.1 43.4 39.7 38.2 0.067 6.62 47.5 0.034 11.43 57.0 0.184 0.478 0.978 1.65 2.44 7.01 5.98 5.15 4.49 4.35 48.6 45.5 42.4 39.48 38.78 0.2 0.5 4.80 4.40 40.9 39.0 0.8 3.93 36.5 0.46 6.19 46.2 0.040 0.088 0.184 0.478 0.478 0.978 1.65 7.55 6.13 5.67 4.78 4.68 4.05 3.68 50.0 46.0 44.4 40.8 40.4 37.1 35.0 0.96 5.0 41.9 2.42 4.20 38.0 0.46 0.96 2.42 6.19 5.0 4.20 46.2 41.9 38.0	dεv dε 1 -0.03 (1) 0.01 0.06 0.03 0.00 0.03 dεv dε 1 dεv dε 1 0.19 0.55 dεv dε 1 0.31 0.90 0.27 -0.13 0.06 0.01 0.00 0.02 0.01 -0.03 -0.10 -0.08 0.31 0 -0.09 0.30 -0.10 0.60 0.25 0.22 1.00 0.16 0.24 0.19 0.30 0.22 0.27 -0.11 1.00 0.26 0.75 0.50 0.36 0.46 0.19 0.50 0.25 0.26 0.31 0.11 1.11 0.30 0.69 0.80 -1.39 0.53 0.55 0.04 0.00 0.71 0.47 0.17 0.40 0.05 0.00 -0.10 0.29 0 -0.08 -0.08 0.17 -0.07 -0.07 -0.07 0.42 0.01 -0.08 -0.08 0.27 0.37 -0.02 0.68 0.83 0.32 0.22 0.71 0.09 0.06 0.59 0.20 -0.02 -0.03 -0.00 -0.06 -0.05 -0.05 -0.04 -0.01 -0.01 -0.01 0.74 -0.01 0.29 0.18 0.07 0.55 0.19 0.08 0.01 0.71 1.33 0.77 0.67 0.60 0.05 0.08 0.024 0.00 -0.2463 0.17 0.71 0.65 0.55 0.15 0.00 0.68 0.02 0.00 -0.05 0.17 0.00 -0.05	σ p (2) [M P a] 136 σ (1) p [M P a] σ (1) p [M P a] σ (1) p [M P a] 136 -------------110 -------200 ---14.2-31.9 -100 73	n (3) [%] 30 n (2) [%] n (2) [%] 30 26 n (2) [%] 30 30 30 29 27 25 29 29 27 24 25 31 32 28 29 21 21 24 30 19 30 18 17 30 25 20 17 25 16 26 39 ≈ 20 25 32 25 32 31 28 28 31 30 27 26 30 30 27 24 29 28 27 25 29 29 41 38 35 35 38 34 31 34 33 15 18 18 46 16 15 15 18 18 17 16 41 35 24 24 24 25 26 38 38 39 28 31 38 38 42 38 35 34 25 26 25 26 28 26	S (4) 1 S (3) S (3) S (3) 1 1 1 0 1 1 0.7 0.5 0.9 0.55 1 1 1 1 0 1 0 1 1 0.0 0 0 1 0 1 1	ρ (5) s [ kN m 3 ] ≈ 26.4 ρ (4) s [ kN m 3 ] ρ (4) s [ kN m 3 ] ρ (4) s [ kN m 3 ] ≈ 26.4 28.1 28.1 26.8 26.4 28.4 25.5 27.6 26.6 27.5 25.7 25.7 25.7 28.1 26.8 26.3 22.6 28.4 28.4 24.5 22.6 -27 26.7 27.2 27.2	B [%] B [%] B [%] 2.1 7.5 B [%] -4.5 3.4 1.3 ----31.3 -9.9 21.6 27.2 8.5 14.8 27.6 -2.1 -8.1 8.9 -3.3 5.8 8.0 1.7 5.2 8.6 -33.8 12.3 12.8 -2.1 -2.7 ---18.5 5 7.2 --6.9 11.0 9.5 -15.3 2.7 -17.5 6.7 12.7 22.3 29.6 31.9 37.2 13.6 21.3 36.4 42.3 5.2 10.4 18.2 --11.9 20.5 ---13.3 15.3 14.3 10.6 --9.6 7.8 5.8 9.9 8.9 6.6 5.0 8.6 7.5 5.0 17 4.0 10.8 6.3 10.1 13.6 2.2 6.4 14 18 ---22.8 29.8 33.6 6.3 13.6 4.5 11.4 -12.6 10.8 8.8 -14 10.5 16.8 19.0 23.9 21.2 18.1
				2.44	3.91	36.3	0.19		32		-
	Continued on Next Page. . .

Table B .

 B 1: Classification of the Elastoplastic model parameters Directly measured Not-Directly measured Elastic K ref , G ref , n e , p ref Critical State and Plasticity φ pp , β, p co , d b Flow Rule and Isotropic hard.

Table B .

 B 4: Parameters of the model proposed by[START_REF] Salim | A new elastoplastic constitutive model for coarse granular aggregates incorporating particle breakage[END_REF]

	Elastic Behavior
	Shear elastic modulus	G
	Swelling-recompression parameter	κ
	Plasticity
	Critical state line	M, Γ, λ cs
	Breakage Factors
	Rate of energy consumption by particle breakage	β
	Plastic-distortional strains	θ, υ
	Breakage rate evolution	χ, µ
	Fitting parameter for initial stiffness	α
	B.6 Kohgo et al. (2007a) Model

Table B .

 B 5: List of parameters in Kohgo model

	Elasticity				Subload Plasticity
	κ	K i	G i γ p	α h	φ	φ cs φ pt R
	State Surface
	λ	Γ	e 0 01	a * 1	λ * f 1	s * m	e 0 02	a * 2	λ * f 2
	Soil water retention
	S re	S rm S rf s e	s m	s f c m c f
	Permeability
	k s	m p	n p

"Conventions Industrielles de Formation par la REcherche"

meters above sea level

the diameter corresponding to the 15% of material passing

in French school it is also called perfect plasticity (plasticité parfaite)

both related by a retention curve considered as unique in wetting or drying paths

number of contacts per particle

The validity of equation 2.20 is discussed in[START_REF] Evesque | Eléments de mécanique quasi-statique des milieux granulaires mouillés ou secs[END_REF], notably because Rowe's relation[START_REF] Rowe | The stress-dilatancy relation for static equilibrium of an assembly of particles in contact[END_REF] is not verified.

see definition in Appendix

the term "fractal" is employed here to characterize the repartition of the number of particles by size as in the theory of[START_REF] Mcdowell | The fractal crushing of granular materials[END_REF]

as defined in Appendix C

by "highly" we mean a size particles range of at least 1 order of magnitude

used currently in engineering practice

[START_REF] Rowe | The stress-dilatancy relation for static equilibrium of an assembly of particles in contact[END_REF] the increase of volume is taken positive (dilatancy)

relative with respect to the sample

International Commission for Large Dams

It is so because the variation of the grain-size distribution curve has not been linked yet mechanically to a strain due to crushing

is the sieve diameter for which 15% of the total mass passes

The original proposition of[START_REF] Wood | Strain softening and state parameter for sand modeling[END_REF] wasMp = Mc -k • ψe

Remerciements

Appendices A.1 Elementary frictional system

In section 1.4 the measured friction angle is supposed to be composed of two parts [START_REF] Terzagui | Soil Mechanics in Engineering Practice[END_REF]. The first component corresponds to the inter-particles friction (φ µ ) given by the friction between two flat surfaces and attributed mainly to mineral composition of the surfaces. The second component is given by a geometric interference (φ g ) component. Here, we would like to show this relation for an elementary frictional system. The elementary frictional system (Figure A.1) is composed of a rigid block in contact with a surface, at limit equilibrium is it required:

Looking at the microscopic level of the contact, the surface of contact is not completely smooth, but it presents some geometrical imperfections. Let us assume that these imperfections are represented as regular "tooth", characterized by an angle φ g (geometrical interference component). The surface of each teeth is characterized by a friction angle φ µ (inter-particles contact). At limit equilibrium each teeth verifies.

where i represents the i th contact. Assuming that there are k tooth along the interface and that there is no rotation, the limit equilibrium equations for the system lead to :

To solve the system let us suppose that all the angles and the forces are the same for all the tooth. i.e.f 1 n = f i n = ... = f k n = f n . And also that all the tooth are geometrically similar, therefore φ 1 g = φ i g = ... = φ k g = φ g . Using equation A.2, the limit equilibrium equations can be rewritten as:

Therefore, the macroscopic friction angle at limit equilibrium can be written as:

which leads to:

The interactions of the contacts in other systems, such as granular materials, can be much more complex. Therefore, a simple expression between the "macroscopic" friction angle and the microscopic characteristics is difficult to find.

Appendix B

Constitutive Models

B.1 Hujeux's Elastoplastic model

The ECP's elastoplastic multi-mechanism model [START_REF] Aubry | A double memory model with multiple mechanisms for cyclic soil behaviour[END_REF]; [START_REF] Hujeux | Une loi de comportement pour le chargement cyclique des sols. Genie Parasismique[END_REF], commonly called Hujeux model is used to represent the soil behaviour. This model can take into account the soil behaviour in a large range of deformations. The model is written in terms of effective stress. The representation of all irreversible phenomena is made by four coupled elementary plastic mechanisms : three plane-strain deviatoric plastic deformation mechanisms in three orthogonal planes and an isotropic one. The projection of any tensors (t) on the deviatoric plane (k) including e i and e j axes concerns only the t ii , t jj and t ij components, such that, we can write the corresponding mean and deviatoric values of strain and stress tensors as:

2)

The model uses a Coulomb type failure criterion and the critical state concept. The evolution of hardening is based on the plastic strain (deviatoric and volumetric strain for the deviatoric mechanisms and volumetric strain for the isotropic one). To take into account the cyclic behaviour a kinematical hardening based on the state variables at the last load reversal is used.

The model is written in the framework of the incremental plasticity, which assumes the decomposition of the total strain increment in two, elastic and plastic, parts. In what follows, a brief overview of the essential aspects of the constitutive model for primary loading paths is given.

The elastic part is supposed to obey a non-linear elasticity behaviour, where the bulk (K) and the shear (G) moduli are functions of the mean effective stress (p ) : where σ is the stress far from the crack influence zone, a is the length of the crack, and β is a dimensionless factor that is, the shape and the ratio a/D, where D is a generalized size parameter.

On the other hand, the theory of stress corrosion establishes that the atomic links at the tip of a crack are more vulnerable to a corrosive agent that those located far away from the tip. As part of their research the authors suggest that water has a corrosive effect and as consequence produces the softening of the material. This softening will develop cracks for lower values of Kc. The speed of cracks propagation is obtained from a thermodynamic formulation applied to the corrosion theory of constraints. The result is the following equation, which is represented in figure B.3.

Where, HR : Relative humidity ; K : stress intensity factor, T : temperature, R : universal constants of gas ; V 0 , E + and b : constants determined by experimental tests.

The thermodynamic balance at the environment implies that the liquid and gaseous phases held the same chemical potential, and then they produce the same corrosion effect constraints in other words the same speed of propagation of cracks. It means that whether the corrosive agent (water) is in a liquid or gas state it will produce the same corrosion. Therefore the only parameter controlling the speed of propagation of cracks are the relative humidity (HR). 3 and assuming that we can determine the intensity factor constraints K for a particle the authors conclude:

• if K K 0 , the crack doesn't grow at all.

• if K 0 < K < K c , crack grows at a defined velocity.

• if K K c , a catastrophic propagation of cracks occur.

We can see in this conceptual model that increasing the water content K 0 decreases and consequently reducing the area of non-cracking. From the oedometric tests containing a control in water content [START_REF] Oldecop | A model for rockfill compressibility[END_REF] have developed a model of compressibility in which they take into account the effect of suction (or water content) on the compressibility curves. The relative humidity in the test, is controlled by a vessel solution (NaCl) and measured with a hygrometer in the circuit. The psycrometric relationship between relative humidity and suction is given by:

where, s + π is the total suction ν is the molar volume of liquid water. The material used in the tests was a slate of angular particles for which the main characteristics are shown in Table B The results of tests are shown in figure B.4. The Test 1 corresponds to the saturated state. The fact of saturate the sample initiates a slight expansion of the specimen. Test 2 starts at dry state (corresponding to 100 MPa in total suction) then the vertical stress is increased until 1 MPa just before unload. After unloading a reload is made up to 0.2 MPa of vertical stress, at this state water content is increased without highlighting any distortion, finally the reload stage is continued until a vertical stress of 1 MPa showing an effect of water content on deformation. The sample of Test 3 also starts at dry state, as Test 2, but in this opportunity the wetting stage is performed when the vertical load reaches 0.6 MPa. In this occasion a "collapse deformation" has been observed while wetting, and afterwards a similar behavior to test 2 is observed. The last test has two wetting stages, the first one when vertical stress is near 0.4MPa, at that point a reduction in suction from 100 to 5 MPa is imposed in the test, and after that the load continues until reach a vertical stress of 0.9 MPa.Then the reduction in the suction was imposed again. As general observation collapse deformations has been observed in both cases.

The main reflexions of the authors are worth to be mentioned:

• It may be accepted that an only compression line exists for a given level of suction. • The oedometric compression stiffness depends on the suction level.

• An elastic domain can be defined by neglecting small hysteretic effects.

• The consequence of wetting (reducing suction) depends on the level of confining stress. For low values an increase in volume may be experienced while collapse deformations are experienced above a given threshold of confinement stress. This is consistent with the theory of cracks propagation summarized in figure B.3.

• Just a water content corresponding to the pores saturation is needed to cause collapse deformation.

Moreover, [START_REF] Oldecop | Theoretical investigation of the time-dependent behaviour of rockfill[END_REF] believe that the long term behavior of rockfills is not negligible. The Compressibility model proposed by [START_REF] Oldecop | A model for rockfill compressibility[END_REF] is based on two mechanisms of plastic deformation, called

• Instantaneous deformation mechanism (IDM) which is independent of the action of water and present at any stress level, and

• Time-dependent deformation mechanism (TDM), which depends on the action of water and is active beyond a threshold of the total stress value σ y . This parameter represents the threshold which defines the beginning of particle crushing.

The action of water is represented eather by the water content (ω) or by the total suction (ψ), these two quantities are linked together by a retention curve supposed unique. The increments of stress for the IDM and TDM are described by the following relations:

The slope of the line in the stress-strain space (σ verticalε vertical ) is λ i when the mechanism IDM is active (dry state). When the TDM mechanism is active the slope is λ i + λ d (w). The exponents i and d correspond to the mechanisms IDM and TDM.

The suction in the model is correlated with the water content inside the particles. The stress-strain relationship for the unload path is given by:

The slope κ is supposed to be independent of the water content. In order to represent the swelling deformation as consequence of the increase in water content at low stress levels [START_REF] Oldecop | A model for rockfill compressibility[END_REF] propose the following expression.

The previous equation can be expressed in terms of suction as:

The yield surface in the σw space is given by:

Where σ y is the clastic yield stress and σ * 0 is the yield stress of the dry state rockfill but also a strain hardening parameter given by:

By the consistency condition, the flow rule is given by:

Then the compressibility coefficient is made a function of the water content as: The formulation of the plastic potential supposes a non-associated rule for the deviatoric part and an associated for the isotropic part. These potentials are given by:

In this work [START_REF] Chávez | A constitutive model for crushed granular aggregates which incluyes suction effects[END_REF] used the proposition of [START_REF] Wan | A simple constitutive model for granular soils: Modified stress dilatancy approach[END_REF] in order to take into account the effect of initial density over the volumetric behavior by using: In the last relation ψ m es the mobilized dilatancy, φ m is the mobilized friction angle, φ cr is the friction angle at critical state, e is the current void ratio, e cr is the critical void ratio and α is a constant of the model. The results of the model and test are showed in figure B.6.

We can observe through equations B.36 to B.44 that the number of parameter (summarized in table B.3 ) to be determined is significant.

B.5 Elastoplastic Constitutive Model Incorporating Particle Breakage

An elastoplastic model that takes into account the degradation of particles during triaxial tests was proposed by [START_REF] Salim | A new elastoplastic constitutive model for coarse granular aggregates incorporating particle breakage[END_REF]. This model is based on the analysis at the particle level and considering an energy dissipation by breakage. The model is developed based on triaxial results on crushed basalt. The main equations are shown hereafter, however no deduction is done. The yield equation is given by

(B.49) where φ f is the angle of friction corrected by dilatancy as suggested by [START_REF] Rowe | The stress-dilatancy relation for static equilibrium of an assembly of particles in contact[END_REF]. If the breakage of particles is ignored at critical state equation B.49 reduces to the critical state relationships as showed in [START_REF] Indraratna | Modelling of particle breakage of coarse aggregates incorporating strength and dilatancy[END_REF].

The flow rule is given by:

where dE B is the energy consumption due to particle breakage per unit volume. The previous equation is rewritten after finding the correlation of dE B with breakage. This passage make to appear several parameters as shown in figures B.7. finally the new flow rule is expressed as: where

As the model follows the plasticity theory, a plastic potential is needed. Actually only its derivatives are required and in that sense the authors propose:

To formulate the yield loci and the hardening function [START_REF] Salim | A new elastoplastic constitutive model for coarse granular aggregates incorporating particle breakage[END_REF] postulate two statements with reference to granular aggregates.

• Coarse granular aggregates experience plastic strains when, and only when, there is a change in stress ratio η.

• The undrained stress paths are parabolic in the p-q plane and are expressed by a the relationship

1-p 0 /pcs . Finally [START_REF] Salim | A new elastoplastic constitutive model for coarse granular aggregates incorporating particle breakage[END_REF] shows the application of their model and the comparison with the experimental data as showed in figure B.8 . The results of the model seems satisfactory, something to remark is that the difference between the results with and without taking the breakage into account are not very different in the planes qε 1 and ε vε 1 however there is a good assessment in the evolution of breakage.

The parameters of the model are summarized in table B.4

Here we are interested in the second moment of the retained curve in the final gradation, it means the J 2u . An issue rices here for the estimation of the equation C.1 in a discrete curve. We have adopted to make an approximation as that used to evaluate integrals or areas. Then the procedure is resumed as follows:

• In the grain size distribution F u convert the x scale to an arithmetic scale • to derive numerically the F u and find the retained mass curve p(D).

• to apply

The method just mentioned have been compared by using a fractal distribution, for which the analytical value of J 2u = 16.66mm is known. We found that the discretization level (number of division along the x axis) is important to find an acceptable response, figure C.1 shows the approximation to the theoretical value by increasing the number of divisions over the experimental curve.

Error Estimation In order to compare the different variations of the parameters over the model response let us define an error parameter estimation given by the quadratic error (equation C.3).

where E is the error y exp is the experimental value of y variable, y model is the estimated value an n is the number of observations.

Appendix to Chapter 5

In chapter 5 a representation of the relation between the input and output power in a soil's sample has been proposed. The particularity of this plot is that it allows to identify the tendency of the material to store or dissipate the energy. In this sense some basic behaviors should be identified in order to compare against the experimental results. This appendix is focused on the theoretical description of a linear elastic and a pure frictional mechanisms on the power in-power out plot.

D.1 Power in -Power Out in linear elastic behavior

Let place us in the coordinate system of three mutual orthogonal axis representing the principal direction of the stress and strain tensors. The elastic behavior is represented by:

where λ and µ are the Lamé's coefficients. In the case of triaxial test σ ij = 0∀i = j and σ 22 = σ 33 = σ 0 , ε 22 = ε 33 , and ε 33 ε 11 = -ν. Therefore the relation between the total input power and the total output power can be expressed by: Appendix E

Appendix to Chapter 6

In soil's mechanics the mechanical characteristics of materials are deduced from laboratory or in-situ tests. However, sometimes the dimensions of the particles composing the material are very important, turning laboratory tests not feasible and in-situ test of difficult interpretation. In this case it has been accepted that parallel reduced gradations (e.g. [START_REF] Fumagalli | Test on cohesionless materials for rockfill dams[END_REF]) can be tested under laboratory conditions to represent the behavior of a larger material. However, this extrapolation can be affected by the level of particles crushing including a size-scale effect. It has been shown in chapter 6 how the size-scale effect can be taken into account to extrapolate in a more reliable way the mechanical behavior of materials with large crushable particles. The theoretical framework has been presented in chapter 6, and applied to a particular case named Pyramid dam material. This material has been tested by [START_REF] Marachi | Strenght and deformation characteristics of rocfill materials[END_REF], who also tested two additional materials, a crushed basalt material, and another material coming from the dam of Oroville. As the results of Pyramid dam material has already been shown, this appendix shows the application of the size-scale effect to the remaining materials.

E.1 Crushed Basalt Material