Keywords: Rockfill, stone, hydraulic structure, seismic, FDM, DEM, continuum-discrete

Firstly, I'm truly grateful to my academic supervisor in LTDS-ECL, Eric Vincens. I was really lucky to work with a supervisor who offered me the flexibility to explore on my own while still providing the support I needed to recover when my steps slipped. He has been an incredible guide, teacher, and supporter throughout my research journey. His knowledge and dedication helped me grow as a researcher.

I also want to thank my co-supervisor at ITASCA Consultants, Fabian Dedecker. His advice and expertise significantly improved the quality of this thesis. His feedback was essential in shaping my research and solving the technical and theoretical problems that we faced during the Ph.D. I would like to sincerely thank Fabrice Emeriault and Fernando Lopez-Caballero for accepting to be the reviewers of the manuscript. Their comments have improved the quality of its content. I also want to thank the examiners, Pierre Breul who is also the jury president, Paulo Lourenco and Claire Silvani. I'm thankful to all my colleagues at ITASCA Company for their support and the maintenance of all my work requirements. I would like to thank the engineers who helped me to overcome the technical problems during my PhD including Marco Camusso, Huy Tran, Sacha Emam, and Rima Ghazal. A special thanks also to the actual director, Lauriane Bouzeran and the previous director, Danial Billaux for their support. I extend my appreciation to my colleagues at LTDS laboratory in ECL and to it itself. Our shared knowledge, and friendly atmosphere made our research environment productive and enjoyable.

I also want to thank EDF Company for the funding of the project and for the essential provided data. I would like to thank our contacts in EDF including Roland Plassart, Philippe Kolmayer and Jeremie Moreau. They collaborated with me, supported me, and provided valuable data.

To my friends, you were a main part of my journey. Your encouragement and support played a crucial role in my success. I am so lucky to have you by my side.

Lastly and most importantly, my family members deserve a special mention and thanks. Their unconditional love, constant encouragement, and constant belief in my abilities have been my bedrock throughout this journey. Their sacrifices and support have been my guiding light, and I am forever indebted for their unwavering presence in my life.

Thank you all for being a part of my journey.

Abstract

Rockfill dams with dry-stone pitching are specific hydraulic structures generally built in the early XX th century. Their specificity lies in the existence of upward and downward faces made of stones that were hand-placed without mortar in a way similar to what can be found in dry-stone retaining walls. Few studies have been carried out to quantify the specific and critical role of the pitching in the stability of such dams. The main objectives of this work are first to set up a numerical model able to retrieve the mechanical behavior of such type of dams with reasonable computational times and secondly to better understand the mechanisms that explain their resistance. Gaining experience from past numerical studies on dry-stone retaining walls, a mixed FDM-DEM (continuum-discrete) approach for the modeling of rockfill dams with dry-stone pitching is herein proposed. Considering the overall isotropic nature of the dam body, constructed from dumped rockfill, this latter is modeled as a continuum medium. In contrast, the individual stones comprising the pitching are individually modeled to accurately reflect the discrete and anisotropic nature of these peculiar dry-stone structures.

A validation of the numerical tool is first proposed on the basis of results from the literature involving experiments on scaled-down rockfill dams with dry-stone pitching. The effect of the pitching characteristics (dimensions, weight, material, stone shape, among others) on their behavior is consequently studied.

Secondly, a case study is proposed on a specific rockfill dam with dry-stone pitching located in France that was studied in the past using different approaches than that proposed herein. First, the modeling of the construction and water-filling phases allowed us to confront the numerical results with available historical measures and to validate the chosen set of model parameters. Then, both the static and pseudo-static resistances of this dam are evaluated through the design of different safety factors. The quantification of the role of the different elements constituting the dam to its overall resistance is discussed in detail.

Finally, the seismic resistance of this specific dam is studied through true dynamic computations and different proposals to quantify it are proposed. The chosen set of model parameters was validated by comparing the obtained crest settlement at the end of the shaking with the postearthquake crest settlements of several reviewed dams.

Résumé

Les barrages en enrochement avec perré en pierre sèche sont des ouvrages hydrauliques bien particuliers, généralement construits au début du XX me siècle. Leur spécificité réside dans l'existence de faces amont et aval constituées d'une ou plusieurs rangées de pierres placées à la main sans mortier, à l'instar de ce que l'on peut trouver dans les murs de soutènement en pierre sèche. Peu d'études ont été menées pour quantifier le rôle spécifique du perré dans la stabilité de ces barrages.

Les principaux objectifs de ce travail sont d'une part de mettre en place un modèle numérique capable de retrouver le comportement mécanique de ce type de barrages, en des temps de calcul raisonnable, et d'autre part de mieux comprendre les mécanismes qui expliquent leur capacité de résistance.

Grâce à l'expérience acquise lors d'études numériques précédentes sur les murs de soutènement en pierre sèche, une approche mixte FDM-DEM (continue-discrète) pour la modélisation des barrages en enrochement, avec perré en pierre sèche, est ainsi proposée. Compte tenu de la nature globalement isotrope du corps du barrage, construit à partir d'un enrochement déversé, ce dernier est modélisé comme un milieu continu. En revanche, les pierres composant le perré sont modélisées de manière individuelle pour refléter précisément la nature discrète et anisotrope de ces structures très particulières.

Une validation de l'outil numérique est d'abord proposée sur la base des résultats de la littérature, comprenant des expériences sur des barrages en enrochement à échelle réduite avec perré en pierre sèche. L'effet des caractéristiques du perré (dimensions, poids, matériau, forme de la pierre, entre autres) sur son comportement est ensuite étudié.

Puis, une étude de cas est proposée, en considérant un barrage en enrochement avec perré en pierre sèche, situé dans les Hautes-Pyrénées (France), et étudié par le passé en utilisant des approches différentes de celles proposées ici. Tout d'abord, la modélisation des phases de construction et de remplissage nous a permis de confronter les résultats numériques aux mesures historiquement disponibles puis de valider l'ensemble des paramètres du modèle. Les résistances statiques et pseudo-statiques de ce barrage sont évaluées selon différents concepts de mesure du facteur de sécurité. La quantification du rôle des différents éléments constituant le barrage dans sa résistance globale est également discutée en détail.

Pour finir, la résistance sismique de ce barrage spécifique est étudiée par des calculs dynamiques réels et différentes approches pour la quantifier sont proposées. Le jeu choisi de paramètres du modèle a été validé en comparant le tassement de la crête obtenu à la fin du tremblement avec les tassements de la crête après séisme de plusieurs barrages examinés.

General Introduction

A rockfill dam with dry-stone pitching is a dam composed of decametric blocks forming its body which is protected by a layer of hand-placed stones. These dams are characterized by a very steep slope of about 1/1. Several dams of such type were built in France in the early 20 th of which about 10 are still operated by the stakeholder " Électricité de France" (EDF).

The mechanical behavior of this type of dam is not well understood with very few past studies. Past modeling using a 2 dimensional fully DEM approach showed that the pitching plays a significant role in the dam's stability allowing it to resist seismic motions up to 0.5g [START_REF] Deluzarche | Modélisation discrète des enrochements: Application aux barrages[END_REF][START_REF] Vincens | Dry Stone Retaining Structures: DEM Modeling[END_REF]. Such findings needed confirmation which motivated EDF to fund scaled-down experiments on such hydraulic structures within the research project PEDRA (2011-2014) funded by the French Ministry of Ecological Transition. Pseudo-static tests were carried out by progressively inclining 1/10 scaled-down dams up to failure.

The PhD thesis project "Rock'hy" was then born in 2020 and is presented in this manuscript.

It is carried out within the framework of a CIFRE contract between the stakeholder "EDF" and "ITASCA Consultants SAS" in collaboration with the LTDS laboratory at Ecole Centrale de Lyon.

The project aims to introduce and validate a new mixed FDM-DEM approach for the modeling of rockfill dams with dry-stone pitching which would then help to better understand both their static and dynamic behavior and the influence of each of their composing elements on their safety and behavior.

This thesis is divided into four chapters:

The first chapter presents a general introduction and literature review of the rockfill dams with dry-stone pitching. The historical overview and development of such dams are first shown with a focus on one of the dams that is still operated by EDF (Escoubous Dam). The mechanical behavior of the rockfill including its strength and deformability is then reviewed. In addition, a description of the PEDRA experiments on scaled-down rockfill dams with dry-stone pitching is presented. Finally, the previous modeling attempts and approaches of Escoubous dam and other dry-stone structures are shown.

The second chapter suggests a mixed continuum-discrete (FDM-DEM) approach for modeling rockfill dams with dry-stone pitching. The dam body was modeled as a continuum medium using an elastoplastic constitutive model and the pitching stones were modeled discretely using a linear contact model. The validation of this approach was carried out by the modeling of the PEDRA pseudo-static experiments on scaled-down rockfill dams with dry-stone pitching.

In the third chapter, the modeling of a full-scale rockfill dam with dry-stone pitching, called Escoubous and located in the Pyrenees Mountains in France is presented. The safety of the dam has been evaluated using two suggested approaches. Additionally, the effect of each of the elements composing the dam on its safety and behavior is presented.

The fourth chapter deals with the true dynamic analysis of Escoubous dam. The dynamic model with its adapted boundary conditions and parameters is presented. A Mohr-Coulomb (MC) model which is coupled to the hysteretic damping model has been used for the dam body. The obtained crest settlement at the end of the shaking is compared to the post-earthquake settlements of several reviewed dams which validated the chosen set of model parameters. Finally, two dynamic safety assessment approaches are suggested and used.

Chapter I

Literature Review

Introduction

The construction of rockfill dams is considered one of the earliest hydraulic engineering methods; evidence of this may be found as far back as ancient Mesopotamia. These dams are constructed for the purposes of renewable energy that participate to meet the electricity demand and for the retention of water in reservoirs allowing irrigation of fields and domestic use. Some of them also play an important role in flood control. The most basic technology of dam construction has advanced throughout time, and the largest works of this sort now reach heights of more than 200 m.

The design of a rockfill dam has to consider two main features that define its resistance on both the short and long term. The first one is related to its ability to resist to static and seismic loads. This stability is ensured by an appropriate geometry together with the strength and resistance of the dam materials (rockfill, pitching, etc.) which strongly depends on its mechanical characteristics.

In addition to the stability, the dam has to be well-sealed to retain the water and avoid any leakage which would alter its stability and would lead to loss of the retained water. This sealing is ensured either by a clay or earth core, or a waterproof mask on the dam upstream slope.

Rockfill dams became popular in the 20th century, and a lot of projects using this technique popped up, projects that would be qualified as ambitious and even reckless. Thus, despite current slope fruit standards of 1.3-1.4, many dams, particularly in France and Italy, show fruits lower than 1. These constructions typically include a pitching on the slopes, which is a layer of dry stones that have been manually positioned to protect the dam body from climatic hazards and animals. In fact, rockfills constituting the dam body alone would not withstand such slopes lower than 1 since they are generally characterized by a repose angle lower than 45°. The pitching stabilizes the rockfill surface probably by its own weight and by the resistance induced by the friction between its stones [START_REF] Deluzarche | Modélisation discrète des enrochements: Application aux barrages[END_REF][START_REF] Thomas | Review of shearing strength of rockfill[END_REF].

In this literature review chapter, the development of such dams is shown in addition to the damages that they suffered throughout their serving period with a focus on a particular dam called Escoubous. In addition, the particularities, and characteristics of the rockfill and the dry-stone pitching composing such dams are shown and discussed. Finally, past experiments on such dams and associated numerical modeling are presented [START_REF] Deluzarche | Modélisation discrète des enrochements: Application aux barrages[END_REF][START_REF] Chen | Modélisation du comportement mécanique des grands CFRD: Identification des caractéristiques des enrochements et comportement du masque d'étanchéité amont[END_REF].

History Overview

2.1. Dry-stone building sector and historical usage

Rockfill dams with dry-stone pitching belong to the class of dry-stone structures and the very dry-stone pitching to the class of dry-stone retaining structures. Dry-stone structures generally characterize structures where rubble stone blocks were carefully hand-placed without mortar in such a way to provide a long-term resistance to the structure. The dry-stone building technique was already known by the human beings of the Neolithic, near the end of the Stone Age, when civilization progressively transitioned from hunting and gathering to an agricultural one. Then, dry-stone structures are one of the first known and executed projects by human beings in civil engineering history [START_REF] Kerisel | History of retaining wall design[END_REF]. Since then, a lot of cases of such structures can be found worldwide as soon as the resource was available. Materials from the area are always utilized, whether they are dug out in a nearby quarry, as was the case for Macchu Pichu [START_REF] Bingham | Inca land: explorations in the highlands of Peru[END_REF], or more directly from a nearby river, as was the case for the Borobudur temple [START_REF] Soekmono | A monument of mankind[END_REF]. Most famous and known examples of this heritage are "Stonehenge" in England, the "Iron Age Brochs" in Scotland ( [8]), the "Great Zimbabwe" in Zimbabwe, and more recently "Machu Picchu" retaining walls and dwellings in Peru. • Rockfill dams with dry-stone pitching: stone blocks placed without mortar constitute the faces of the dam; the pitching also contributes to the stability of those dams,

Some examples of dry

• Dry-stone walls retaining slopes and roads: The weight of the very wall allows it to resist the backfill loading and to the traffic load on the road,

• Dry-stone bridges: A vault made of cut or rubble stones without mortar transfer the loads to the abutments or the foundations.

In this dry-stone heritage, the rockfill dams with dry-stone pitching on their faces are one of the most remarkable and specific structures. Their inner backfill is usually made of decametric blocks which are quarried from the dam site. The faces of those dams are covered with a layer of carefully hand-placed rubble or cut stones.

Development of Rockfill Dams with dry stone pitching

Rockfill dams with dry-stone pitching on both faces became popular in the early 20th century. However, this building technique was used as far as in the pyramid age (about 2600 BC) in Egypt [START_REF] Nicholas | A history of dams: the useful pyramids[END_REF]. Indeed, EL-Kafara dam in Egypt is the world's oldest known water-retaining construction of this size (Figure I.2). Since then, dams' design has varied from country to country, with dropped rockfill zones for backfill occasionally mixed in with hand-placed rockfill zones [START_REF] Léon | Leçons sur les grands barrages[END_REF]. In the panel of newer rockfill dams constructed in the 20th century, two different backfill designs coexist: one typical of Italy and the other primarily created in the United States (Figure I.3). For example, because of the availability and low cost of qualified labor in Italy, the rockfill of the backfill could be completely hand-placed. In the USA, due to the shortage of qualified masons at that time, rockfill was typically hand-placed solely under the upstream face. the American design but with certain Italian quirks such as the steep downstream fruit and faces covered by dry-stone pitching. The majority of those dams are still run by the French stakeholder "Electricity of France" (EDF). They are mainly composed of a backfill made of either dumped or hand-placed rockfill. A layer of cut or rubble stones is arranged on the faces of the dam forming a rigid and dense pitching of blocks. In the case of a dumped rockfill in the dam body, a transition layer of hand-placed rockfill can be sometimes found under the pitching to form a planar surface on which the pitching stones are laid. It allowed to smooth the surface and to avoid local stress concentrations [START_REF] Vincens | Dry Stone Retaining Structures: DEM Modeling[END_REF]. Some of those dams include a berm or more on its downstream surface. A concrete mask is generally casted over the upstream pitching to form a sealing and to prevent water leakage. The international construction design at the time for dumped rockfill had three main features: selection of a quarry of hard and unaltered rock, dumping blocks of various dimensions in the dam body, and sluicing large quantities of water in the rockfill. Nevertheless, it is not guaranteed that those conditions were respected during the construction of French dams. Disrespecting these steps would lead to high settlements after strong rainfalls. The rockfill type was mainly granodiorite, granite, gneiss and schist having an average compressive strength of about 100 MPa. The average size of rockfill is about 30 cm; the material is generally narrowly graded with an average coefficient of uniformity between 2 and 3. The average porosity of the dumped rockfill is around 0.4 and that of the hand placed stones is around 0.2 [START_REF] Vincens | Dry Stone Retaining Structures: DEM Modeling[END_REF].

Different types of stone pitching can be found two of which are dry without mortar: sized blocks without mortar ( ). A probable cause of this dismantling could be the arching effect through the dam body between the upstream and downstream pitchings as mentioned by [START_REF] Léon | Leçons sur les grands barrages[END_REF] and [START_REF] Vincens | Dry Stone Retaining Structures: DEM Modeling[END_REF].

• Bulging or buckling of several rows of the drystone pitching (Figure I.8d). This is usually reported on pitching made of rectangular sized blocks (Figure I.7a). This kind of pitching is very stiff. These local bulging occurred when the pitching tried to adapt to the deformed shape of the dam after settlement. The rockfill that constitutes the body of the dam was progressively dumped during construction without being compacted and therefore corresponds to a loose state. Its properties are defined by a porosity (n) estimated at 0.45, a D 50 of 15 cm and a uniformity coefficient C u < 3. The pitching is made of 30 cm square granite blocks. It is laid on a layer of rockfill arranged by hand allowing the construction of the pitching in good conditions: the thickness of this sub-layer was identified before rebuilding deformed zones of the pitching (during maintenance). This sub-layer is 90 cm thick on the downstream side and 20cm thick on the upstream side.

Historical Data

There is almost no historical data related to the tracking of the dam deformations during its construction and impoundment in 1954. However, early measurements (1955-1959 period when the reservoir water level elevation rose 15 m) indicated a maximum downstream displacement of 50 mm during this process.

By extrapolation of the theory established on structures of this type (Gréziolles dam in particular, which is operated by EDF and has been much more auscultated throughout time), it should be expected at the dam crest during the first impoundment:

• settlements of 0.7% of the dam, i.e. 11.2-12.6 cm;

• a horizontal downstream displacement of 0.35% of the dam height, i.e. 5.6-6.4 cm.

Plots E, F, G and H in

Mechanical Behavior of rockfill

Rockfill is a type of granular material that is characterized by a size larger than 8 cm [START_REF] Guyer | An Introduction to Identification and Classification of Soil and Rock[END_REF]. All of the granular material types share some common mechanical behavior traits but in some cases, the behavior of the rockfill differs and has its unique features.

General characteristics

The important size of the rockfill makes this material peculiar within the class of granular materials since their mechanical characterization is difficult to provide. However, like other granular materials, some of their main geometric properties include the shape of the rockfill (elongation ratio and angularity), the size of the rockfill, its grading and the initial porosity state.

The shape of the rockfill can be characterized by two variables: the particle form (elongation or aspect ratio) and the particle angularity (roundness). Elongation represents the general form of the particle whether it is elongated (the maximum dimension of the particle is much higher than its minimum dimension) or uniform (not elongated: the maximum dimension of the particle is similar to its minimum dimension). The elongation of the particle is measured by the elongation ratio ER which is the ratio of the smallest dimension to the largest dimension of the particle (Equation I.1).

ER = d min d max (I.1)
The particle angularity characterizes the particle shape locally. It depends on the form of the vertices forming the particle. The rockfill particles of fluvial origin are normally rounded particles whereas the quarried rockfill is usually angular. The angularity can be measured using the variable r v (Equation I.2) as proposed by Marsal ([14]). This variable is equal to 1 in the case of a sphere.

r v = 6 V π. d3 (I.2)
where • V is the particle volume;

• d 1 , d 2 and d 3 are the particle dimensions in three orthogonal axes;

• d is the average dimension of the particle ( d = (

d 1 + d 2 + d 3 )/3)
Several charts were proposed to classify the granular particles as a function of their form as the one proposed by Barton ([15]). The chart is shown in Figure I.14 As mentioned by Chen [START_REF] Chen | Modélisation du comportement mécanique des grands CFRD: Identification des caractéristiques des enrochements et comportement du masque d'étanchéité amont[END_REF], the form of the rockfill particles has an effect on the porosity: the more angular the blocks, the higher the porosity obtained. Higher density can be obtained with rounded particles. Fewer contacts between particles are found with higher angularity and elongation (not totally accurate with ordered material) which increases the local stresses at the contacts and the possibility for the crushing of the rockfill. Marachi ([16]) has shown that the internal friction angle of rockfill rounded particles is 3 to 4°higher than rockfill angular particles crushed under pressures higher than 100 psi (690 kPa). The opposite trend is observed on sand.

The rockfill grading is also a key factor affecting its mechanical behavior. It includes the dimension of the rockfill material and the heterogeneity of the distribution of the rockfill block sizes. From this distribution, the uniformity coefficient C u = d 60 /d 10 can be evaluated: the higher this coefficient: the more the heterogeneity and spread are found in the rockfill block sizes. According to Leps ( [START_REF] Thomas | Review of shearing strength of rockfill[END_REF]), several conclusions can be derived in relation to the effect of this distribution on the rockfill behavior. The lower the uniformity coefficient of the rockfill, the higher the rockfill porosity is (for same compaction energy) and the higher the block breakage is (due to lower contact points and then higher local stresses), and the lower the friction angle is.

The porosity of the rockfill reflects its arrangement and state and it has an important effect on its mechanical behavior and resistance. Indeed, dense granular material (with low porosity) have higher shear strength [START_REF] Marsal | Mechanical properties of rockfill[END_REF]. An increase of porosity leads to a lower internal friction angle of the granular material as shown by several authors [START_REF] Barton | Shear strength of rockfill[END_REF][START_REF] Marschi | Evaluation of properties of rockfill materials[END_REF][START_REF] Nakayama | Selection and quality control of materials for rockfill dam of pumped storage project in phyllocrystalline schistose area fourteenth congres on large dams[END_REF]. From the triaxial tests on Monterey sand shown in figure I. [START_REF] Barton | Shear strength of rockfill[END_REF], we can see that the dense material shows a resistance peak followed by a softening phase. The intensity of this peak decreases remarkably in the case of loose material.

The dense material has a larger internal friction angle associated to dilative deformations when the confinement is low. On the other hand, when confinement is high, dilatancy is severely constrained, which causes the peak friction angle to drop. This behavior can be generalized for all granular materials including rockfill.

Shear strength criteria

The shear strength of any type of soil including the rockfill is an important characteristic which defines the stability in case of geotechnical structures. The main parameter that characterizes the rockfill's shear resistance is the maximum mobilized friction angle or internal friction angle.

The shear strength is normally measured using laboratory tests and mainly using a cylindrical triaxial apparatus. The problem with rockfill is there big dimensions which could reach 1m of average diameter. No triaxial cells are that big and can hold a Representative Elementary Volume with such rockfill. There are several approaches in the literature that addresses this subject in order to estimate the rockfill shear resistance. It is proven that the internal friction angle for rockfill increases with the decrease of confining pressure [START_REF] Vincens | Dry Stone Retaining Structures: DEM Modeling[END_REF][START_REF] Thomas | Review of shearing strength of rockfill[END_REF][START_REF] Barton | Shear strength of rockfill[END_REF][START_REF] Frossard | Comportement macroscopique des matériaux granulaires mis en oeuvre dans les barrages[END_REF]. Consequently, the rockfill failure surface that delimits the domain of acceptable stresses able to be born by the material is not a straight cone but flattens at high pressures. In addition, rockfill is often characterized by an axial deformation of about 10% to 20% at peak which is higher than that of the sand [START_REF] Chen | Modélisation du comportement mécanique des grands CFRD: Identification des caractéristiques des enrochements et comportement du masque d'étanchéité amont[END_REF].

De Mello shear criterion

De Mello suggested a non linear equation (Equation I.3) that describes the variation of the rockfill shear resistance τ versus the normal stress σ n (normal to the failure plane) applied on the rockfill [START_REF] Victor | Reflections on design decisions of practical significance to embankment dams[END_REF].

Leps diagram

The shearing strength of rockfill was reviewed by Leps [START_REF] Thomas | Review of shearing strength of rockfill[END_REF] by collecting experimental data for triaxial tests available for rockfill [22][START_REF] Taylor | Shearing properties of ottawa standard sand as determined by the mit strain-controlled direct shearing machine[END_REF][START_REF] Wesley | Triaxial shear tests on pervious gravelly soils[END_REF][START_REF] Hall | Triaxial testing with large-scale high pressure equipment. Laboratory shear testing of soils[END_REF][START_REF]Research on the behavior of granular materials and rockfill samples[END_REF][START_REF] Raúl | Large scale testing of rockfill materials[END_REF]. The shear test results have been re-interpreted when necessary to plot the internal friction angles as a function of normal pressures across the failure plane, as deduced from the use of the Mohr diagram [START_REF] Thomas | Review of shearing strength of rockfill[END_REF]. The obtained data are plotted in Figure I.17 where the shear strength is given by the value of the internal friction angle against the normal stress. Three trend lines have been deduced from the collected data and plotted on the same figure: the top one represents the case of dense, well graded and strong rockfill, the bottom one represents the loose, poorly graded and weak rockfill and that which is located between the two others the average rockfill.

From these data, we can notice how the shear strength increases remarkably with the decrease of the pressure. The friction angle can reach values up to 60°at low pressures for dense rockfill. In addition, from the plotted data of angular sand and Ottawa sand in Figure I.17, one can note that rockfill is characterized by a much higher internal friction angle than the sand at the same pressure.

Barton

The friction angle is expressed as a function of the stress state and state factors including porosity, grain shape, grain strength, and grain size in Barton and Kjaernsli's ( [START_REF] Barton | Shear strength of rockfill[END_REF]) strength criteria. The formula for the internal friction angle reads I.5:

ϕ = R.log( S σ ′ n ) + ϕ b (I.5)
where:

• R: Equivalent roughness depending on the shape of the blocks (angular or rounded) and the porosity ( • ϕ b : Basic friction angle which corresponds to the contact friction angle between two polished rock faces. It can be measured by a "tilt test". The principle of the test is: two rock plates with a flat surface (at least 1775 mm 2 ) are installed on a tilting table. The mock-up is rotated slowly until the upper block starts to slide. The angle existing between the flat surface with the horizontal measured at failure defines the basic friction angle ( [START_REF] Usbr | Procedure for determining the angle of basic friction (static) using a tilting table test materials engineering and research laboratory[END_REF]). This value generally ranges between 25°and 35°. • ψ: critical angle or non-dilative angle of friction (ϕ crit );

tan 2 ( π 4 + ϕ 2 ) = [1 - dε v dε 1 ]tan 2 ( π 4 + ψ 2 ) (I.6) ψ = φ b + θ M (I.7) ϕ -ψ = A.[I D (10 -ln(p ′ ) - 1
• φ b : basic friction ratio of a sawn surface of rock;

• θ M : Basic friction angle increase due to the particle roughness [START_REF] Md Bolton | The strength and dilatancy of sands[END_REF];

• p ′ : Mean effective stress in kPa;

Deformability and Compressibility

The deformability of rockfill material depends on several characteristics including the elongation and angularity, the grading, and its initial porosity. It has been deduced experimentally that the axial and volumetric deformations in a triaxial test are more important for angular particles than for rounded ones [START_REF] Marsal | Mechanical properties of rockfill[END_REF][START_REF] Marschi | Evaluation of properties of rockfill materials[END_REF]. Marachi ([17]) studied the deformation of rockfill used in Oroville dam. At a confinement pressure of 2.41 kPa, he obtained an axial deformation of 13% and volumetric deformation equal to 5.5% at peak for the angular elongated material. This material is able to highly deform before reaching failure. On the other hand, the rounded material reached an axial deformation of 6.5% at failure and a volumetric deformation of 1.5%. This could be due to the ability of the rounded particles to have denser initial state than the angular particles where the peak of shear resistance is obtained at lower deformation levels. These previous remarks and conclusions have been also studied and verified to be true in most of the configurations by Barton [15] and Miural [START_REF] Miura | Physical characteristics of sands with different primary properties[END_REF] (Figure I.20) with some exceptions at certain porosity values [START_REF] Chen | Modélisation du comportement mécanique des grands CFRD: Identification des caractéristiques des enrochements et comportement du masque d'étanchéité amont[END_REF]. In Figure I.20, it is noticed how the secant strain modulus decreases with the increase of particle angularity.

Fig. I.20. Correlation between particle angularity (A 2D

) and axial strain (ε ac ) and strain modulus (E ′ 50 ) [START_REF] Miura | Physical characteristics of sands with different primary properties[END_REF] Laigle illustrated the effect of the particle angularity and state throughout triaxial tests [START_REF] Laigle | Proposition de modèles de comportement pour les enrochements et justification sur le cas du barrage de cogoti[END_REF]. In Figure I.21, an angular rockfill with high porosity has a high contractive behavior accompanied with a critical state that is reached at higher deformations.. On the other hand, an angular or rounded rockfill with lower porosity has a resistance peak with smaller deformations and a higher rate of dilation. It is also noticed that the angular particles gave higher shear resistance peak than the rounded particles which was expected. In addition to the angularity, grading represented by the coefficient of uniformity has a remarkable effect on the rockfill compressibility and deformability. Figure I.22 gives an estimate of the minimum and maximum void ratios as a function of the uniformity coefficient and the block angularity [START_REF] Ropers | Contribution à l étude du compactage[END_REF][START_REF] Biarez | Influence de la granulométrie et de son évolution par ruptures de grains sur le comportement mécanique de matériaux granulaires[END_REF]. We can see that as the particle angularity increases its minimum and maximum void ratio increase. Rounded material can reach denser states than the angular ones: angular material can reach looser state than the rounded ones. The increase of the uniformity coefficient leads to a decrease in the maximum and minimum void ratios. A well graded rockfill is able to reach lower void ratios and denser state than a poorly graded material.

E = K b .p a ( σ 3 p a ) m (I.9)
where K b = 4/n 3 , p a = 100 kPa is the atmospheric pressure, σ 3 is the confining pressure and E is the Young's modulus.

According to Equation III.5, rockfill with a porosity equal to 0.2 (typical porosity in dry stone retaining walls or in the pitching of rockfill dams) has a Young's modulus which is 11 times higher than a rockfill with a porosity equal to 0.45 (case of rockfill dumped in the dam body).

Cyclic behavior of rockfill

Earthquakes and seismic motion are major factors that are considered in the design and stability study of rockfill dams. To understand the rockfill dams behavior under seismic loading, the cyclic behavior of rockfill material should be reviewed and understood.

No cyclic tests have been performed directly on rockfill but the case of sands can illustrate what is expected in that case [START_REF] Leslie | Compaction of sands by repeated shear straining[END_REF]. In Figure I.23, cyclic deformations induces irreversible densification which rate is high for the first cycles and which tends to harden quickly with cycles. [START_REF] Bolton Seed | Seismic design of concrete faced rockfill dams[END_REF][START_REF] Gazetas | Seismic analysis and design of rockfill dams: state-of-the-art[END_REF]. The rockfill is characterized by a higher secant modulus degradation with cyclic deformation than what is observed for sand [START_REF] Bolton Seed | Moduli and damping factors for dynamic analyses of cohesionless soils[END_REF]. Accordingly, higher damping than for sand is reached for rockfill at high strains. The maximum initial shear modulus G max depends on the state of the rockfill and on its average confining pressure. It can be estimated using several relationships from the literature. For example, Seed [START_REF] Bolton Seed | Moduli and damping factors for dynamic analyses of cohesionless soils[END_REF] proposed equation I.11 which gives the elastic shear modulus G max of a granular material as a function of the confining pressure σ ′ m and parameter k 2 which depends on the soil state.

G max = 220.K 2 .(σ ′ m ) 0.5 (I.11)
In addition, Kokusho [START_REF] Kokusho | Cyclic triaxial test on sands and coarse materials[END_REF] proposed Equation I.12 for rockfill with elastic shear modulus depending on the void ratio e and the average confining pressure σ ′ m .

G max = 13000 (2.17e) 2 1 + e .(σ ′ m ) 0.55 (I.12)

PEDRA experiments on scaled down rockfill dam

The lack of a strong database to validate numerical models for rockfill dams with dry stone pitching, such as the one produced by Deluzarche [START_REF] Deluzarche | Discrete modelling of rock-fill: Application to dams; modelisation discrete des enrochements: Application aux barrages[END_REF][1] motivated a specific experimental campaign. Scaled-down dams were built by dry stone masons and tilted. This test also denoted pseudo-static test is a simplified method for determining the seismic resistance of structures in which the horizontal inertial forces caused by a seismic motion to the structures are approximated by constant forces (e.g., European standard EUROCODE 8).

Studied cases

The physical models' shape and dimensions are inspired by a particular dam operated by EDF, with a down-scale of 1/10. The models were 2 m high, 4.2 m in length and 2.25 m wide with a 45°s lope on both faces (Figure I.26). These models were created within a crane truck bucket. The pseudo-static tests were conducted by rotating the bucket with the truck crane until dam failure, which was defined as the collapse of both the pitching and the upper part of the dam body in the bucket. The bucket was rotated by elevating one side of the bucket with the truck crane's arm while keeping the other side of the bucket on the ground floor I.27a). The bucket's location upon failure is determined, and the angle formed by the bucket floor with the horizontal is calculated. The rotation of the bucket increases the destabilizing forces (horizontal forces) on the pitching and dam body (Figure I.27b). Then, a larger rotation angle of the bucket at failure indicates that the investigated structure can withstand more destabilizing forces. It qualifies as a structure that is more resistant to the seismic motion. Steel plates were welded to the bucket's two internal sides and then coated with a polyane sheet to reduce friction between the structure of the dam and the side walls. Two steel rods were connected to the truck's base to support the first row of stone pitching at both face toes (upstream and downstream faces). I.2). Furthermore, an experiment using a scaled-down rockfill dam without pitching (1 m in height) was carried out to validate the importance of dry-stone pitching in dam stability.

To track dam displacements, six sensors were installed on the downstream dry stone pitching. (Figure I.26). Furthermore, photogrammetry with four cameras and a 3D laser scanner was employed to locate stone deformations throughout the tilting process. In addition, a high-speed camera was also used to better track the pitching kinematics toward failure. Table I.3 shows the obtained tilting angles at failure for the four cases. Herein, the displacements and pitching deformations of the sensors are only presented for case 4 (Figure I.30), yet the results are representative of those found in the other cases (Experimental Results for cases 1,2, and 3 are found in Appendix A). Despite the careful design of the system, failure did not occur in a perfect plane strain mode due to movement constraints at the lateral end sides. Similar characteristics were identified in scaled-down experiments involv-ing dry stone retaining walls, however, repeatability testing revealed that they have no substantial influence in the quantitative parts of the problem [START_REF] Savalle | Experimental and numerical studies on scaled-down dry-joint retaining walls: Pseudo-static approach to quantify the resistance of a dry-joint brick retaining wall[END_REF].

Main results

By comparing the results of the four studied cases, several conclusions can be drawn:

• Due to the equal pitching weight in cases 1 and 2, it can be noticed that a larger friction angle (29°for granite versus 26.5°for schist) results in higher dam resistance;

• The comparison of cases 1 and 4 (the pitching is composed of a layer of granite in case 1 and two layers in case 4) demonstrates that pitching weight is a crucial parameter influencing dam resistance;

• Case 3 demonstrates how pitching anchoring in the backfill affects dam resistance. 1. The central area (in white) consists of the assembly of approximately square domains, composed of about 700 blocks, already in equilibrium before being installed to create the dam. This area is modeled by the duplication of this small square sample representing a loose state.

Former numerical modelings

3. The hand-laid pitching (in black) is composed of an assembly of blocks of type n°3 of 50 cm long, regularly laid on the external limits of the dam. For the downstream pitching, the hand-placed sub-layer on which was placed the pitching was considered as an integral part of the pitching, by modeling a double-row of blocks for this latter.

Once generated, the stability of the dam was evaluated. The role of the pitching on the stability of the dam was evaluated: first, the dam (with unbreakable blocks) was modeled without the pitching, which led to the complete failure of the structure (Figure I.33). The dam could not be stable without the upstream and downstream pitching. In a second step, the dam with pitching was modeled, and a stable state for the dam at the end of its construction was obtained which showed the key role of the dry stone pitching on the stability of the structure. Then, the filling of the reservoir was simulated under two conditions, when the blocks are unbreakable and then breakable. When considering breakable blocks, the behavior was different: about 6% of the blocks broke during the filling of the reservoir, resulting in displacements three times higher than those measured with unbreakable blocks. These displacements were mainly horizontal but were in this case partly irreversible: they did not disappear completely when the reservoir was emptied (Figure I.36). These additional displacements also resulted in a slip in the lower zone of the downstream face, linked to the presence of the berm discontinuity. This phenomenon was also observed in actual dams, although it is generally less significant. However, this slippage did not seem to affect the overall stability of the dam; thus, it seemed that the stabilizing role of the pitching is maintained.

A study has been done to evaluate the static stability level of the dam by reducing the friction coefficient between the blocks. In a first step, the coefficient was decreased from its initial value 0.7 down to 0.5: the dam remained stable. This reduction corresponds to a rockfill static safety factor (= tan(initial rockfill friction angle)/ tan(reduced rockfill friction angle)) of 1.85 and a pitching static safety factor (= initial pitching stones friction coefficient)/ reduced pitching stones friction coefficient)) of 1.4 as found by Deluzarche in his procedure [START_REF] Deluzarche | Modélisation discrète des enrochements: Application aux barrages[END_REF]. In a second step, the coefficient was decreased to a value of 0.1 which caused the failure of the dam. This corresponds to a static safety factor of 2.8 for the rockfill and 7 for the pitching. Therefore, the rockfill static safety factor was found between 1.85 and 2.8, and the pitching static safety factor was found between 1.4 and 7. Consequently, in this work, the range of the static safety factors is wide and does not indicate an accurate value.

Another approach was also used to study the stability level of the rockfill dam in which pseudostatic tests were carried out. Rotation was rather imposed on gravity (which is equivalent) and results were just saved for two values of the tilting rotation, due to the great amount of information to store: 5°and then 10°. The dam remained stable for a rotation of 5°and failed for a rotation of 10°. The seismic safety factor (against seismic motion) of the dam was then estimated between 1.2 and 1.4.

Dynamic computation on Escoubous dam has been also carried out by Deluzarche using a fully DEM model. In a first computation, the reservoir is supposed to be empty. Three values for the PGA (peak ground acceleration) of the applied signal have been tested: 0.1g, 0.3g, and 0.5g. The dam remained stable in those three studied cases which shows that such dam can resist to seismic forces up to 0.5g. The maximum displacement at the crest is 9 cm with a PGA of 0.1g, 16 cm with a PGA of 0.3g and 23 cm with a PGAof 0.5g. It was found that the post-seismic displacements are mainly associated with the blocks crushing. Without the possibility for crushing, the post-seismic displacements are found 7 times smaller.

With a full reservoir, the post-seismic horizontal deformations of the dam were found to be higher.

There is a global downstream deformation of the dam body. Moreover, in this case (full reservoir), the settlements were found lower than when the reservoir is empty. This is because the upstream face is confined by the hydrostatic pressure which prevents the deformations on this side and decreases the crest settlement.

The DEM model of Deluzarche had several limitations: first of all, he used a 2D model (excluding 3D behavior and effects) which in essence tends to overestimate dilative deformations (blocks forced to move in the plane). Moreover, this model was complex to build because of the creation of boxes full of particles, of different and sometimes complicated shapes, which were assembled to constitute the body of the dam. Furthermore, the contact of the pitching on the body of the dam with sufficient contact points between the two entities remained an issue that was difficult to handle. Finally, the number of circular particles used to build the body of the dam was important, which induced long computation times that nowadays would probably be reduced. The dumped rockfill forming the body of the dam was modeled by LK-Enroch constitutive model (noted LKE, detailed description of this model is found in Appendix C), developed by EDF, while the upstream mask was modeled by an elastic model and all other parts of the model using a Mohr-Coulomb model. The used parameters were extrapolated from previous studies on similar structures, in particular Gréziolles dam [START_REF] Laigle | Réévaluation de la stabilité du barrage de gréziolles : étude du comportement mécanique sous sollicitations statique et sismique[END_REF][START_REF] Silvestre | Development of a constitutive model for rockfill application for cfrd stability[END_REF].

The filling of the dam with water was simulated and the vertical and horizontal displacements were computed. They are respectively given in ) and was able to retrieve the damages observed on the concrete mask in a peculiar actual. He also studied the cracks found on the upstream concrete mask in which his numerical model using LK Enroch was able to reproduce them qualitatively. This also proved the ability of LKE model to well reproduce the rockfill behavior. For Escoubous dam, the safety factor of the dam was also estimated using the parameters degradation method. The displacements of different points at the downstream face were traced as a function of the reduction factors. The safety factor corresponded to the margin existing from the case where the slope stability of one of the dam slopes was not warranty anymore. A static safety factor equal to 1.5 was deduced.

Although this fully continuum model was on the overall relevant and acceptable i terms of computation time, it showed some limitations with respect to the local behavior of the dry-stone pitching. In that work, the pitching is assumed to be isotropic and it is modeled together with the hand-placed rockfill sub-layer as a single continuum element. In addition, no interface was used between the stone pitching and the dam body which prevents possible differential displacements and sliding between them. Finally, the local behavior of the pitching cannot be correctly reproduced due to the missing of its discrete nature. This fact limits the study of the high deformations and buckling (existence of bellies) in the downstream pitching that were described earlier.

Dry-stone retaining walls: Fully DEM vs mixed FDM-DEM approach

Experimental and numerical studies performed on Dry Stone Retaining Walls have helped to calibrate experiments for scaled down dams and motivated the choice of the numerical approach that will be taken for this PhD thesis. *

To better understand their mechanical behavior, dry-stone structures have undergone a number of experimental campaigns in the past. The earliest known research on dry stone retaining walls (DSRWs) were conducted in the nineteenth century by Anon [START_REF] Anon | Experiments carried on at chatham by the late lieutenant hope, royal engineers, on the pressure of earth against revetments, and the best form of retaining walls[END_REF] in 1845 and Burgoyne [START_REF] Burgoyne | Revetments of retaining walls[END_REF] [START_REF] Mundell | Behaviour of drystone retaining structures[END_REF], [START_REF] Mundell | Limit-equilibrium assessment of drystone retaining structures[END_REF] carried out large-scale studies in Great Britain to determine the cause of the bulging phenomena that is observed in situ and contributes to the failure of the majority of ancient DSRWs. In addition, scaled-down pseudostatic experiments were performed by Savalle (2020) [START_REF] Savalle | Experimental and numerical studies on scaled-down dry-joint retaining walls: Pseudo-static approach to quantify the resistance of a dry-joint brick retaining wall[END_REF] to investigate the seismic behavior of DSRWs. He studied the influence of the wall slenderness on the wall's failure tilting angle and on the failure mode. The failure angle decreases with more slender walls. In addition, sliding failure mode was obtained with less slender walls whereas more slender walls failed by toppling.

All those experiments on dry-stone structures allowed to constitute an extensive database for the validation of analytical design tools including the homogenized yield theory (Colas et al. 2018 [START_REF] Terrade | Upper bound limit analysis of masonry retaining walls using piv velocity fields[END_REF]). It also allowed to validate several numerical models including fully 2D DEM model [START_REF] Oetomo | Modeling the two-dimensional failure of dry-stone retaining wall[END_REF] and a mixed continuum DEM approach [START_REF] Deluzarche | Modélisation discrète des enrochements: Application aux barrages[END_REF][START_REF] James | Modeling the 2d behavior of dry-stone retaining walls by a fully discrete element method[END_REF]. For example, Savalle [START_REF] Savalle | Experimental and numerical studies on scaled-down dry-joint retaining walls: Pseudo-static approach to quantify the resistance of a dry-joint brick retaining wall[END_REF], [START_REF] Savalle | Seismic behaviour of dry stone retaining walls: Experimental and numerical pseudo-static studies[END_REF] validated a numerical model by modeling down-scaled tilting tests. Past studies proved that a mixed continuum -DEM approach (continuum approach to model the backfill, DEM approach for the retaining wall) was able to retrieve results obtained throughout experiments at a much more efficient computation time than a fully DEM approach. In 2013, Oetomo modeled a dry-stone retaining wall that was studied experimentally by Colas [START_REF] Colas | Assessing the two-dimensional behaviour of drystone retaining walls by full-scale experiments and yield design simulation[END_REF] using both approaches [START_REF] Oetomo | Modeling the two-dimensional failure of dry-stone retaining wall[END_REF][START_REF] James | Modeling the 2d behavior of dry-stone retaining walls by a fully discrete element method[END_REF][START_REF] Oetomo | Comportement des murs de soutènement en pierre sèche: une modélisation par approche discrète[END_REF]: fully DEM and mixed FDM-DEM (Figure I.41). The height of the backfill was increased gradually until failure was obtained (The backfill height at failure is called the critical height). Both approaches gave excellent precision in retrieving the mode of failure and the critical backfill height (Table I.4) with an advantage for the mixed FDM-DEM approach that needed much lower computational time than the fully DEM approach. However, these authors did not perform 3D computations. Considering the numerical studies performed in the past on Escoubous dam and on dry-stone retaining walls, the following key points can be pointed out

• A global safety factor in the range of 1.3 to 1.5 was found for Escoubous dam whether using a static or pseudo-static approach;

• Although the fully DEM model described well the discrete nature of the materials, the approach was complicated to process and results were obtained at the time with high computational cost;

• The fully continuum model created by Laigle in EDF showed that it can predict the dam's displacements after water filling using the calibrated elasto-plastic LK-Enroch model, but had difficulties to take into account precisely the discrete nature of the pitching and mobility between its blocks;

• Previous studies on a dry-stone retaining wall have shown that the mixed FDM-DEM approach can successfully reproduce the model behavior with results very close to the ones obtained thoughout experiments. This mixed approach can be considered as a compromise between the fully DEM which is time consuming and the mixed FDM-DEM approach: it produces accurate results and allows to take into account in a more rigorous way the discrete nature of some key elements of the dam, namely the pitching while being satisfactory in terms of computational time.

Conclusion

Rockfill dams with dry-stone pitching became popular in France in the 20th century but this technology was abandoned a decade after World War II. Such dams are still operated in France. They generally have steep slopes contrary to the actual designs which are much milder. They are usually made of rockfill which can be dumped or hand-placed. This rockfill is protected by a stone pitching which is hand-laid. A transition sub-base layer is found at the rockfill face under the pitching to help lay the pitching blocks. It has been shown in several previous studies that this pitching plays a major role in stabilizing the steep slopes of those dams. Those dams have been subject to several damages throughout their operation period such as bulging and cracks.

In the past, a fully 2D DEM approach was used to model the behavior of one of these dams. The obtained results were found comparable with a plain strain continuum approach but with a high computation cost. More details could be obtained in that case for the pitching which is a key element of these types of dams. Numerical studies involving dry-stone retaining structures which also involve a backfill with a wall made of individual blocks showed that another approach can be used to address the behavior of such structures in a more efficient way. It consists of describing the backfill/infill as a continuum material and the wall/pitching as a discrete system composed of individual blocks.

Herein, a mixed FDM-DEM approach is proposed to evaluate the mechanical behavior of rockfill dams with a dry-stone pitching with a focus on one of them: Escoubous dam.

Chapter II

A mixed FDM-DEM approach for rockfill dams with dry-stone pitching

Introduction

The suggested mixed FDM-DEM approach is used to model the scaled-down physical models of rockfill dams with dry-stone pitching (PEDRA experimental campaign's cases) described in paragraph 4 in chapter I.

The main objective is to develop and validate a mixed FDM-DEM numerical approach to study the mechanical behavior of such structures. The second objective is to study the effect of stone pitching characteristics (kind of material, density, assemblage) on the resistance of the dam. Finally, the effect of geometrical imperfections in the blocks on the global dam resistance and on the pitching behavior is studied. Indeed, real stones forming the pitching do not have perfect parallelepiped geometry which decreases the surface of contact between the blocks. This could have an effect on the pitching resistance and its deformability.

Numerical model 2.1. Model description

Dams such as those built in PEDRA experimental campaigns are composed of two distinct elements: the dam body represented by an assembly of boulders/gravel (the rockfill) and an arranged and well-organized assembly of stones (stone pitching) forming the skin of the dam. It is therefore easy to differentiate these two types of elements: the dam body (rockfill) is modeled as a continuous medium using FLAC3D (FDM) while the dry-stone pitching is modeled as a discrete system using PFC3D (DEM).

On the one hand, the engineering mechanics computing software FLAC3D uses a three-dimensional explicit Lagrangian finite-volume (hence referred to as FDM) model where materials are repre-sented by polyhedral elements (tetrahedral or hexahedral) on a three-dimensional continuous grid.

On the other hand, PFC3D is a distinct-element modeling software in which the medium is represented by independent bodies (balls, clumps, rigid blocks and walls) that interact with one another according to a built-in contact law. Both commercial software are developed by ITASCA [START_REF] Cg Itasca | Pfc-2d, particle flow code in 2 dimensions[END_REF] and may easily be communicated to each other to share pieces of information: information is shared at the interface between the backfill and the pitching. The description of the theory behind those software is presented in Appendix B. This mixed FDM-DEM approach is then used to simulate the behavior of PEDRA experiments.

The blocks forming the pitching have the same geometric characteristics as those used in the PEDRA experiments (Figure II.1). The rockfill represented by a continuous medium (mesh) is studied in section 2.4. For this system, it is necessary to study the influence of the mesh size to guarantee the validity of the numerical results.

Boundary conditions

In order to replicate the two lateral sides of the bucket, which are two smooth plates covered with a polyane sheet (to reduce friction), the normal displacement (y-direction in Figure II.3) of the backfill vertical planes was fixed, and a no-friction condition was imposed. A frictional interface was used to link the dam base to its foundation, and the dam foundation was fixed in all

Modeling procedure

The calculation was divided up into a number of subsequent steps. The dam body geometry and the mesh for the backfill were created within FLAC3D as the initial phase. Gravity was applied until equilibrium was established: for such structures, equilibrium is meant to be reached when the unbalanced forces ratio is less than 1.10 -5 (value usually used in static stability in FLAC3D).

It is defined as the average unbalanced mechanical force magnitude for all gridpoints in the model (or at contacts for infinitely rigid bodies as stated herein in PFC3D) divided by the average applied mechanical force magnitude for all gridpoints (or on a body in PFC3D). The dry-stone pitching was generated within PFC3D on the dam's upstream and downstream faces. At this point, the coupling between the two softwares is processed by creating an interface wall between the stone pitching and the dam body. To reach a new mechanical equilibrium, the whole dam model is once more run under gravity. The following section 3 discusses all mechanical parameters of the dam as well as the pitching-backfill interface.

Mesh generation and study

The dam body mesh (backfill) was developed by investigating the influence of mesh size on the failure angle and maximum pitching displacement at a tilting angle of 20°. The major goal was to identify the biggest possible mesh size (lowest gridpoint number) that had no effect on the outcomes. The evolution of the failure angle and maximum displacement as a function of the gridpoint number on the backfill face in contact with the stone pitching is displayed in II.5. Both curves reach a "plateau" after a specific number of high gridpoints, as indicated by the red dotted line. The mesh corresponding to this threshold number was chosen as an optimum mesh size that ensured mesh independence of outcomes while also limiting computation time. In Case 1, for example, each stone is in contact with at least one gridpoint on the backfill surface. As a result, the total number of gridpoints on the backfill surface that come into contact with the pitching is more than the entire number of stones that make up the pitching. In all four modeled cases, the final resulting mesh for the dam body FLAC3D model is the same.

Equilibrium and failure criteria

The dynamics equations are solved using an explicit time technique in FLAC3D and PFC3D.

A criterion must be specified in order to determine when a static equilibrium is reached. The average unbalanced force ratio is generally monitored. The targeted ratio is specified by the user, however ITASCA recommends a common value of 10 -5 when accurate results are expected. This criterion, however, may be relatively long to achieve and even never be reached indicating that a failure is in fact on the way. As a result, another criterion was created to determine whether failure occurred prior to a fully visible collapse and the fully stability of the entire system. In case of a detected failure, this additional criterion should enable the simulation to be stopped, thus saving computational time. In such a quasi-static system, monitoring the kinetic energy of the system or a subset of the system can supplement the preceding criterion: if it exceeds a user-defined threshold, a continuing collapse is validated, and the computation can be stopped.

The developed criteria were comparable to those applied by Oetomo or Savalle in their respective works [START_REF] Oetomo | Modeling the two-dimensional failure of dry-stone retaining wall[END_REF][START_REF] James | Modeling the 2d behavior of dry-stone retaining walls by a fully discrete element method[END_REF][START_REF] Oetomo | Comportement des murs de soutènement en pierre sèche: une modélisation par approche discrète[END_REF][START_REF] Savalle | Dem modelling of dry stone retaining walls[END_REF]. The model is run (gridpoint locations and stone placements are changed) through thousands of iterations for a particular tilting angle. If the model reaches the threshold unbalanced force ratio during these iterations, the model is considered to be in equilibrium, and another increment of tilting angle can be processed. Otherwise, a second criterion based on pitching kinetic energy is used to determine whether the dam has reached failure. Failure is assumed to have occurred if the kinetic energy is greater than 1.10 -3 J for 10 consecutive observations with 500 iterations between each. If not, it suggests that the model may tend to a stable state and that the entire system has a further chance of achieving equilibrium (few thousands more iterations). During a preliminary investigation in which pitching kinetic energy was tracked, the kinetic energy threshold value was determined. It means that the calculated threshold value is only applicable to the studied system. Figure II.6 shows the full scheme.

Model calibration

In this study, a sophisticated constitutive model known as LK-Enroch (LKE) developed at EDF for rockfill is used. Additionally, this advanced model is compared with a more typical constitutive model Mohr-Coulomb (MC). Data from the PEDRA experimental campaign was used to calibrate the model parameters. When information was unavailable, calibration was carried out using data from the literature. This section describes in detail the calibration process for each component of the system, including backfill, stone pitching, and the backfill-pitching interfaces.

Backfill

The dam body is composed of a backfill made of angular granite rockfill with a grading ranging from 10 mm to 80 mm with a uniformity coefficient of 2 and a D 50 of 40 mm (Figure II.7 ).

Porosity is equal to 0.45.

Mohr-Coulomb

Experimental triaxial tests conducted for confining pressures of 100 kPa and 200 kPa were used to calibrate the model parameters for the material forming the dam body. For the granular material under consideration, a Representative Elementary Volume (REV) was obtained using a 1-m diameter cell. The mean pressures present in the scaled-down system, which range from 1 to 10 kPa at mid-height from the dam crest to the base, are not comparable to the confining pressures that are typically used in the experiments. It means that there will be a two-stage calibration process for the model parameters. The first one will include calibrating the model parameters using the results of the two available triaxial tests. In a second stage, the set will be adapted to be valid for the range of mean pressures existing in the scaled-down system.

In the first stage for the identification of the model parameters being completed (Table II.1), the internal friction angles at 100kPa and 200kPa throughout the triaxial tests using a large cell were found equal to 42°and 39.5°respectively. According to LEP's chart (Figure II.8), these values are compatible with a rockfill at low density.

In the second stage, these values were adapted to low confining pressure, which corresponds to the average mean pressure in the PEDRA dam, using LEP's chart [START_REF] Thomas | Review of shearing strength of rockfill[END_REF]. The internal friction angle of the rockfill was estimated to be 52°considering a low density for a poorly graded material and an average pressure ranging from 1 to 10 kPa (Table II.1). Bolton's approach was used to calculate the dilatancy angle [START_REF] Md Bolton | The strength and dilatancy of sands[END_REF]. According to the author, in plane strain conditions and a relative dilatancy index of 0 < I R < 4 (I R = I D (10-ln(p')) -1; where I D is the relative density and p' is the confining pressure), the dilatancy angle ψ max can be determined using Equation II.1, which matches the rockfill features described herein. In Equation II.1, φ peak represents the backfill internal friction angle, while φ crit represents the backfill friction angle at the critical state.

First, the dilatancy angle for the rockfill and for confining pressures of 100 and 200 kPa was found equal to 6°and 2°, respectively (Table II.1). For the range of confining pressures of the studied scaled-down system (1-10 kPa), Equation II.1 gave a dilatancy angle close to 18°. Finally, based on Equation II.2 [START_REF] Hicher | Micro-macro correlations for granular media. application to the modelling of sands[END_REF], which considers the typical dependency of the Young's modulus on the mean pressure and the reference pressure P re f of 100 kPa [START_REF] Hicher | Micro-macro correlations for granular media. application to the modelling of sands[END_REF], E 0 was determined to be 50 MPa throughout the triaxial tests with the large cell. The Young's modulus was then calculated at P=5 kPa (average pressure in the scaled-down numerical model computed at mid-height and at the dam's center). Finally, Poisson's ratio ν was maintained as independent of mean pressure, and cohesion C was set to zero. Table II.1 gives all of the values for the MC model. The other model parameters were identified throughout a trial-and-error technique with the help of simulations of triaxial tests. Then, considering available actual triaxial tests at 100 kPa and 200 kPa, both elasticity threshold parameters a 0 and m 0 were identified.

E = E 0 [ P P re f ] 0.5 (II.2)
The peak threshold parameters which mainly control the value of the deviatoric stress peak are namely a peak , m peak and σ c . First, the compressive strength σ c is assumed to be equal to 100 MPa [START_REF] Ozturk | Strength classification of rock material based on textural properties[END_REF], then a peak and m peak values were calibrated (Table II.3) to reach the experimental deviatoric stress peak (Table II.2).

Concerning the hardening parameters, γ˙p eak was chosen by fitting the position of the deviatoric stress peak with respect to the axial strain. γ˙r es and η were chosen to fit the post-peak softening phase whereas X ams was chosen to fit the curvature of the pre-peak curve.

The volumetric deformations parameters were set with the help of the triaxial experimental volumetric curves. First, the residual friction angle φ res equal to 38°was deduced from the actual triaxial tests. Then, the dilatancy angle ψ 0 was chosen by fitting the volumetric curves.

Finally, the isotropic mechanism parameters p c0 and β were estimated using Hicher and Rahma approach [START_REF] Hicher | Micro-macro correlations for granular media. application to the modelling of sands[END_REF] (the main equations of this approach are given in Appendix D). While p c0 was identified equal to 1 MPa throughout the simulations of triaxial tests, it was taken equal to 15 kPa for the modeling of the scaled-down dam (average mean pressure of 5 kPa). The shape form parameter γ that controls the deviatoric failure envelope dissymetry is typically taken equal to 0.87. The final values for LKE parameters are listed in Table II The LKE volumetric strain curve for 100 kPa is quite close to the one produced throughout the experiments which is not exactly the case for 200 kPa. Triaxial repeatability experiments (not illustrated herein) revealed that the volumetric deformations curve is sensitive to flaws in the tested samples: replicating the observed trends contributes to confirm the quality of the calibration obtained.

Sophisticated elastoplastic models are often valid for the mean pressure range involved with the calibration process (normally hundreds of kPa). The mean pressure in the scaled-down model was substantially different (1-10 kPa), implying that the one found may not be appropriate for the examined scaled-down system. An approach similar to that utilized for identifying MC model parameters was used. The adaptation is intended to improve the prediction of peak strength in particular. The peak threshold parameters a peak and m peak were then adjusted to achieve a target value for the internal friction angle of 52°(LEP's approach [START_REF] Thomas | Review of shearing strength of rockfill[END_REF]) as used for MC model calibration.

Dam-foundation interface

A linear frictional Coulomb model was used to model the dam-foundation interface. A friction angle of 52°(equivalent to the internal friction angle of the rockfill) is applied. The normal and shear stiffnesses K n and K s were defined using the FLAC3D manual's suggestions. Stiffnesses are advised to be ten times the equivalent stiffness K eq of the most rigid zone near the interface defined by EquationII.3, where K and G are the bulk and shear moduli, respectively, and ∆z min is the minimum dimension (height) of volume element meshes at the interface. Then, the values were identified to K n = K s ≈ 1.10 10 Pa/m.

K eq = max[ (K + 4 3 G) ∆z min ] (II.3)

Stone pitching

Using PFC3D, the stones forming the pitching are modeled as infinitely rigid bodies with deformable contacts. The description of PFC3D software is found in appendix B. For contact between blocks, a linear frictional law is used with multiple variables such as the local friction coefficient, shear stiffness, and normal stiffness (Appendix B). The shear and normal stiffness inputs can be substituted by the effective modulus E * and stiffness ratio K * , where K * = K n /K s , with K n signifying normal stiffness and K s denoting shear stiffness.

For the pitching, whether granite or schist, the friction coefficient µ (tangent of the friction angle of the stone-stone contact) is used to define the contact's friction. The granite stones in cases 1 and 4 have µ granite = tan(29°)= 0.55, and the schist stones in cases 2 and 3 have µ schist = tan(26.5°)= 0.5. In the PEDRA campaign, the local friction angles of granite and schist stones were measured experimentally. Slide slip tests were performed by stacking two stones on top of each other and rotating the entire system until the top block slipped on the bottom one. For the granite and schist stones, an average value was calculated based on many experiments (Figure II.12).

The effective modulus and stiffness ratio used are 5.10 7 N/m 2 and 2, respectively. They were obtained as a compromise between actual stiffness and computation times through a parametric study. High stiffness, in fact, causes a small time increment for the resolution of the dynamics equations, penalizing time computation. However, the contact stiffness can be reduced without affecting the quality of the results. A global damping of 0.7, which is commonly utilized in the context of quasi-static calculations was considered [START_REF] Peter | A discrete numerical model for granular assemblies[END_REF].

Pitching-backfill interface

The pitching-backfill interface is activated by the coupling PFC3D and FLAC3D software. A coupling mechanism is processed on the downstream and upstream faces within a wall created in PFC3D whose facets coincide with the facets of the dam body in FLAC3D (Figure II.13). The coupling logic works by taking the contact forces and moments from the infinitely rigid bodies in PFC3D (pitching stones) with wall facets and determining an equivalent force system at the facet vertices. These forces are passed to the gridpoints along with stiffness contributions to the facets of the dam body in FLAC3D.

A purely frictional linear elastic model, using the same set of parameters as the one used for the and4 is based on the stone-stone friction coefficient, according to Colas [START_REF] Colas | Upper-bound solution for the stability of stone-facing embankments[END_REF]. Thus, for cases 1 and 4 (granite stones), the contact friction coefficient was µ granite = tan( 29°)= 0.55, but for case 2, the interface friction coefficient was µ schist = tan(26.5°)= 0.5. In the peculiar case 3, the pitching was anchored in the backfill. This anchorage was not physically modeled; instead, it was assumed that the presence of anchors caused the shearing failure surface at the interface to shift into the backfill (Figure II.14). This phenomenon can be observed at a variety of rough soilstructure interfaces when the interface friction angle is assumed to be equal to the backfill friction angle. Different authors successfully used this statement when designing dry stone retaining walls with a configuration similar to that of this study (e.g . [START_REF] Savalle | Dynamic numerical simulations of dry-stone retaining walls: identification of the seismic behaviour factor[END_REF]). This statement describes a rough interface between the concrete structure and the surrounding soil that results from in-situ concrete casting, as stated in the European Code for Geotechnical Works (EUROCODE 7). Therefore, the friction angle of the rockfill is chosen for the interface friction angle of case 3 which is µ rock f ill = tan( 52°)= 1.28. The percentage of the difference between the numerical (MC and LKE) and experimental results are also calculated and shown in table II.4: very small error (less than 5 %) is observed in the predicted failure angles in most of the cases except for case 2 where we have an error of about 9 % when using LKE model which is still considered satisfactory. In fact, the physical models have several imperfections contrary to the perfectly plain strain numerical model. This leads to the conclusion that the numerical model developed herein gave an excellent prediction of the dam's resistance in the 4 cases. In most of the sensors locations, MC underestimated the displacements as a function of the tilting angles: this model is characterized by a large elastic domain which generates a rigid behavior.

On the contrary, LKE model is able to generate further (irreversible) deformations from small deformations due to a limited initial elastic domain that expands according to a non-linear isotropic evolution law. In consequence, LKE model was able to retrieve the evolution of the displacements with a very good accuracy: the displacement at each sensor increased gradually with the increase of the tilting angle which is the case in the experiments.

To illustrate the difference between MC and LKE prediction of the displacement, the displacement field of the backfill and of the dry-stone pitching for case 4 (MC and LKE models) for a tilting angle of 27°(just before failure) are given and compared to the corresponding experimental results in Figure II.17 parametric analysis was conducted on Case 1 taking into account 3 distinct interface friction angles: 30°, 40°, and 50°to better understand the influence of the interface friction angle. These three values produced failure positions of 0.37*H, 0.39*H, and 0.42*H, respectively. This study demonstrates a small effect of the interface friction angle on the failure position of the rockfill dams, but it does not account on its own for the significant variation seen in case 3. Instead, a second factor, such as the pitching block dimensions, can explain why the maximum relative displacement position shifts upward.

Safety factors

An important step in the process of justifying and designing any structure, including dams, is the evaluation of the safety factors. Two different safety factors based on two different approaches have been computed: the static safety factor and the pseudo-static safety factor [START_REF] Deluzarche | Modélisation discrète des enrochements: Application aux barrages[END_REF].

The first safety factor F φ is linked to the system's static resistance . The internal friction angle of the material forming up the dam body, for instance, must be at least equal to both the angles of the upstream and downstream slopes for the system to remain stable -for a given earth dam without mask or pitching on both the upstream and downstream sides. According to this framework, the first safety factor is calculated as the ratio of the actual internal friction angle of the dam body material and the reduced internal friction angle required to cause failure. This first factor is calculated using Equation II.4 by gradually decreasing the backfill friction angle φ until no stable condition can be identified.

F φ = tan(φ initial ) tan(φ reduced ) (II.4)
The safety factors were computed on the basis of the rockfill dam using LKE model. In order to calculate the reduced value of the friction angle φ reduced in Equation (II.4), three LKE parameters a peak , m peak and σ c allow the reduction of dam body resistance. Since a peak is the parameter that most affect the variation of the peak resistance, we chose to gradually reduce its value while maintaining constant m peak and σ c . Finally, from the reduced LKE set of parameters, the equivalent Mohr-Coulomb friction angle φ equivalent at failure is determined. It is determined by the trial and error method where the friction angle (using MC), giving the same triaxial test peak resistance as the LKE reduced parameter, is chosen as the reduced friction angle φ reduced .

The second safety factor F θ is associated with the system's resistance to horizontal inertial forces.

It is a common calculation used to examine a system's stability against earthquakes using a simplified pseudo-static approach (for example, slope stability). This safety factor is defined as the ratio of the maximum possible inclination of the dam slope (which is directly related to the maximum potential horizontal inertial force the system can support) to the dam slope at rest. Equation II.5 defines this factor, where α is the dam slope angle (here, 45°) and θ f ailure is the tilting angle at failure.

F θ = tan(α + θ f ailure ) tan(α) (II.5)
The safety factors were determined for the four cases examined. An additional calculation was performed in the situation where the pitching was removed. This is a reference test that will show how much pitching helps to improve dam stability. In the case of F φ , calculations were performed taking into account the pressure due to water in the reservoir (F φ (ww) , with water) on the pitching and without pressure (F φ (wow) , without water).

Table II.5 presents all of the calculated safety factors. Firstly, the results demonstrate that without pitching, the rockfill dam is very close to failure based on both safety parameters. Secondly, according to both safety parameters, pitching leads to a very good safety margin. They are multiplied by a minimum of two in comparison to the reference case (without pitching). Cases 3 and 4 result in higher levels of safety than Cases 1 and 2. This is attributed primarily to a heavier pitching weight. However, the weight of case 3 pitching is less than that of case 4 pitching, but the safety considerations for the former are greater. As expected, the backfill-pitching interface friction angle plays a crucial role: such a trait is also found true for retaining walls [57] [62].

Because of the higher resistance (contact friction angle) of granite blocks in the pitching compared to schist blocks in Case 2 (Table II.5), Case 1 is slightly more resistant than Case 2. However, no firm conclusion can be taken at this stage regarding the significance of its role because no computation was performed with a very different type of stone with a considerably greater or lower contact friction angle. Colas and al. [START_REF] Colas | Upper-bound solution for the stability of stone-facing embankments[END_REF] conducted research on the resistance of rockfill dams with dry stone pitching. The yield design method was used to determine failure in a parametric study on the block-block friction angle. According to their findings, increasing these parameters increases the dam resistance.

Finally, by comparing the static safety factors with and without water pressure, it is noticed that the reservoir water has little influence on the dam resistance: the failure is generally triggered in the one-third upper part of the dam where the hydrostatic pressure on the upstream face is very low. This phenomenon is shown in figure II.20 which presents the displacement increments, measured at the last stage before failure in the static safety factor calculation, for Case 1 with and without water pressure. Without water pressure, the dam displacements before failure are symmetrical with 2 failure surfaces starting from the crest and ending at the bottom part of both upstream and downstream pitchings. Conversely, with water pressure, we have only one failure surface on the downstream side which also starts from the upper part of the dam near the crest where the water pressure is low. On the upstream side, the water pressure plays a role in higher stability and lower displacement even with reduced parameters.

Influence of pitching block geometry

The initial numerical model described until now considers a pitching which is modeled using perfectly parallelepipedic blocks which is not the case in real situations. The stones forming the pitching hold a disparity in sizes and shapes inducing that the contacts between stones involve just a part of the possible contact surface area with the surrounding stones which changes the possibility for relative movements. This fact was the motivation for a numerical study of irregular blocks with imperfections on their faces. This parametric study is only carried out for case 1. After the construction of the 3 numerical models including geometrical imperfections in their pitching blocks, they were tilted until failure to study the effect of such defects on the dam resistance. The failure angles of all the studied imperfection cases are shown in table II.6. We can first notice that all the imperfection cases in the blocks lead to a decrease in the failure angle (dam resistance). In fact, this was expected since the pitching is now less rigid with more flexibility in the blocks movements with respect to each other. In addition, less contact surface is found between the blocks which decreases the frictional resistance. Furthermore, by comparing the failure angle obtained in the 3 imperfection cases, the case with imperfections on the bottom (internal) side of the blocks has the lowest failure angle (21°) compared to 23°obtained in the 2 other cases: the pitching mainly deforms towards the outside of the dam body which makes the pitching with blocks with imperfection on the bottom less resistant to such deformations.

The failure position of the downstream pitching obtained after the dam rotation is now compared between the 3 imperfection cases. It is measured as a percentage of the dam height (Table II.6).

Starting with the case of blocks with bottom imperfection, the obtained failure position is at 43% of the dam height. This value decreases to 35% for the blocks with bottom and upper imperfections then to about 30% for the blocks with just upper imperfections. We can notice a trend of the failure position decreasing when the imperfection passes from the bottom side to the upper side. The normal displacement of the downstream pitching at a tilting angle of 20°is shown in In the second step, since the real contact surface between the blocks is unknown, the influence of its contact surface length is studied: the initial contact surface (Coef parameter) fixed to 20% of the total block length is now increased to 30%, 40% and 60%. We can see that the failure angle decreases gradually from 25°for the full contact surface length case (1*d corresponding to regular and perfect parallelepipedic blocks) as the contact surface length decreases. A failure angle of 22°is found with a 30% contact surface length and reaches a value of 21°for a 20% contact surface length (Table II.7). The decrease in resistance is due to the reduction of the frictional resistance between the blocks as a result of the reduced contact surface. In addition, with a shorter contact surface length, the blocks have more rotational freedom that increases the pitching ability to deform locally and decreases its resistance. This is illustrated in In conclusion, the pitching mechanical behavior and resistance are mainly influenced by the quality and contact surface area between blocks in the sub-horizontal plane (contacts between one row of blocks and another): contacts between blocks along their vertical plane (vertical imperfection) have little influence on the resistance and deformability of the dam and its pitching.

Conclusion

The numerical modeling of rockfill dams with dry-stone pitching was validated in this chapter. It was based on the simulation of actual scaled-down pseudo-static tests in which different pitching features were assessed. Four cases from PEDRA campaign were analysed where differences lied in the kind of stones for the pitching, weight of pitching and the existence of anchored pitching blocks within the dam body. For the dam body, two constitutive models were used: the Mohr-Coulomb (MC) model and an advanced elasto-plastic model known as LK-Enroch (LKE).

Regardless of the dam body's constitutive law, the tilting angle leading to collapse was obtained with high precision throughout the PEDRA campaign's tests. However, because of its inherently stiff nature, the MC model was unable to retrieve the correct displacement field while the LKE model had successfully predicted the displacement field. It is due to the presence of a small initial elastic domain in the LKE model, which allows for the generation of irreversible deformations for small deviatoric stresses.

The relative displacements study of the pitching throughout simulations revealed that the backfillpitching interface friction angle played a significant role in the failure position. Furthermore, the effect of the pitching weight on the dam's rigidity was observed: increasing the dam weight tends to reduce pitching deformations (normalized by pitching thickness) at failure, resulting in stiffer dam behavior.

To assess the dams' safety factors, two approaches were used. Both approaches quantitatively justified the crucial role of the pitching in dam stability. The dam is close to failure without the pitching, with a safety factory slightly larger than 1. For all cases involving pitching, the safety factors increased to values greater than 2.5. Furthermore, by comparing the four cases, the pitching weight and interface friction angle were numerically justified to have an important impact on the dam's resistance. It has been also shown that the water pressure on the upstream side has negligible effect on the dam resistance since the failure is triggered from the upper part of the dam were the water pressure is very low.

The effect of the geometrical imperfections in the pitching stones on the global behavior of the dams was studied. Several geometrical imperfections have been introduced to the pitching of Case 1. The resistance of the dam decreases when imperfections are added in the sub-horizontal contact planes between blocks whereas it does not highly affected when added in the vertical contact planes. The case with the internal imperfections (bottom) has the highest effect on the dam's resistance. In addition, the imperfections lead to a higher ability to have displacement peaks and bulging effects.

Finally, it should be noted that the results and conclusions obtained in this study can be different from actual real-scale dams. Very different features and conditions can be at stake on actual structures. If the main qualitative conclusions related to the role of the pitching may hold true for actual dams, scale effect induced by far larger mean pressures in the dam body may greatly modify the quantities computed herein.

Chapter III

Static and pseudo-static behavior of Escoubous dam

Introduction

After the validation of the mixed FDM-DEM approach for modeling the mechanical behavior of scaled-down rockfill dams with dry stone pitching, a case study on an actual dam is proposed. The studied case is Escoubous dam that was previously studied by Deluzarche [START_REF] Deluzarche | Modélisation discrète des enrochements: Application aux barrages[END_REF][START_REF] Deluzarche | Discrete modelling of rock-fill: Application to dams; modelisation discrete des enrochements: Application aux barrages[END_REF] and Laigle [START_REF] Laigle | Réévaluation de la stabilité du barrage de gréziolles : étude du comportement mécanique sous sollicitations statique et sismique[END_REF]. This study aims is to position the mixed FDM-DEM approach compared to previous ones and to strengthen the robustness of the proposed approach as a standard tool for engineering practice.

The 18 m height dam is composed of a dumped loose rockfill forming its body, a hand-placed rockfill sub-layer on the two faces of the dam and a dry-stone pitching which is laid by skilled workers on the sub-layer. The sub-layer is used in such dams to form straight dam faces and then to facilitate the placement of the pitching stones. The dam rests on a moraine layer which itself lays on a rocky massif. A detailed description of the dam elements and dimensions is found in Chapter 1.

The dam body has been modeled as a continuum medium using LKE constitutive model in FLAC3D whereas the stone pitching has been modeled using DEM in PFC3D. LKE model parameters for the dam body were extrapolated from a previous study on a similar dam studied by EDF. More precisely, an inverse method has been used from the measurements of the displacements of Gréziolles dam throughout its first reservoir water filling. The hand-rearranged rockfill sub-layer on which the pitching lays and the foundation elements have been modeled using the Mohr-Coulomb model and their parameters were all obtained from the literature as average values for the involved materials.

First, the construction of the dam was modeled under several stages and then the dam reservoir was filled. The first objective of this study is to validate the model and its used parameters by comparing the obtained displacements to the historical and/or expected ones.

The second objective of this study is to evaluate different safety factors for this dam using two approaches: the first one, denoted static safety factor is related to the safety margin against slides of both downstream and upstream slopes. The second safety factor denoted as a seismic factor is related to the ability to resist horizontal inertia forces provided by earthquakes.

The third objective of this chapter is to understand and quantify the contribution of different elements of the dam to its overall stability. The study starts from the most simple possible dam model and then other elements are successively added until reaching the final model. Several supplementary cases are presented also.

The final objective is to study the large deformations that were observed on the dam downstream pitching during its operational period and to investigate their effect on the dam's stability. They include bulges found on the bottom part of the pitching and cracks under the berm. Blocks of the pitching with perfect geometry were used but another study was provided using geometrical imperfections for them. It may be closer to the actual block shape and provide the pitching a further ability to deform. Then, the presence of the chamber valve that exists at the bottom of the downstream pitching was also simulated.

Safety factors

The evaluation of safety factors is an important step in the justification and design process of any structure, including dams. For a stability study, it is necessary to quantify the safety margin against the risk of overall mechanical instability of the structure. Given the complexity of the rheological behavior of rockfill dams, which are heterogeneous and highly deformable, the used modeling methods do not allow us to derive directly a slip safety coefficient, usually defined by simplified failure analyses such as slip circles. The safety factor selected is the one for which calculations show a sudden increase in displacements and consequently in the kinetic energy of the pitching: the failure criteria described in section 2.5 in Chapter 2 have been used with an adaptation of the critical kinetic energy.

Static safety factor

The static stability of the dam is related to the stability of both the upstream and downstream slopes and was estimated by degrading the properties of the backfill and pitching, and determining under which property conditions the dam ceases to be stable. For the rockfill, modeled by the LK-Enroch law, the safety factor (F φ ) is defined as follows.

If σ c-max denotes the maximal strength of the material, the maximum deviator evaluated under triaxial loading is equal to

Q max = σ c-max m pic σ min σ c-max a pic
, m pic and a pic being the parameters of the maximum strength criterion and σ min the confining stress. For a given rock strength σ c , the associated deviator Q is then equal to Q = σ c m pic

σ min σ c a pic .
The safety factor (F φ ) is then defined by the ratio

Q max Q = σ c-max σ c 1-a pic .
In the logic of the concepts introduced in the LK-Enroch model, the residual friction angle φ res also evolves as a function of F φ . Indeed, physically, a degradation of material maximum strength is related to a degradation of intrinsic properties and then is bound to be reflected by a decrease in the residual friction angle. Numerically, the maximum resistance and residual friction angle are respectively calculated from the following formulas:

σ c = σ c-max F 1/(1-a pic ) φ (III.1) tan(φ res ) = tan(φ res-max ) F φ (III.2)
Concretely, a reduction in the safety margin of 20 % is associated with a reduction of the internal friction angle of the rockfill (backfill) of almost 5 % and effectively a reduction in the maximum deviator Q max of 20 %.

In parallel, the reduction in strength of the stone pitching is taken into account at the same time as that of the backfill, by modifying the angle of friction of the block-block contact φ P , based on the following expression:

tan(φ P ) = tan(φ P-max ) F φ (III.3)

Seismic safety factor

A second safety factor F θ is constructed and obtained by global rotation of the dam, similar to that used for pseudo-static tests in Chapter II. This safety factor is therefore associated with resistance to horizontal inertial forces. It is expressed in equation III.4 in which α is the slope of the studied dam face (equal to 45°in the case studied) and θ f ailure is the angle of rotation at failure [START_REF] Deluzarche | Modélisation discrète des enrochements: Application aux barrages[END_REF] (see section 3.2 for the definition of the failure criterion).

F θ = tan(α + θ f ailure ) tan(α) (III.4)

Numerical model

Geometry and boundary conditions

The dam body and the foundation were modeled as a continuous medium in FLAC3D, while the dry-stone pitching was modeled as a discrete system in PFC3D. The advantage of such an approach is that realistic results, as previously justified [START_REF] James | Modeling the 2d behavior of dry-stone retaining walls by a fully discrete element method[END_REF], can be achieved in less time than The lateral width of the model was set at 9 m, representing a total of 22 blocks (the same number as in the PEDRA models) on a lateral line of pitching blocks. This lateral number of blocks was chosen to provide a sufficiently wide pitching to allow an observation of the dam behavior under a plane strain deformation in the middle part of the system (lateral boundary conditions having very little impact on the dam's behavior) without penalizing computation time. The pitching blocks have dimensions of 20 × 30 × 40 cm (height × thickness × width): these dimensions were estimated from photographs of the pitching provided by EDF.

The mesh size of the backfill was chosen to ensure sufficient grid points in contact with the pitching blocks. Then, the mesh dimensions are of the order of 20 cm. To maintain the same mesh size ratios according to backfill height, the number of meshes decreases as approaching the crest (Figure III.1).

The number of meshes in the model is thus 120,960, and the total number of blocks forming the The following boundary conditions were applied:

• The two lateral sides of the modeled dam + foundation (vertical sides along the length of the dam) have been fixed in the direction normal to the vertical face (y-direction) and free to move in the other directions;

• At the level of the foundation on which the dam rests, the lateral edges of the model have been fixed in the horizontal direction (x-direction), free in the other directions. The base of the model has been fixed in the vertical direction and is free to move in the other directions;

• For convenience, crest blocks were not explicitly simulated, but a vertical stress equal to the surcharge produced by these blocks was applied at the crest, i.e. a stress of 6.24 kPa (crest block thickness × density × g = 0.3 m × 2,080 kg/m 3 × 10 m/s 2 );

• As far as the pitching is concerned, the 2 (down the slope) stone columns of the pitching located at the lateral limits of the system studied ("vertical walls") have been fixed in the direction normal to the wall (y-direction) and free to move in the other directions. The first row of stones in contact with the foundation soil has been completely fixed.

Failure criterion

The failure criteria used in this chapter concerning Escoubous dam are the same as the ones used in Chapter 2 for the scaled-down dams. The first criterion concerns the unbalanced forces ratio. The targeted value is equal to 10 -5 which is generally associated to a computation with a very good precision. After several thousands of cycles, if the first criterion is not met (10 -5 not yet reached), the other criterion concerning the pitching kinetic energy is then tested: if the pitching kinetic energy is greater than a critical value for 10 consecutive measurements, performed every 500 iterations, failure is assumed to have been reached. Otherwise, the model can tend towards a new stable state, and the whole system is given another chance to reach equilibrium (a few thousand more iterations). The critical kinetic energy value is adapted to the real-scale studied dam: a devoted study showed that a compromise between a very low value and a value corresponding to the loss of stability of the dam can be found using 10 J. This value was chosen by tracing and quantifying the pitching kinetic energy in a preliminary study and by comparing the weight of the real-scale pitching with the scaled-down pitching (the critical kinetic energy of the scaled-down model was 10 -3 J).

Model parameters 3.3.1. Dam body: rockfill

The dam body, which consists of simply dumped rockfill, has been modeled using the LK-Enroch (LKE) model. The model parameters are listed in Table III.1. These parameters were extrapolated from a previous study carried out by EDF on a similar structure with similar dumped rockfill characteristics: Gréziolles dam [START_REF] Laigle | Réévaluation de la stabilité du barrage de gréziolles : étude du comportement mécanique sous sollicitations statique et sismique[END_REF][START_REF] Silvestre | Development of a constitutive model for rockfill application for cfrd stability[END_REF]. This dam has been the subject of more detailed monitoring and experimental testing, which has enabled the parameters of the numerical model to be calibrated.

A triaxial test at a confining pressure of 100 kPa was performed using these LKE model parameters (Figure III.3). The observed volumetric behavior is purely of contractive type, consistent with the fact that the rockfill forming the dam body is loose. Furthermore, with regard to deviatoric stress behavior, we note that the peak is reached at high axial strains of about 12%. It should be noted that Young's modulus E, identified throughout triaxial testing, is generally a pseudo-elastic modulus determined in a range of deformation where irreversible deformations are observable (10 -3 ). Then, it is not a "true" modulus of elasticity that must be determined in the deformation ranges of 10 -6 -10 -5 .

Table III

Hand-rearranged rockfill sub-layer

The hand-rearranged rockfill layer denoted herein as sub-layer (for the pitching) on both sides of the dam is modeled using the classical Mohr-coulomb (MC) model. The parameters were defined on the basis of literature data (Table III.2):

• Rockfill was manually rearranged in a denser state, which makes it possible to maintain an assumption of isotropic behavior for the material. The arranged nature of the rockfill, • The internal friction angle has been identified to 58°, in accordance with LEP's diagram [START_REF] Thomas | Review of shearing strength of rockfill[END_REF] for a dense material and low confining pressure;

• The dilation angle has been accordingly identified to 18°in correspondence with its dense state using Bolton studies [START_REF] Md Bolton | The strength and dilatancy of sands[END_REF];

• The Poisson's ratio has been taken equal to 0.25 and cohesion is equal to 0 by definition for a coarse granular material;

• Young's modulus E has been taken equal to 11 times the value of the rockfill dumped into the dam body. This choice was made according to the equation below (III.5) developed by François Laigle [2][42]. Since the porosity n of the backfill is 0.45 and that of the handplaced rockfill 0.2 [START_REF] Vincens | Dry Stone Retaining Structures: DEM Modeling[END_REF], the pseudo-elastic Young's modulus of the hand-placed rockfill is approximately 11 times higher.

E = K b .p a . σ 3 p a m (III.5)
where K b = 4/n 3 , p a =100 kPa is the atmospheric pressure and σ 3 is the confining pressure. The parallelepipedic stones were modeled as infinitely rigid bodies with deformable contacts. The contact is linear elastic involving normal and tangential stiffnesses and the transmitted tangential force is limited by Coulomb's law of friction. A friction coefficient of 0.7 (tangent of the friction angle of the block-block contact) was used for the block-block contact of the pitching: this was obtained from a previous study (slip tests) for granite stones with similar dimensions [START_REF] Lr Alejano | Stability of granite drystone masonry retaining walls: Ii. relevant parameters and analytical and numerical studies of real walls[END_REF]. The effective modulus and stiffness ratio used were set at 5.10 7 N/m 2 and 2, respectively, which are the same as the ones used and justified in the scaled-down model in Chapter 2. An overall damping of 0.7, classically used in quasi-static calculations [START_REF] Peter | A discrete numerical model for granular assemblies[END_REF], was considered.

The 20 cm thick upstream concrete mask was considered in our model, by increasing the weight of the upstream blocks and adding cemented bonds between the blocks, considering a tensile strength of 3.9 MPa and a cohesion of 3 MPa: these are average values obtained from Brazilian concrete tests, taken from a previous study [START_REF] Lelovic | Determination of mohr-coulomb parameters for modelling of concrete[END_REF].

Interface model parameters

Pitching-backfill interface

The coupling between PFC3D and FLAC3D at the interface between the pitching and the backfill is activated automatically. A purely frictional linear elastic model, using the same set of parameters as for the block-block contacts (see section 3.3.4), was used to simulate the mechanical behavior of the interface.

Dam-foundation interface

For the modeling of the dam-foundation interface, as for the pitching-backfill interface, a linear elastic model with Coulomb friction was used. An interface friction angle of 42°(equal to the internal friction angle of the backfill) was imposed. It is typical of very rough interfaces involving a stiff material in contact with a granular material. The normal and shear stiffnesses K n =5.10 9 Pa/m and K s =5.10 9 Pa/m are defined in accordance with the recommendations of the FLAC3D manual.

Dam construction phase

The The displacement field is not symmetrical between the upstream and downstream sides, with displacements being greater on the downstream. The dissymmetry is mainly due to the presence of a berm on the downstream side with a zone where shear deformations are concentrated. It is the consequence of the transfer of a part of the weight of the upper pitching towards the backfill with a punching effect of this latter. The dislocation of the berm with its detachment seems to confirm this hypothesis.

In general, there is a gap in the continuity of movement between the pitching and the backfill, indicating sliding between these two subsystems, made possible by the presence of a frictional interface. Friction at this interface contained the displacement field in the backfill close to it.

While the upper pitching shows displacements of the same order of magnitude as the backfill, this is not the case for the lower pitching, which shows small displacements on the overall. As the backfill in this area has settled, it seems that isotropic deformations have essentially been generated in the backfill. The discontinuity in the displacement field between the backfill and the pitching is very marked here, resulting in strong sliding between the backfill and the pitching in this zone.

The berm thus seems to have cut off all the continuity in the transmission of forces within the pitching, and the upper and lower downstream pitchings appear to operate relatively independently. As a result, the berm induces a weakness in the mechanical behavior of the pitching. The distribution of normal forces confirms the independent behavior of the upper and lower downstream pitchings. Indeed, the cumulative effect of the weight of the upper pitching that increases the value of the normal forces ends at the berm. A new cumulative effect of the weight of the lower pitching starts at the berm and accumulates till the dam toe. Joints under the berm generate more tangential forces: this is linked to the de-structuring of the joints in the area where the upper pitching punches the backfill. One can note that the ratio of mobilized friction in the block-block contact at the end of construction is very low with a ratio of about 20% except close to the berm which is not critical.

The forces at the backfill-pitching interface are given in Figure III.8. There are two main points to note. The maximum friction at the interface is mobilized on the lower downstream pitching, confirming sliding between the pitching and the dam body as a result of the significant settlement undergone by the backfill under the berm. In addition, in the upper and lower zones of the berm, there are losses of contact between the backfill and the pitching (transparent zones -not affected by color), which should be compared with the zones where the mobilized friction between the blocks of the pitching is greater (Figure III.7), which is consistent. This finding confirms the role played by the berm as both a point of mechanical discontinuity in the pitching and a weak point constraining the movement of lower pitching blocks in its vicinity.

Reservoir water filling phase

After construction, the water level was gradually increased with 1 m increments until the maximum level was reached (15 m above the upstream foundation). The resulting hydrostatic pressure was simulated by applying normal forces to the external faces of the upstream pitching stones.

The force applied to each stone depends on its position relative to the water level at any given time. The distribution of forces applied to the upstream stones at the end of water filling is shown in Figure III.9. As the upstream face of the dam is impermeable due to the presence of a concrete mask, no water flow through the dam body is simulated and thus there is no water pressure in the dam body. Thus, only the mechanical action has been modeled in this work. Water pressure leads to displacements throughout the dam body, reaching all the way to the downstream part of the dam, including the lower part. One can note that there is a central zone of the dam body that experiences small displacements. This central zone is then supposed to act as a screen for the downstream part of the dam and damp out horizontal movements which is not the case. We can observe that larger vertical displacements are generated close to the upper downstream pitching. It may be the consequence of the rearrangement of the upper pitching which effect propagated through the denser sub-layer towards the dam toe. Then, the rockfill located close to the lower downstream pitching may have been more disturbed than expected.

The horizontal and vertical displacements at the crest are 6 cm (downstream direction) and 10 cm respectively, which is very close to the extrapolated historical values described in Section 2.3.2 in Chapter 1 (around 6 cm horizontal and 11 cm vertical displacement at the crest). This validity, which was expected since the set of LKE model parameters was calibrated by Laigle [START_REF] Laigle | Réévaluation de la stabilité du barrage de gréziolles : étude du comportement mécanique sous sollicitations statique et sismique[END_REF] for this purpose, is confirmed in this work. Concerning the pitching, significant forces were only created at the base of the upstream pitching which was expected (Figure III.12). In the downstream pitching, the filling did not lead to any change in the distribution of normal forces. Shear forces changed mainly at the base of the pitching. This lower part seems to take up horizontal forces induced by hydrostatic thrust and not absorbed by the backfill. This aspect is accentuated by the boundary condition at the base of the downstream pitching, where the first layer of blocks is fixed in the model. It may also be noted that some contacts are lost during filling at the bottom of the downstream pitching (Figure III.13).

Figure III.14 gives the intensity and direction of displacements, and within the dam body, after construction (left) and after water filling (right). It can be seen that the lower downstream part of the dam body moves parallel to the pitching, which justifies the very low displacements observed on the lower part of the downstream pitching and the high frictional mobility at the interface in this region.

Safety factors computation

The two approaches described in Section 2 have been used to assess the dam stability. In those approaches, the detection of failure is done using the developed criteria which is described in paragraph 3.2.

The first safety assessment was carried out with and without a water filling phase. The strength parameters of the dam body and pitching were progressively reduced until failure is reached. The resulting static safety factor was found equal to 1.4 with and without water pressure which is very close to what was found by Laigle in his full-continuum study of the dam ( see section 5.1.2 in Chapter I). Filling with water had little or no influence on the strength and safety of the dam: this The seismic safety factor was obtained by applying a rotation to the gravity vector until dam failure is reached in the direction where the weaker slope is involved. The dam failed at a rotation of 8°, giving a safety factor of about 1.33. This safety factor is consistent with that obtained by It may be noted that large shear deformations are initially initiated at the foot of the downstream pitching, in the most heavily loaded zone, and then propagate towards the crest, almost parallel to the pitching.

Influence of different elements in the dam stability

In this section, the different stages of the complete dam construction are presented to highlight the separated contribution of different elements constituting the dam in its stability.

The simplest model was first studied and tested, then the complementary elements and conditions

Initial model (Backfill+homogeneous foundation)

In a very first step, a very simple model was built: the simplest model considers only the "backfill" part forming the body of the dam, without taking into account the pitching, the hand-rearranged rockfill sub-layer, the concrete mask, and the berm. The dam foundation was considered to be a homogeneous granodiorite rock mass.

This first model was built in four stages, but the dam failed as soon as the 1 st construction layer was run (Figure III.19). This behavior was expected since the backfill has an equivalent friction angle of around 42°, which is less than the slope angle of the two faces of the dam.

In conclusion, this initial model is not stable, leading to a safety factor of less than 1. Finally, the safety factor of the dam was studied quantitatively using the two approaches described above. The static safety factor is found equal to 1.4. Next, pseudo-static tests are carried out by tilting the dam to failure. The dam fails at a rotation angle of 8°, giving a seismic safety factor of 1.33. In conclusion, the role of the pitching in dam stability has been quantitatively justified. The pitching increased the dam safety from less than 1.0 to more than 1.3.

Berm

A dam berm is a horizontal zone that breaks the continuity of the dam slope. It is widespread in rockfill dams with pitching in France, mainly on the downstream side. There are rockfill dams with one, two or even several berms, depending on the height of the dam. The berms of these dams are generally 1 to 2 m wide, according to several examples provided by EDF [START_REF] Vincens | Dry Stone Retaining Structures: DEM Modeling[END_REF].

To better understand its effect, the berm was added on the downstream side and the dam behavior was again studied and compared with the previous case (without berm). Displacements at the end of dam construction with the berm are compared with the initial case without berm (Figure III.22).

The presence of the berm tends to significantly increase crest displacement (factor 2). The upper part of the backfill moves sharply downstream. This displacement is linked to the initiation of a "landslide" at the foot of the upper downstream pitching. The fracture surface opens up beneath the berm, causing it to rotate upstream. This would carry away the first blocks of the upper part of the lower downstream pitching, which would explain the horizontal "crack" observed on site. Shear strains on this sliding surface are of the order of 30% (Figure III.22), and the angle of friction mobilized at this level is that corresponding to the residual state. The factor of safety obtained by reducing the properties decreases from 1.4 to 1.3 for the dam with berm. Pseudo-static rotation tests were also carried out. The dam with the berm fails at a rotation angle of 6°(corresponding to a safety factor F θ = 1.24, i.e. 2°less than the dam without the berm (F θ =1.33).

Hand-rearranged rockfill layer under the pitching

Underneath the pitching, the rockfill has been more densely arranged by hand to ensure that the pitching can be laid in the right conditions. This layer is 90 cm thick on the downstream side and 20 cm thick on the upstream side. The aim of this section is to study the influence of this layer on the behavior of the dam. The displacements and shear strains of the dam at the end of construction are shown in Figure III.25 (with and without the sub-layer). The presence of the hand-rearranged layer reduces the negative effect of the berm, diminishes the punching effect of upper downstream pitching on the dam body and makes the displacement field more symmetrical. It also limits the intensity of shear strains emerging from the downstream face. The safety factor of the dam was studied with the hand-rearranged layer to quantify its influence on its stability. The static safety factor increases from 1.3 to 1.4. This represents an increase in the dam safety of around 7%, justifying the important role of this layer in the dam's stability. The dam also fails at an angle of rotation of 8°, resulting in a seismic safety factor of 1.33 (initially 1.24 without this layer and with the berm).

Concrete mask

The last element to be integrated into the final model is the concrete mask on the upstream face of the dam. To study its effect on the dam behavior, it is now added to the previous model (model with berm, hand-rearranged rockfill sub-layer and non-homogeneous foundation). As mentioned in paragraph 3.1, the 20 cm thick upstream concrete mask has not been modeled explicitly but has been integrated into the upstream mask, by increasing its weight and adding cohesion and tensile strength between its blocks compatible with those of the concrete material. From a safety point of view, the removal of the downstream pitching has considerably reduced the dam safety factor. The static safety factor was found of 1.15, representing an 18% reduction in dam safety. After a pseudo-static rotation test, the dam failed very quickly, after a rotation angle of only 3°: this is equivalent to a safety factor of 1.11. These results conclude that the dam without the downstream pitching is stable, but very close to its failure limit. In this part of the study, a pessimistic case was considered and studied in order to predict the behavior and safety of the dam under more unfavorable conditions.

The conditions and parameters of the hand-rearranged layer were modified. Porosity is now considered equal to 0.3 instead of 0.2, which changes its density from 2000 kg/m 3 to 1750 kg/m 3 . In addition, its elastic parameters were modified accordingly: the modulus of elasticity is now taken equal to 4 x 75 MPa instead of 11 x 75 MPa (see section 3.3.2). The angle of friction was modified from 58°to 55°, which is an average value for the considered porosity and confining pressure. Displacements after water filling are also compared (Figure III.36): horizontal and vertical displacements after water filling are very similar overall. The only difference is slightly higher displacements in the upper downstream region, due to lower stiffness.

The static safety factor decreased from 1.4 to 1.35 for the modified pessimistic case. In addition, pseudo-static rotation tests were carried out. The dam, in the modified pessimistic case, fails at a rotation angle of 7°(corresponding to a seismic safety coefficient F θ = 1.28, i.e. 1°instead of 1.33 in the initial studied case). In conclusion, even with the pessimistic parameters for the handplaced layer, no major change is observed in the displacements and the safety factors for the dam stability which reassures the used model parameters for the sub-layer.

Summary of the previously studied cases

The results of all the cases studied in Section 4 are summarized and presented in Table III.4. They include the vertical and horizontal displacements at crest and mid-height (berm region), and the safety factors for each case. Note that displacements are given in cm.

Modeling of high deformations in the downstream pitching

Several large deformations have been observed on the studied dam, during its operational period since the 1950s. The main observed features include bulges or buckling developed in the lower downstream pitching and a longitudinal crack under the downstream berm (Chapter I, section 2.3). The aim of this part of the study is to identify and understand the causes of these deformations by trying to reproduce them numerically and to determine whether they have an adverse effect on the dam's stability (reduction in its safety factor).

Geometry imperfection in the pitching blocks

First, a geometric imperfection was introduced in the shape of the blocks, which is considered not to be perfectly parallelepipedic. A slightly trapezoidal shape was then used with a 3 mm opening on the bottom (internal) side of the block in contact with the backfill. The shape defect is therefore quite small (Figure III.37). The imperfection reduces the block's contact area with the backfill (19.4 cm (20 -0.3 x 2) x 40 cm instead of 20 cm x 40 cm). It also reduces the surface area in contact between the blocks (6 cm (= 0.2 x 30) x 40 cm instead of 30 cm x 40 cm), giving the pitching greater flexibility. This shape was used to model the PEDRA test campaign (Chapter II), and this configuration for the geometrical imperfection is the one that led to maximum internal mobility for the pitching. At the end of the construction phase, no significant differences are observed in displacements and contacts. Differences appear after the water filling phase: more concentrated displacements are observed in the downstream pitching just below the berm, with the use of trapezoidal blocks which could explain the longitudinal crack found on-site at this location. In addition, buckling with loss of contact with the dam body was observed at the base of the upper downstream pitching just above the berm. This is due to the greater rotational flexibility of the trapezoidal blocks, in addition to the compressive forces developed in this part of the pitching. In the lower part of the downstream pitching, after water filling, a greater loss of contact is also noted with the use of trapezoidal blocks. In addition, the use of trapezoidal blocks leads to greater shear strains on the lower downstream face of the backfill.

The safety of the dam was again estimated with the use of trapezoidal blocks. The factor of safety, which is estimated by parameter reduction, did not change with the use of trapezoidal blocks (F φ = 1.4). However, the dam failed for a rotation angle of 7°(rotational safety factor F θ equal to approximately 1.3) which is slightly lower than the case with parallelepipedic blocks (failure angle of 8°: F θ = 1.33).

Influence of the valve chamber

The main deformations observed on the Escoubous dam are the bulges or bellies on the lower part of the downstream pitching. More specifically, these bulges are concentrated in the area around the valve chamber, which is assumed to be a rigid structure (Figure III.41a). This structure could have a probable influence on the development of these bulges because of the restrictive displacement conditions it induces on the pitching.

To simulate the presence of this structure, some of the stones at the bottom of the downstream pitching were fixed (Figure III.41b). The dimensions of the structure (1.6 m × 2.25 m) were estimated visually from photos supplied by EDF. The width of the model was doubled from 9 m to 18 m in order to provide a sufficient number of stones in the horizontal direction. This allows the boundary conditions to be moved away from this rigid structure. In this model, blocks with geometrical imperfections are used (as the ones used in the previous section).

To evaluate numerically the effect of the rigid structure on the dam behavior and bulges formation, the normal displacements at the dam downstream pitching face and interface contacts (represented by disks) are compared for the case with and without the rigid structure (Figure III.42).

The concentrated displacements were found just above and on both sides of the rigid structure. These displacements are higher than in the immediately surrounding area. In addition, contacts are lost in the zone where the largest normal displacements are concentrated. This leads to the conclusion that bellies develop around the rigid structure with a loss of contact with the body of the dam. However, the bellies obtained numerically after water filling are less severe than the actual bellies on site. These numerical bellies could evolve further when modeling the creep behavior of the dam body together with successive water filling/emptying phases, which are outside the scope of this study. A 3D effect linked to the shape of the valley could also contribute to an increase in deformations.

Another concentrated area of normal displacement had already been noted just below and beneath the berm. The rigid upstream rotation of the berm had been considered as the origin of these longitudinal "cracks" found on site.

The normal and shear forces and the ratio of mobilized friction for each block of the pitching are also analyzed, at the end of construction and after filling with water (Figure III.43). The previous conclusion is reinforced when friction is strongly mobilized in the area around the rigid structure, in addition to the part under the berm. This could be justified by the high resistance shear forces observed in these areas. The static safety factor was also evaluated and found equal to 1.45 without water pressure, which is slightly higher than the 1.4 found for the case without a rigid structure and blocks with a geometrical imperfection. This could be justified by a mobility restriction zone played by the presence of the rigid structure. With the addition of water pressure, which leads to the formation of bellies around the rigid structure, the reduction safety factor is now 1.4, which is slightly lower than the case without water pressure.

Secondly, a pseudo-static rotation test was carried out: the dam failed at a rotation angle of 7°, giving a seismic safety factor identical to the case without a rigid structure (using blocks with geometrical imperfections). Two sections of the pitching are given after construction and at the last step before failure (tilting angle = 6°): the first section is at exactly the middle of the model and the second one is just 1 m far from the dam lateral boundary (Figure III.44). In section 1, we can see the formed belly above the rigid structure whereas we can see the belly beside the rigid structure in section 2 which is lower in position than the first one. A supplementary case with pitching blocks having geometrical imperfections on the bottom internal side but with 6 mm opening instead of 3 mm was tested also. Identical results and conclusions

Conclusion

A rockfill dam with dry-stone pitching, located in the Pyrenees mountains in France and built in the early 1950s, was studied using a mixed FDM-DEM numerical approach. This approach was already validated in a previous study on similar structures.

The dam was first constructed in several stages and then the water reservoir was filled. The set of model parameters was validated by comparing the simulated horizontal and vertical displacements after water filling with the expected ones taken from reconstructed historical data and from EDF experience on the displacement trends for such structures. We have seen that the upper and lower part of the downstream pitching (which are divided by the berm) behave independently. The transmission of forces toward the foundation for the upper downstream pitching is then cut and transferred to the denser rockfill sub-layer in the dam body. This tends to create a punching of the dam body material below the berm. Sliding between the lower downstream pitching and the backfill is found due to the isotropic deformations of the downstream lower part of the backfill below the berm which is also coherent with the mobilization of the maximum friction at the backfillpitching interface of this part. After water filling, a shear surface emerging on the downstream face appeared with shear strains up to 7%. The shear forces at the bottom of the downstream pitching increased after water filling which is the result of the resistance of the forces induced by the water pressure and not absorbed by the dam body in addition to the fixity of the first row of blocks.

Two approaches were used to evaluate the safety of the dam, one linked with the static safety margin and the second one with the equivalent seismic safety margin. Safety factors in the range of 1.3-1.4 were obtained in both approaches. In addition, in the last stages before failure, the pitching resistance role (retaining of the emerging surfaces) was illustrated by the strong mobilization of the friction angle between its blocks.

The effect of each of the dam elements on its mechanical behavior and stability was also studied in detail. From this study, it has been shown that the dam is not stable without the stone pitching and the denser rockfill sub-layer. The addition of the pitching stabilized the dam. By adding the berm on the downstream side, it has been justified that it represents a weak point in the dam and decreases its safety. This berm decreased the stability role of the downstream pitching by dividing it into two independent parts with no force transmission. The addition of the denser rockfill sublayer increased the dam safety and decreased the negative effect of the berm. Afterward, the usage of a non-homogeneous foundation increased slightly the displacements due to the lower rigidity of the foundation that increases the ability of the dam's overall movement. Finally, the integration of the concrete mask increased the displacements after construction due to the added weight but decreased them after water filling due to the stiffness increase of the upstream pitching.

Supplementary studied cases have shown the complete dam elements but without the stone pitching is not stable which was expected since it is characterized by a friction angle lower than the slope angle. The dam failed at the upstream part even before the completion of construction. This is because the sub-layer is smaller in width on the upstream side. Therefore, another test was done by considering the complete dam element but without only the downstream pitching (the upstream pitching is modeled). The dam was stable after construction but with a static safety factor close to 1 which means that the dam is very close to failure without its downstream pitching. Finally, a pessimistic case was considered with a lower Young's modulus and a lower friction angle for the hand-placed sub-layer. This case gave very similar displacement fields and safety factors and no major changes have been observed which comfort the used model parameters for this specific layer of the dam body.

Finally, the bulges found on the bottom part of the downstream pitching were investigated. As a first step, the pitching was modeled using blocks with geometrical imperfections. Such blocks added more flexibility to the pitching which leads to higher local concentrated deformations in the pitching and more ability to show bulges. The second step represented the consideration of the valve chamber found in the middle of the base of the downstream pitching which restricts its displacement around it. The normal displacement of the downstream pitching after water filling has shown local displacement concentration in the area around this structure but with lower intensity than the actual ones observed on site. The safety factor of the dam was found close to the initial case without the rigid structure.

Chapter IV

True dynamic behavior of Escoubous dam

Introduction

An earthquake is a natural hazard that occurs all over the world and is due to the movements of the earth's crust (tectonic plates). It releases huge energy that could damage human-made structures such as buildings, bridges, dams, etc. Its power can vary from time to time and from one place to another depending on the presence of nearby tectonic faults, their intensity and the geological composition of the subsoil. An earthquake's magnitude is normally measured by several definitions that depend mainly on the released energy and sometimes on the distance from the epicenter. The earthquake signal is also characterized, at a certain ground point, by a PGA (peak ground acceleration), PGV (peak ground velocity) and PGD (Peak ground displacement) which are the maximum values of the acceleration, velocity and displacement in time, respectively, throughout the earthquake duration. The intensity of the signal (PGA, PGV and PGD) depends on the neighboring environment of the measured point and the bottom geological soil layers which changes the intensity of the upcoming signals. The signal is also characterized by a spectrum that allows us to determine the dominant frequencies (or periods) of the signal, that is to say, the range of frequencies where most of the seismic energy is concentrated.

To simplify the study of the behavior of the civil engineering structures subjected to earthquake loading, the measured signals data are used to form a seismic risk zone distribution for each country. With this aim of simplification, the regulatory frameworks concerning seismic risk have defined a reference acceleration a gR for the different seismic zones, such that the probability of occurrence of a higher level earthquake occurring at that location is less than a certain percentage.

In France, as a reference case, this probability of exceedance percentage is conventionally taken to be equal to 10 % in fifty years (earthquake return period of 475 years) [START_REF] Pecker | Overview of Seismic Regulations for French Industrial Facilities[END_REF][START_REF] Loudiere | Risque sismique et sécurité des ouvrages hydrauliques[END_REF]. Figure IV.1 and Table IV.1 give the reference accelerations a gR to be considered for each French seismic zone. As a result, since 2011 and according to article D.563-8-1 of the French Environment Code [START_REF]Article d.563-8-1 of the french environment code[END_REF], the seismic risk must be taken into account a priori for any structure design. Exceptions related to certain structures, seismic zones and the coefficient of importance of the structures are explicitly given when this calculation is not expected [START_REF] Légifrance | Arrêté du 22 octobre 2010 relatif à la classification et aux règles de construction parasismique applicables aux bâtiments de la classe dite « à risque normal[END_REF].

In this chapter, the time-dependent seismic behavior of Escoubous dam is studied. The dam is located in zone 4 seismic area with an a gR = 1.6 m/s 2 ≃ 0.16 g. The main objective of this chapter The Antilles 3 is to retrieve the dam's response to a probable earthquake in its region and to quantify its safety and stability in case of such an event kicks the structures. Other objectives include the study of the effect of the horizontal incident motion acting solely and the cumulative effect of horizontal and vertical incident motions on the structure. The effect of the soil layers under the dam and of the reservoir water pressure on its behavior and deformation under a seismic motion are also addressed.

The numerical model with its adapted dynamic boundary conditions is first presented. Then, the constitutive models used in this dynamic study are shown and their parameters are estimated and calibrated including the elastic, strength and damping parameters. The damping model used in the dam body and foundation was calibrated by fitting the shear modulus degradation and damping ratio curves to suggested curves from the literature. Then, the artificial signals used in this study are presented. They were generated based on a target acceleration spectrum (damping of 5%). The results obtained in several cases are then presented and discussed.

Numerical model and its dynamic boundary conditions

The numerical model of Escoubous dam described in Chapter III is used for the dynamic modeling but with adapted dynamic boundary conditions. The mesh of the numerical model was verified to ensure the accurate transmission of seismic waves throughout the model. Indeed, Kuhlmeyer and Lysmer (1973) [START_REF] Roger | Finite element method accuracy for wave propagation problems[END_REF] showed that for an accurate representation of the wave transmission through the soil model, the spacial element size ∆L, must be smaller than approximately one-tenth to oneeighth of the wavelength λ associated with the highest frequency component of the input wave (Equation IV.1).

∆L ≤ λ min 10 (IV.1)
Where λ min is the wavelength associated with the highest frequency component f max = 30 Hz (cutoff frequency) that contains appreciable energy: λ min = V wave f max where V wave is the wave propagation speed.

This analysis leads to Equation IV.2. This condition was verified for all the zones in the numerical model.

∆L ≤

V wave 10 × 30 (IV.2) Furthermore, the following boundary conditions (Figure IV.2) were applied to the numerical model during the dynamic loading:

Quiet boundary

A quiet boundary (also called absorbing or compliant boundary) was applied to the foundation base to represent an infinite extension medium for the foundation. It mainly prevents the reflection of the outward propagating waves back into the model by absorbing them. In FLAC3D, the viscous boundary created by Lysmer and Kuhlemeyerin 1969 is used [START_REF] Lysmer | Finite dynamic model for infinite media[END_REF]. Quiet boundaries are composed of independent dashpots in the normal and shear directions. When body waves approach the boundary at angles of incidence greater than 30°, the method is almost completely effective in absorbing them. There is still energy absorption for lower angles of incidence or for surface waves, but it is imperfect. The advantage of this scheme is its time frame action. In both finiteelement and finite-difference models, its effectiveness has been proven [START_REF] Rr Kunar | A tested soil-structure model for surface structures[END_REF]. The provided viscous normal and shear stresses are given by equations IV.3 and IV.4. In FLAC3D, the stresses t n and t s are calculated and applied at every timestep in the same way boundary loads are applied. Where ν n and ν s are the normal and shear components of the velocity at the boundary, ρ is the mass density and V p and V s are the p-and s-wave velocities.

t n = -ρV p ν n (IV.3)

Free-field boundary

A free-field boundary was applied to the vertical sides of the system (sides normal to the x-axis) on which the dam lays. Quiet boundaries cannot be used on those sides when the dynamic signal is applied at the model base since the energy will leak out of those sides. On site, there are no true lateral boundaries but a semi-infinite medium. In a numerical system with definite vertical boundaries, waves generated in the system radiating outwards will be reflected on the boundary towards the inner system. The free field boundary solves this issue by retaining the non-reflecting properties of the sides (i.e., outward waves originating from the structure are properly absorbed). The continuum finite-volume code NESSI (Cundall et al. 1980 [START_REF] Pa Cundall | Nessi, soil structure interaction program for dynamic and static problems[END_REF]) used this strategy. For FLAC3D, a method of this kind that involves running free-field calculations in parallel with the main-grid analysis was developed. Two plane free-field grids make up the free-field model side boundaries and corners. To ensure that there is a one-to-one correspondence between gridpoints in the free field and in the main grid, the plane grids are generated to match the main-grid zones on the side boundaries. The calculations used for the plane free-field grids are two-dimensional and assume infinite extension in the plane's normal direction. The plane free-field grids consist of standard FLAC3D elements, which have gridpoints constrained in such a way as to achieve the infinite extension assumption. The lateral boundaries of the main grid are coupled to the freefield grid by viscous dashpots to simulate a quiet boundary and the unbalanced forces from the free-field grid are applied to the main-grid boundary.

Dynamic signal

The dynamic signal is applied at the base of the foundation. The goal here is to obtain the target outcrop motion on the downstream surface of the foundation (Figure IV.2). Since the downstream part of the foundation is characterized by a linear elastic behavior (granodiorite), the deconvolution analysis is not necessary. The outcrop signal can directly be used to apply the dynamic load at the foundation base.

Since a quiet boundary is applied at the foundation base, the acceleration time history of the dynamic load cannot be directly applied as mentioned in FLAC3D manual. Instead, the signal at the base is applied by converting the velocity history of the dynamic load into normal and shear stresses and then applying those stresses to the foundation base. The stresses are applied according to equations IV.5 and IV.6. The "0.5" factor is found in those equations because the surface signal (outcrop) is double the information held by the incident signal. Indeed, the measured surface signal is the sum of the incident and reflected waves [START_REF] Lay | Chapter 4-surface waves and free oscillations[END_REF] at the boundary. In addition, the "2" factor must also be added to the equations since half of the applied signal at the foundation base is absorbed by the quiet absorbing condition at the base. The applied vertical and horizontal components of the signal are described in detail in section 4.

σ n = 2(ρ.C p ).ν n × 0.5 (IV.5) σ s = 2(ρ.C s ).ν s × 0.5 (IV.6)
where:

• σ n and σ s are the applied normal and shear stress respectively;

• ρ is the mass density;

• C p and C s are the speed of p-wave s-wave propagation in the medium, respectively;

• ν s is the input shear velocity: ν s = V (1749 signal);

• ν n is the input normal velocity:

ν n = 2 3 x V (1751 signal).
Other conditions:

• The two lateral sides of the modeled dam + foundation (vertical sides with a normal colinear with y-axis in Figure IV.2) are fixed in their normal y-direction and free to move in the other directions (x and z directions)

• The two columns of the pitching at the lateral limits are fixed in the y direction and free to move and rotate in other directions. The velocity of the first row of blocks in contact with the ground surface is updated each 100 cycles to be equal to the velocity of the foundation gridpoints in contact with it. This is done to ensure the signal transmission from the foundation to the pitching that rests on it.

Constitutive model

In this chapter, due to the lack of robustness of LKE model in dynamic computations, we used an incremental calibrated hysteretic model (see section 3.2) coupled to a Mohr-Coulomb (MC) model that limited the range of stresses to those with a physical meaning. Plastic deformations will be created accordingly when MC failure surface will be reached. This is not detrimental to the study carried out, as MC model provided an accurate estimate of the dam's resistance and the incremental hysteric model is known to hold the advantage of retrieving correctly non linear phenomena observed throughout cyclic loadings for a range of cyclic deformations smaller than 10 -3 . Moreover, one must remind that the principal aim of this work is to estimate both the static and seismic resistance and safety factor of the dam which are properties driven by stresses and less by deformations (in the studied case).

The static state of the model, after construction or after water filling but before any dynamic computations, is firstly simulated using LKE model to get an initial stress and deformation state closer to reality. The final static deformations and stress states of the dam are recovered and then used as the model's initial state of the dynamic simulations using MC model coupled with the hysteretic damping model. The strength and elastic parameters of MC model are given and discussed in section 3.1

MC model parameters

The definition of the elastic parameters (mainly the shear modulus) of the studied soils is a very important step in the dynamic study to correctly retrieve their behavior under seismic loading. Indeed, each soil is characterized by an initial maximum elastic shear modulus G max that degrades with the increase of the cyclic shear strain.

The maximum shear modulus of each of the moraine and alluvium layers in the foundation has been estimated using the equation proposed by Seed [START_REF] Bolton Seed | Moduli and damping factors for dynamic analyses of cohesionless soils[END_REF] and described in Chapter I. K 2 was taken equal to 90 as per the study of Seed [START_REF] Bolton Seed | Moduli and damping factors for dynamic analyses of cohesionless soils[END_REF] and as in the dynamic study done by Rezaee [START_REF] Rezaee | Seismic analysis of rockfill dam, a case study of ghoocham dam[END_REF]. For the moraine layer, a maximum shear modulus of about 300 MPa was obtained for an average pressure of 200 kPa (average value in the moraine layer), and for the alluvium layer, a shear modulus of about 200 MPa was obtained for an average pressure of about 105 kPa (average value in the alluvium layer). The shear modulus for the granodiorite layer and for the concrete barrier were taken as in the static model in Chapter 3 due to their high rigidity.

Concerning the dam body, the maximum shear modulus was estimated using the equation proposed by Kokusho [39] for the rockfill material and presented in Chapter 1. A void ratio of about 0.818 (porosity n = 0.45) corresponding to the rockfill loose state has been used in the equation.

The equation becomes G max = 13.07.(σ ′ m ) 0.55 MPa for this void ratio. The maximum shear modulus is then calculated for each zone depending on its final confining pressure (after construction or after water filling, the value was retrieved numerically) before the beginning of the dynamic loading. The obtained G max and V s profiles along the vertical column under the dam crest is plotted and shown in Figure IV.3. It is also compared with the approach of Gazetas [START_REF] Gazetas | A new dynamic model for earth dams evaluated through case histories[END_REF][START_REF] Dakoulas | A class of inhomogeneous shear models for seismic response of dams and embankments[END_REF] who estimated the G max profile using equation IV. 7 (Figure IV.4). G max,b corresponds to the maximum shear modulus at the dam base which was estimated using Kokusho's equation [START_REF] Kokusho | Cyclic triaxial test on sands and coarse materials[END_REF] at a confining ). The small difference may be due to the asymmetrical shape of Escoubous dam which includes a berm in contrast on the downstream slope with the symmetrical shape assumed by Gazetas.

Kokusho's equation was also used to estimate the maximum shear modulus of the hand-rearranged rockfill sub-layer on the 2 faces of the dam. A void ratio of about 0.25 (porosity n=0.2) corresponding to the dense state of this layer was used. For a mean pressure estimated to be about 50 kPa (mean pressure in the sub-layer at the final state before the beginning of the dynamic loading), a maximum shear modulus of about 330 MPa was obtained and used.

The strength parameters including the friction angle and cohesion and the dilatancy angle are all the same as the ones used in Chapter II for the different foundation elements and the hand-placed sub-layer. For the dam body, the internal friction angle was estimated to be equal to 42°using

Hysteretic damping model

A crucial step in the dynamic modeling of soils is the definition of damping for the involved materials. It is important to capture an accurate realistic soil response and its energy dissipation during dynamic loading. Several numerical damping formulations are found including the Rayleigh damping and the hysteretic damping which are the most commonly used. Rayleigh damping is a viscous damping that is proportional to a linear combination of mass and stiffness mobilized in the shaking. It can be a simple and effective way to model the soil damping response but in an explicit FDM scheme, it decreases the required time-step drastically which may increase excessively the computational time. In addition, it does not capture the frequency-independence of damping in soils in a correct way.

On the other hand, the hysteretic damping has no effect on the time step and the computational time in FLAC3D. In fact, damping is provided throughout cyclic movements by a time increment constitutive model where the non-linear elastic response of the soil can be accounted for. Indeed, with such a model, the secant shear modulus can be degraded according to an increase in cyclic deformations. Accordingly, the damping which is independent of the frequency is naturally created by the creation of hysteresis loops.

A hysteretic model available in FLAC3D was chosen to capture the damping behavior in both the rockfill and the foundation materials (Equation IV.8). The sigmoidal-4 model is defined through the value of its normalized secant modulus M s (G/G max ) which varies as a function of the cyclic deviatoric γ:

M s = y 0 + a 1 + e -L-x 0 b (IV.8)
where a, b, x 0 , and y 0 are model parameters and L = log 10 (γγ 0 ). γ 0 is the initial deformation value of the moment when there is a reversal of movement direction (peak of response). The The same hysteretic damping formulation was used to apply damping in the moraine and alluvium elements in the foundation but with different model parameters. There is not a detailed knowledge of the component of those two foundation elements and it was assumed by EDF that they are composed of a mix of sand and gravel. Degradation and damping curves of sand given by Idriss and Seed [START_REF] Idriss | Seismic response of soil deposits[END_REF] and of gravel by Rollins [START_REF] Kyle M Rollins | Shear modulus and damping relationships for gravels[END_REF] are used to calibrate the parameters of the hysteretic model (Figure IV.6). The numerical hysteretic curves were successfully fitted by a trial-and-error method to the theoretical average curves of the sand provided by Idriss and Seed. Actually, this curve is also very close to degradation curves derived for gravel by Rollins [START_REF] Kyle M Rollins | Shear modulus and damping relationships for gravels[END_REF]. No damping was added to the granodiorite layer in the foundation due to its high rigidity.

Finally, in addition to the hysteretic damping, a very small percentage of Rayleigh damping (0.05 %) at a central frequency of 3 Hz was added to the model to remove high-frequency noise and to compensate the very low damping provided by the hysteretic model for cyclic deviatoric deformations smaller than 10 -4 . The central frequency was chosen by defining the range of the predominant frequencies in the dam body during the applied dynamic shaking as suggested by the FLAC3D manual. More details are shown in a preliminary elastic study in section 5.1.

A damping of 5% is used for the pitching rigid blocks. Savalle [START_REF] Savalle | Dynamic numerical simulations of dry-stone retaining walls: identification of the seismic behaviour factor[END_REF] proposed such a value in his dynamic numerical simulation of a dry-stone retaining wall. He tested two values of damping: 5% and 15%. He concluded that a value of 5% is sufficient and gives identical results as the 15% case.

Input dynamic signals

The adopted target earthquake is provided by EDF from an internal report. It is characterized by a magnitude of 6.1 at a depth of 18 km and an epicenter that is 4 km away from the site. The ground conditions used to define this spectrum correspond to rock with an average shear wave speed greater than 800 m/s in the first 30 m. For Escoubous dam case, the spectrum is not plotted in the internal report of EDF: only the Peak Ground Acceleration (PGA) of 0.16 g is given, i.e. slightly more than that given by the attenuation law used in the RFS (Règles fondamentales de sûreté) 2001-01 recommendations using Berge-Thierry law [START_REF] Berge-Thierry | New empirical response spectral attenuation laws for moderate european earthquakes[END_REF]. It was decided by EDF to increase the spectrum according to Berge-Thierry law [START_REF] Berge-Thierry | New empirical response spectral attenuation laws for moderate european earthquakes[END_REF] in order to have a maximum acceleration of 0.16 g. The final target spectrum that was used to generate the accelerograms is shown in There are an infinite number of accelerograms corresponding to a given spectrum. A methodology has therefore been adopted by EDF to generate movements consistent with the target spectrum, but also with other parameters known as "harmfulness parameters" whose influence on the modeling can be significant, particularly for multimodal structures and non-linear behavior. This whole methodology is outlined and was validated by Humbert in previous studies [START_REF] Humbert | Guide de génération des signaux sismiques -utilisation de powerspec[END_REF][START_REF] Humbert | Time histories for seismic analysis of structures-pros and cons of available methods[END_REF]. The used signals were cut at 30 Hz to avoid high-frequency components which allows us to limit the refinement of the meshing. In fact, the vast majority of the signal energy lies in the range of frequencies lower than 30 Hz. As a consequence, the use of this cut-off will also not require calibrating again the PGA to 0.16g (Figure IV.10.

In our model, the velocity history of the 1749 signal is used for the horizontal component of the seismic load whereas two-thirds [START_REF] Loudière | Risque sismique et sécurité des ouvrages hydrauliques[END_REF] of the velocity history of the 1751 signal (which is initially a horizontal motion) is used for the vertical component of the seismic motion. Seismic building codes, such as Uniform Building Code or Eurocode 8, generally assume that the peak vertical acceleration is simply a fraction of the peak horizontal acceleration. A value of two-thirds is most often used as the maximum effective ratio between vertical and horizontal acceleration [START_REF] Nathan | Earthquake spectra and design[END_REF]. The usage of other signals than the 5 generated ones for the application of the motion's horizontal and vertical components will be also tested in another study for comparison.

Studied cases 5.1. Purely elastic model without damping

As a first step, after adapting the model to the dynamic conditions, a dynamic simulation was carried out using a purely elastic model for all the materials and without any damping and by applying only the dynamic load in the horizontal direction at the foundation base (signal 1749). This simulation is carried out to provide an estimation of the maximum cyclic shear strain level and the range of predominant frequencies in the dam body. It will also allow us to estimate the representative shear modulus degradation and damping that would be found at such strain values. The maximum cyclic shear strain throughout the dam body was found equal to about 0.05% (Figure IV.11). The corresponding effective cyclic shear strain which is equal to 75% of the maximum cyclic shear strain was found equal to about 0.037%. At this strain level, the shear modulus would be approximately equal to 0.4 of the value obtained at small strains and the damping ratio would reach 12% (Figure IV.12). The range of predominant frequencies in the dam is found between 2 and 6 Hz. The predominant frequencies are neither the input frequencies nor the natural modes of the system, but actually a combination of both. The numerical outcrop motion (horizontal acceleration downstream at the level ground) was measured and its spectrum was plotted (Figure IV.13). It is compared with the target outcrop spectrum and found to be identical which validates our input procedure at the foundation base.

Mohr-Coulomb model with damping

The dynamic behavior of Escoubous dam and of the site was studied gradually to better isolate the contribution of all the parts of the systems to the overall response of the site and the structure. 4 different cases were studied and their details and differences are shown in table IV.3. In both cases, one can note higher displacements at the downstream part of the dam due to the unfavorable effect of the berm. In addition, by comparing cases A and B, it is clear that the consideration of a heterogeneous foundation in case B has led to higher displacements in the dam body. To illustrate this point, the displacements at the crest were also monitored throughout the dynamic loading (Blue curves for Case A and Green curves for case B in Figure IV.15). For case A, the dam crest settled by about 11 cm and moved about 6 cm to the downstream side (observed at the end of the dynamic loading). These displacements are mainly due to the shear sliding of upstream and downstream slopes. These displacements increased in Case B: the crest settles by 3 cm more (14 cm in total) than in Case A (homogeneous foundation) while the crest horizontal displacement is almost the same.

Figure IV

. [START_REF] Marachi | Strength and deformation characteristics of rockfill materials[END_REF] shows the shear strain filed in the dam body for cases A and B at the end of the dynamic load. Two shear surfaces can be noticed on both sides of the dam with values up to about 8% with small local deformation of about 10% for Case A. In Case B, higher shear strains are observed with values up to 12%. ). Indeed, this was expected since the dam is characterized by a natural frequency of about 6 Hz which lies in the dominant frequencies range of the input signal. The natural frequency of the dam was estimated to be around 6 Hz using two approaches proposed by Dakoulas [START_REF] Dakoulas | A class of inhomogeneous shear models for seismic response of dams and embankments[END_REF] and Mahdizadeh [START_REF] Mahdizadeh | Calculation of natural frequency of earth dams by means of analytical solution[END_REF].

The higher displacements obtained in Case B (heterogeneous foundation case) can be justified by the higher accelerations and motions entering the dam body at its base. The horizontal acceleration at the base in the case of a non-homogeneous foundation (Case B) is shown in Figure IV.19(b) and compared with the outcrop motion and the crest horizontal acceleration. We can see that the acceleration at the dam base is higher than the outcrop motion and the "dam base acceleration in the case of a homogeneous foundation (Case A)". This is due to the signal amplifying when it passes by the moraine layer under the dam. This leads to higher accelerations in the dam body and at the crest which reaches values higher than 4 m/s 2 (Higher than the crest accelerations in the homogeneous foundation case (Case A)). These overall higher accelerations eventually lead to higher displacements and deformations in the dam body and the pitching that were shown earlier.

Different components for the shaking

In this part, previous case B results will be compared to case C where the actual foundation (with its complexity) is submitted to a shaking having both a horizontal and vertical component (respectively signal 1749 and two-thirds of signal 1751). The objective is to analyze the additional contribution of vertical shaking to the seismic response of the dam (with respect to case B taken as a reference here). The shear strain of the dam for case C is almost the same as for case B (Figure IV.22) with a slightly higher shear deformation in the region under the berm. This is due to the higher punching effect of the upper downstream pitching caused by its slightly higher displacements toward the berm.

Finally, the horizontal displacements of the monitored points on the downstream and upstream pitchings have slightly increased for case C compared to case B (Figure IV.23). The residual vertical displacements are higher for case C for the upper downstream pitching (red and green In this part, it has been seen that the vertical component of the shaking increases in just a small amount (but not negligible) the computed displacements of the dam body and the pitching after the shaking. Thus, it was preferred to take it into account in the rest of the studies to have a more realistic earthquake scenario. One can note that the displacements decreased remarkably on the upstream side due to the confining effect of the water pressure. On the downstream side, the displacements slightly increased in the region around and above the berm whereas it decreased at the crest region.

It could be due to the lower displacements at the upstream area of the dam which limits the settlement of the crest region. This is also illustrated in Figure IV.25 which shows the displacements at the crest throughout the dynamic load with and without the water pressure. With the presence of water, the crest settles about 5 cm less than the case without water and displaces horizontally toward the downstream side more than double the horizontal displacement of the case without water. In fact, the water pressure acts and pushes the dam towards the downstream side which could mobilize higher horizontal displacements in that direction during the seismic motion.

The shear deformations have decreased to values close to zero on the upstream region of the dam body with the addition of the hydrostatic pressure (Figure IV.26(b)). Indeed, hydrostatic pressure counterbalances horizontal inertial effects in the upward direction. On the other hand, higher shear strains have been found in the case of water pressure in the downstream region. A sliding surface has developed from the middle of the dam to the bottom of its downstream side with values up to 12%. No hydrodynamics effects due to the water reservoir were taken into account in this calculation which constitutes a limit of this work.

Concerning the monitored points on the upstream and downstream pitching ( From now on, the conditions corresponding to case D which is associated with a more complex and complete situation are considered in the following studies and results: the foundation is heterogeneous, the shaking holds both a horizontal and vertical component and a hydrostatic pressure is applied on the upstream side of the dam.

The horizontal acceleration was traced at four different points in the downstream pitching through- out the shaking (Figure IV.28). High accelerations up to 7g were found in the lower downstream pitching (points 3 and 4) which is very high whereas much lower accelerations (less than 1g) were found in the upper downstream pitching (points 1 and 2). Again, as noticed in chapter III, we see that the upper and lower pitching function as two independent systems. The only difference between them is that the lower pitching has a fixed base point whereas the upper pitching has a base point above the berm that can move. The lower pitching would function as a beam that is virtually independent of the embankment over a large part (with no contact with it), embedded at the base and with imperfect embedment (translation and rotation possible) at its top point below the berm. , huge amplifications are found in the pitching. This lower part of the pitching acts as a beam that receives the signal at its base which is embedded in the foundation and then the signal highly amplifies as it reaches its top point below the berm. Indeed, the top of the lower downstream pitching is not well constrained by the berm above it which increases its freedom of movement. The acceleration in the adjacent sub-layer is much lower than the pitching but a bit higher than the adjacent dam body rockfill. This could be due to reflected surface waves or to the added inward signals generated by the high accelerations of the pitching next to it.

Concerning rblock 1 case representing the upper downstream pitching (Figure IV.29(a)), the pitching accelerations are much lower than the lower downstream pitching. Indeed, the wave transmission in the downstream pitching from the foundation is cut at the berm which leads to a low signal at the base of the upper downstream pitching. The acceleration of the upper pitching is mainly due to the signals transmitted from the dam body to it through the interface. This explains the close accelerations between the pitching and the dam body. We can notice also that no amplification is found in the sub-layer under the upper downstream pitching (same acceleration as the adjacent dam body rockfill). This was expected since no external effect of the pitching on the sub-layer is found in this case.

The normal and shear forces and the mobilized friction between the blocks and at the pitchingbackfill interface were compared before and after the shaking (Figures IV.30 and IV.31). The main post-shaking situation is an increase of shear forces around the berm resulting in more blocks where friction is totally mobilized at their contacts. It can be noted also that some of the interface contacts are lost and others are recovered. This can be explained by a destructuring of the berm that limited the influence of it on the arrangement of the pitching blocks.

Dynamic deformations analysis

The obtained results from the numerical dynamic model studied in this chapter need to be compared with the available historical displacements in order to justify its validity. Swaisgood [START_REF] Swaisgood | Embankment dam deformations caused by earthquakes[END_REF] did an extensive review of about 70 available case histories of embankment dams' behavior and settlements during an earthquake. Several embankment dam types were reviewed including concrete face rockfill dams (CFRD), Earth core rockfill dams (ECRD), Hydraulic fill (HF), and the Earth fill dams. The percentage of the crest settlements of the dams was traced as a function of the shaking PGA they faced. This normalization includes not only the dam height but also the total height of the deformable layers on which lays the dam.

The collected data is shown in Figure IV.32. Overall, there is a general trend for an increase in the crest settlement percentage as a function of the increase in the PGA, which was expected. In addition, the intensity of the damage to the dam after the earthquake is also classified according to its crest settlement percentage.

In this work, the average crest settlement at the end of the applied earthquake (PGA = 0.16g) is about 15 cm (Figure IV.25). In addition, the dam height is 18 m and the average height of the deformable moraine layer under the dam is about 6 m. This gives a crest settlement percentage of about 0.6% for a PGA of 0.16g. This point is put in Figure IV.32 (green pentagon point, "Numerical model") with the other historical data. The percentage obtained herein lies within an acceptable range of historical data with a slight overestimation and is qualified as a moderate damage. It may lie in the choice of a signal that is not an actual recording but a synthesized one that overestimates the damaging power of a possible existing motion on site.

To illustrate the validity of this latter aspect, another simulation was carried out using real recorded signals in the Pyrenees mountains in France. It is characterized by the same PGA (0.16g) and magnitude (about 6), an epicentral distance less than 20 km, and a duration of the strong phase of the motion of about 4 s. A crest settlement of about 4 cm was obtained giving a settlement percentage of about 0.16%. This new point was added with the other historical data in Figure

Dynamic safety factor assessment

The seismic safety factor that was previously defined in Chapter III is intimately related to a simplified approach (pseudo-static) of the seismic behavior of structures as proposed by Eurocode 8 [START_REF]eurocode 8: design of structures for earthquake resistance-part 1: general rules, seismic actions and rules for buildings[END_REF]. The effect of a time-dependent loading on the structures is simplified as an equivalent inertial constant force acting at the center of gravity of the moving masses. This approach is a "low frequency" approach which is known to be conservative [START_REF] Savalle | Experimental and numerical studies on scaled-down dry-joint retaining walls: Pseudo-static approach to quantify the resistance of a dry-joint brick retaining wall[END_REF]. Other approaches to evaluate the safety margin of Escoubous dam against a seismic motion are proposed. Note that the dynamic motion used in this safety assessment study is based on the artificially generated signals in order to have more conservative safety factors.

Parameters degradation

The first approach that was used to evaluate the dam's safety during the dynamic loading of PGA of 0.16g is the parameter's degradation method which is similar to the ones used in the static study. The friction angles of the dam body and at the contacts between the pitching blocks were degraded progressively until reaching failure (Equations IV.9 and IV.10). Then, the safety factor F φ ,dynamic is considered equal to the degradation coefficient at failure.

tan(φ rock f ill ) = tan(φ rock f ill-initial ) F φ ,dynamic (IV.9) tan(φ pitching ) = tan(φ pitching-initial ) F φ ,dynamic (IV.10)
More in detail, the calculation is performed through two steps.

a) first, the static equilibrium is checked for a certain increment for the degradation of mechanical properties of the dam body and of the pitching.

b) secondly, if the static equilibrium is tested, a dynamic computation is performed to check if the system in that situation is able to withstand the shaking. The minimum degradation state associated with a failure throughout the shaking is defined as the dynamic safety factor.

Several degradation coefficients for the dam strength were tested. The dam failed during the dynamic loading with a degradation coefficient of 1. A static post-shaking safety factor can be computed in order to estimate the degradation of the dam's safety after the shaking. It complements the previous one but is truly a static factor. The procedure led to a static post-shaking factor of 1.3 which is slightly smaller than 1.4 obtained before the shaking. The seismic motion earthquake has led to a slight decrease in the static safety of the dam due to the earthquake-induced deformations and damages. However, the state of the system after the shaking is certainly not the one obtained herein. We expect that cyclic deformations induce the densification of the rockfill which is supposed to increase both the hardening of its mechanical behavior (stress-strain curve, lower contracting and appearance of a possible dilative behavior) and its ability to resist shearing (an increase of internal friction angle). In the absence of an elasto-plastic constitutive law for the rockfill valid in a dynamic context, it is difficult to state about this reduction of the static safety factor.

PGA augmentation

The second approach that was used to evaluate the dam's dynamic safety factor involves performing sequentially dynamic calculations with PGA values that increase. This approach may be seen as the one that holds most a physical meaning since it is done by scaling up the applied signals on the foundation base by a certain scaling factor and then running the dynamic modeling of the dam (New PGA = scaling factor × 0.16g) using the scaled-up signal. The safety factor is defined to be equal to the smallest value of the scaling factor inducing the dam failure during the simulated earthquake. Finally, the dynamic safety factor using this approach is in the range between 3.0 and 3.3. This also means that the dam is able to resist seismic motions up to 0.48g (3*0.16g) which is consistent with the findings of Deluzarche [START_REF] Deluzarche | Modélisation discrète des enrochements: Application aux barrages[END_REF]. In addition, from this study, one can note once again the highly conservative results of the pseudo-static approach that predicted a failure acceleration of about 0.14g = tan(8)g (8°is the tilting angle at failure in the pseudo-static study) which is about 3 times lower than what has been found in the true dynamic study.

Conclusion

The dynamic behavior of Escoubous Dam, located in the "zone 4" seismic region, was studied using the mixed FDM-DEM numerical approach. Synthesized seismic motions were used with horizontal (PGA =0.16g) or both horizontal and vertical components.

LKE model was discarded due to its lack of robustness and validity in a dynamic context to the benefit of an incremental hysteric model coupled with the Mohr-Coulomb model for the rockfill. The limits of this approach are that irreversibilities are only generated when the failure criterion is reached and there is no effect of the cycles on the values of the mechanical properties. In particular, there is no possibility for a hardening due to cycles of the mechanical behavior (both stress-strain and densification curve).

Some lessons were taken from the performed studies. The moraine layer in the foundation and under the dam definitely amplifies the motion by about 25% (0.2g signal arriving at the dam base after passing by the moraine layer compared to the 0.16g outcrop signal). Thus, due to this amplification, the presence of this moraine layer (heterogeneous foundation) has led to higher dam deformations when compared to the rocky homogeneous foundation case.

Taking into account the vertical component of the seismic motion amplifies the phenomena obtained using only a horizontal motion. This aspect was expected.

The dam crest settlement was compared to other cases where settlements of actual dams were measured after actual seismic events. The case involved in this work seems to be within the range of observed results but on an upper bound. It may be due to the choice of the synthesized seismic motion which is built on very conservative statements. With the usage of a real recorded signal, deformations are closer to and within the range of the observed historical results. The damage for Escoubous dam in the studied case is considered as moderate compared to what was reported in the literature. A shear surface has developed on the downstream side of the dam stretching from the middle of it to the lower downstream side with local values up to 12%. Note that no hydrodynamic pressure was considered which is a limit of this study (Hydrostatic pressure was used).

The acceleration of the lower downstream pitching was found very high up to 7g. It acts as a beam that is embedded at its base in the foundation where it receives the motion signal and is almost free to move at its top part which is not well constrained by the berm blocks. This leads to high amplifications in this part of the pitching. For the upper downstream pitching, the signal propagating in the lower downstream pitching doesn't reach it since it is cut by the berm. Its acceleration is much lower than that found in the lower downstream pitching and is mainly due to the motion that it receives from the dam body at its interface. Therefore, the 2 parts of the downstream pitching were found independent in the dynamic behavior.

The motion seems to have contributed to the destructuring of the berm and at the same time seems to have smoothened the contacts between the lower downstream pitching and the dam body. Some of the initially lost interface contacts were recovered and other active contacts were lost. In total, the number of active interface contacts in the lower downstream pitching increased by 4%. The shaking has changed its interface contact distribution. In addition, the higher acceleration in the lower downstream pitching could also have a role in this pitching rearrangement.

The post-seism static safety factor was found of 1.3 which is lower than the initial value equal to 1.4. The seismic motion led to the weakening of the dam but to a small extent. The second approach was based on the upscaling of the applied signal and its PGA. A safety factor ranging between 3.0 and 3.3 was obtained.

General conclusion and perspectives

A mixed continuum-discrete (FDM-DEM) approach for the modeling of rockfill dams with drystone pitching was built. The continuum approach has been used to model the rockfill forming the dam body while the discrete method has been used to model the discrete anisotropic nature of the pitching.

The mixed FDM-DEM approach was successfully validated for the modeling of the behavior of such hydraulic structures with the help of available results obtained through a past experimental campaign of pseudo-static tests performed a few years ago on different down-scaled dams holding the particularities of the dams studied herein. Two constitutive models were tested for the continuum modeling of the dam body: the typical Mohr-Coulomb model and an elastoplastic model (LKE) for rockfill developed by EDF. The calibration of the model parameters was performed with a straightforward procedure. The values of the tilting failure angles were retrieved with accuracy using both constitutive models but LKE model was the only constitutive model capable of reproducing and fitting the displacements obtained throughout the experiments. Two different safety factors for the different structures were also estimated and revealed the key role of the pitching in their stability. At this scale, the safety factors were found in the range of 2.1 to 3.7.

From the comparison of the different cases, it has been seen that the weight of the pitching and the friction angle of the pitching-backfill interface have a significant effect on the stability of the dams. In addition, it has been found that the higher the pitching weight, the stiffer the behavior of the dam. Finally, simulations were also performed with blocks for the pitching that do not hold a perfect parallelipipedic shape which condition is closer to actual blocks contrary to ideal shapes that are generally used in such studies. It has been found that the geometric imperfection induces a relative mobility between the blocks of the pitching and decreases the resistance of the dam.

The modeling of a rockfill dam with dry-stone pitching, built in the 1950s and located in the Pyrenees in France, has been carried out using the mixed FDM-DEM approach that was previously validated. The model parameters were extrapolated from previous studies in the literature. The validity of the set of parameters has been confirmed as it gave displacements after water filling which are similar to historical values. The dam's safety was assessed and found to be in the range of 1.3-1.4. To better understand the role of the different elements of the site and the dam in its resistance, different case studies were carried out. The key role of the pitching was here again since the dam was not able to be built without it. The berm was found to play an unfavorable role in the dam's stability. In fact, it cuts the downstream pitching into two independent parts preventing any mechanical transfer of forces from the top one to the lower one. The layer of more dense rockfill under the pitching seems to compensate this aspect however a computation involving a less resistant layer did not prove to decrease significantly the dam safety. The upstream concrete mask and the existence of a heterogeneous foundation had a slight effect on the dam displacements but no effect on the dam's safety.

Finally, an investigation of the downstream pitching deformations and bulges that exist on the actual hydraulic structure was carried out. These deformations were qualitatively retrieved together with their location on the pitching but not quantitatively. The existence of the chamber valve which acts as a hard point in the system increased the obtained deformations in the bottom part of the downstream pitching. A creep behavior together with successive filling/emptying phases could lead to more severe deformations and then to values closer to what is observed on site. In addition, the 3D effect linked to the shape of the valley could also contribute to an increase in the values obtained in this work.

The case study was amplified by studying the behavior of this dam under seismic loadings. First, artificial signals with a PGA of 0.16g, provided by EDF from an internal report, were used. In this part of the study, LKE model was discarded due to its lack of robustness and validity in dynamic computations. Instead, a hysteretic model coupled with the MC model has been used for the rockfill. Several cases were studied, leading to the following conclusions. It has been shown that the moraine layer under the dam amplifies the motion reaching the dam base by about 25% which leads to higher displacements and accelerations in the dam body. The consideration of the vertical component of the seismic motion amplifies the phenomena obtained using only a horizontal motion. In addition, the presence of the hydrostatic pressure on the upstream face of the dam due to the existence of the reservoir tends to play a stabilizing role on this dam side whereas it increases the displacements on the downstream side. The total displacements computed at the crest at the end of the shakings were found to be in the range of the statistics of displacements observed for this kind of dam but on the upper limits. This could be due to the choice of a synthesized signal which is conservative. A further computation using a recorded motion with equivalent peak ground acceleration and epicentral distance led to displacements closer to historical data.

Concerning the downstream pitching, its lower and upper parts clearly behaved independently throughout the shaking due to the cutting effect of the berm. Higher accelerations were found in the lower pitching which behaved like a beam embedded in the foundation, and generally much higher than those found in the upper pitching. The shaking has also led to the rearrangement of the lower downstream pitching with an averaging effect on the interface contacts distribution with the dam body.

On the safety side, the post-seismic safety of the dam was found equal to 1.3 which is just slightly lower than the initial value of 1.4. The seismic motion weakened the dam with damages considered moderate. Finally, the artificial signal was upscaled to study the dam's seismic resistance limits. It failed at a scaling factor of 3.3 and was stable at a value of 3 which concludes that the dam can handle seismic motions up to about 0.5g.

More generally, the work in this thesis has enabled the validation of the mixed FDM-DEM approach in the modeling of rockfill dams with dry-stone pitching and showed its ability to correctly retrieve the behavior of such dams in reasonable computational times. In addition, this newly validated numerical tool adds the possibility to refine the modeling of the pitching including the block shapes and their arrangement among others.

In the dynamic study, a Mohr-coulomb model coupled with a hysteretic damping model was used due to the lack of robustness of LKE model. The limits of the used constitutive models are that they can not generate plastic deformations when unloading after the PGA and do not correctly model phenomena at stake in granular materials throughout cyclic loadings. Moreover, there is no effect of the cycles on the values of the mechanical properties including hardening of the stressstrain curve and the rate of densification in drained conditions. Further work is required to adapt LKE model to seismic computations. It includes a better description of cyclic phenomena in granular materials and robustness.

In general, a better understanding of the dam structure, both in terms of geometry and properties of materials, would be beneficial: a better knowledge of the materials composing the dam can be acquired by doing more in-depth in-situ and lab tests which would allow to better estimate and calibrate the model parameters. Furthermore, the geometry of the pitching blocks and thickness of the rearranged rockfill sub-layers can be investigated on site in order to have more precise data. This would lead to a more accurate numerical model which could pave the way to more precise and robust results.

Finally, a statement of plane strain deformations was taken in this work for the case study involving an actual dam. This statement was also considered in previous work involving this dam. The modeling of the true V-shaped valley may give additional and more precise information. In particular, it may induce a 3D effect in the pitchings and may amplify the detachment of the pitchings creating more deformations of it than those found in this work. a) Nodal force: the resolution process begins by determining the unbalanced forces F i at a given node. They are the sum of the contributions of the various tetrahedra common to this node, each contribution being the sum of three terms: one related to the stress fields on the surfaces of the tetrahedron T i , the other related to the volume forces ρb i and the last characterizing an external loading P i applied to the considered node:

F <n> i = [[ T i 3 + ρb i V 4 ]] <n> + P <n> i (B.1)
The exponent < n > indicates the value of node n in the global numbering of all nodes in the mesh, to be distinguished from the exponent l (l=1,4) which refers to node l in a tetrahedron.

The symbol [[•]] <n> is used to represent the sum of the contributions to node n of all the tetrahedra that meet at this node. This force tends towards zero when the problem reaches a state of equilibrium.

b) Nodal velocity: the new nodal velocity values v i are determined from the indefinite Newton equation which is expressed in nodal terms, leading to an equation of the following form for the node:

F <n> i = M <n> dv i dt <n> (B.2)
where the nodal mass M <n> of the node is equal to a quarter of the tetrahedral mass (ρV /4) of the node.

Considering that the nodal velocity varies linearly during the calculation step ∆t, we can write:

v <n> i (t + ∆t 2 ) = v <n> i (t - ∆t 2 ) + ∆t M <n> F <n> i (B.
3) c) Deformation rate ξ i j : the deformation rate of the various elements is determined using Gauss's theorem and the assumption of velocity linearity and constant deformation rate within each tetrahedron:

ξ i j = - 1 6V 4 ∑ l=1 (v l i n (l) j + v l j n (l) i )S (l) (B.4)
The exponent < l > indicates the value associated with the face l (l=1,4) of each tetrahedron. Note that, by convention, face l is opposed to node l in a tetrahedron normal to face l and n (l) is the normal unit vector of the l face. d) Stress tensor ∆σ i j : the new stress tensor is deduced from the strain rate using the materials law of behavior, whose incremental expression is in the following form : ∆σ i j = H * i j (σ i j , ξ i j , α) (B.5)

The function H * i j characterizes the behavior law. Note that for large-strain calculations, a corrective term linked to the rotation rate tensor is added to the resulting stress to take into account the rotation of the principal stress axes: ∆σ i j = H * i j (σ i j , ξ i j , α) + (ϖ ik σ k jσ ik ϖ k j )∆t (B.6)

where the components of the rate-of-rotation tensor ϖ ik are calculated by:

ϖ ik = - 1 6V 4 ∑ l=1 (v l i n (l) j + v l j n (l) i )S (l) (B.7)
This elementary algorithm is reproduced until an equilibrium state is reached, i.e. a state for which the unbalanced nodal forces are less than a fixed calculation tolerance. The latter is called a "ratio" in FLAC3D, and corresponds to the ratio of the average/maximum unbalanced force on the model by the initial average/maximum unbalanced force. However, this stability criterion leads to a very high number of calculation cycles, making the explicit method not very efficient for linear problems, but has the advantage when studying nonlinear and dynamic problems. Note that in the case of large strain calculations, new nodal coordinates can then be calculated and the mesh is updated.

Mechanical damping

During the first calculation cycles, the system is unstable, but the stress and strain paths are calculated at each time step. To create an operational algorithm, the solid's movements must be damped as rapidly as feasible to achieve a condition with negligible residual imbalance. Damping in FLAC3D is accomplished by exerting a force on the nodes known as the damping force, the modulus of which is proportional to the unbalanced force from the previous computation cycle, and whose direction is such that it always produces dissipative work. Local damping and combined damping are two dampening types that can be used in FLAC3D.

Local damping

The damping force is expressed as a function of the unbalanced force F i and nodal velocity u˙i :

F d = -α|F i |sign(u˙i) (B.8)
where α is the damping constant taken by default equal to 0.8 and the sign function is defined as: sign(x) = -1 if x<0 and sign(x)=1 if x>0.

Combined damping

The combined damping mode gives equal importance to the change of sign of nodal velocities as that of unbalanced forces, in the residual energy dissipation function over time. The damping force is then:

F d = α 2 |F i |[sign( dF i dt
)sign(u˙i)] (B.9)

Timestep (∆t)

In the finite difference scheme, the difference equations (Equation B.3) will not provide valid answers unless the numerical scheme is stable, i.e. the used timestep is less than the critical timestep. In FLAC3D, during mechanical static modeling, an important aspect of the numerical method is that a uniform unit timestep ( ∆t=1) is adopted for the whole system.

At gridpoints in dynamic problems, "real" masses are used instead of the fake masses used in the static solution scheme for optimum convergence. Each tetrahedral subzone provides onequarter of its total mass (determined by zone density and area) to each of the four corresponding gridpoints. In the instance of an eight-noded zone with two overlays, the total gridpoint mass is then split by two. FLAC3D employs lumped masses and a diagonal mass matrix in finite-element language.

The dynamic timestep is then calculated as follows:

∆t crit = min V C p A f max (B.10)
where C p is the p-wave speed, V is the tetrahedral subzone volume, and A f max is the maximum face area associated with the tetrahedral subzone. The min {} function is taken over all zones, and includes contributions from the structural and interface modules. Since Equation B.10 is just an estimation of the critical timestep, a safety factor of 0.5 is taken, i.e. the used timestep ∆t dynamic is half the critical timestep ∆ crit (∆t dynamic = ∆ crit /2).

2. PFC3D: DEM modeling 2.1. Overview of PFC3D PFC3D (Particle Flow Computation) is a software developed by ITASCA that simulates discretely the interactions between rigid finite-sized particles using the distinct element method (DEM).

The model is made up of distinctive particles that move apart from one another while interacting during pair-wise contacts. The mechanical behavior of such a system is characterized in terms of the motion of each particle and the inter-particle forces acting at each contact location where the particles are assumed to be rigid. The essential connection between particle motion and the forces driving it is provided by Newton's equations of motion. The force system may be such that the particles flow, or it may be in static equilibrium (in this case there is no motion). If the particle-interaction law predicts a physical contact between particles, the behavior of the contact is described using a soft contact approach, in which the rigid particles are permitted to overlap in the vicinity of the contact point and a finite normal stiffness is taken to represent the measurable stiffness that exists at a contact. More sophisticated behavior can be reproduced by bonding the particles at their contact locations so that when the inter-particle forces acting at any bond reach the bond strength, the bond is broken. Tensile forces can form between particles as a result of this. One may then simulate the interaction of these bonded "blocks," including the emergence of fissures that can lead blocks to fragment into smaller blocks. Potential energy functions can also be used to derive the particle-interaction law and represent long-range interactions.

DEM Calculation cycle in PFC3D

PFC is an explicit type of code, i.e. the model configuration evolves (the simulation progresses) over time according to the resolution of successive calculation cycles. Calculation cycles stop when a user-defined stopping criterion is reached. The observed behavior is the direct result of the interaction of the various bodies in contact. For each cycle, a series of operations is performed (

Contact models

Several DEM contact models between rigid bodies have been developed and added to PFC3D. The 2 contact models that have been used in this project are presented in the following sections.

Linear frictional model

A linear frictional model is used to model the contacts between the blocks and at the interface. This model is composed of two components that act in parallel with each other: the first one is the linear component which provides the linear elastic frictional behavior with no tension whereas the second is the dashpot component which provides viscous behavior (Equation B.11).

Each component is split into normal and shear forces (Equation B.12). The linear force is calculated using the linear springs with constant normal and shear stiffnesses k n and k s . The dashpot force is calculated using dashpots with normal and shear critical damping ratios β n and β s . The contact can be considered as active or inactive according to the gap between two different bodies.

The surface gap g s is the difference between the contact gaps g c and the reference gap g r . The contact is active only if the surface gap is less than or equal to zero. The contact interface consists of a contact plane with certain position and normal direction. ). The force is updated according to the relative displacement increment. The linear normal force is updated based on the normal displacement (surface gap) increment ∆δ n . Moreover, no tension forces exits in this linear model. The linear shear force is updated based on the relative shear displacement increment ∆δ s . This shear force can increase with the relative shear displacement until a threshold then slipping behavior occurs (when the shear force exceeds the shear strength F u s = -µ.F l n where µ is the defined friction coefficient and F l n is the contact normal force). If the contact changes from being inactive to active during the current timestep, just the amount of this increment that takes place while the surface gap is negative is employed to carry out the incremental update of the normal and shear forces. The dashpot forces are also updated based on the dashpot mode and incremental equations. Detailed description can be found in ITASCA PFC3D documentation [START_REF]CGI Itasca. Pfc3d (particle flow code in 3 dimensions), version 4.0[END_REF]. The other components correspond to the building blocks of the "Parallel Bond" model: normal ( kn ) and tangential ( ks ) stiffness, tensile strength ( σc ), local cohesion ( c) and angle of friction ( φ ). When the bond breaks, the contact model will be purely linear as the one described above. with K e 0 and G e 0 reference values for an effective mean pressure equal to P a = 100kPa. I 1 is the first invariant of the effective stress tensor and n elas is a model parameter.

• Plastic isotropic mechanism

Concerning the plastic isotropic mechanism, the yield surface of the isotropic plastic mechanism is a plane orthogonal to the hydrostatic axis in the space of principal stresses, expressed by:

f i (σ ′′ ; p c ) = I 1 3 -p c (C.6)
The hardening law is considered since the consolidation pressure p c changes according to Equation C.7. p c0 and β are model parameters and the volumetric plastic deformation ε p v is the hardening variable. The hardening variable of the elastoplastic mechanism is the deviatoric plastic deformation γ d defined by: For each parameter a and m involved in equations C.9 and C.10, specific variation laws are defined according to the hardening variable γ d : the elastoplastic surface thus evolves throughout different reference state, namely the peak state and the residual state at large deformations. The different hardening phases are bounded by "threshold" values, noted γ peak and γ res .

γ d = ( 2
In the pre-peak phase, the deviatoric yield surface evolves from the elasticity threshold to the peak threshold, for 0<γ d <γ peak , according to the following laws (X ams is a model parameter): The deviatoric deformations are obtained from the following expression:

ε ḋ i j = λ ḋG i j (C.18)
Where G is the hardening function of the elastoplastic mechanism defined by: [START_REF] Frossard | Comportement macroscopique des matériaux granulaires mis en oeuvre dans les barrages[END_REF])

G i j = ∂ f d ∂ σ i j -( ∂ f d ∂ σ kl n kl )n i j (C.
n i j = β ′ s i j s II -δ i j β ′2 + 3
, where β ′ = √ 6g v 3g v and g v = - Identifying those model parameters usually needs experimental data from triaxial and oedometric tests. The main problem in this process is the lack of experimental data to be based on due to difficulty in carrying out experiments on large samples of rockfill or due to the lack of access to this data. To overcome the lack of data, the literature and previous studies can be helpful in identifying those parameters. Then, the value of p c0 can be deduced from this line at a given relative density D r corresponding to a void ratio e 0 . More details about this approach can be found in its dedicated article [START_REF] Hicher | Micro-macro correlations for granular media. application to the modelling of sands[END_REF].

iii I. 28

 28 Experimental triaxial tests of the rockfill at a confining pressure of 100 kPa and 200 kPa . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . I.29 Dry-stone pitching of the four physical models in PEDRA experimental campaign I.30 Dry stone pitching deformations of case 4 at bucket rotation angle of 27°(just before failure) and sensor time displacements as a function of the bucket rotation angle . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . I.31 Definition of four breakable block shapes (type from 1 to 4, from left to right). . . I.32 (a) Actual dam; and (b) fully 2D DEM model of the dam with the different zones I.33 Failure of the dam during its construction, without pitching . . . . . . . . . . . . I.34 Displacements after filling the reservoir (unbreakable blocks) . . . . . . . . . . . I.35 Displacements after filling the reservoir (breakable blocks) . . . . . . . . . . . . I.36 Horizontal displacement of the crest during filling and emptying of the reservoir . I.37 Elements and geometry of the Escoubous dam model (©EDF-2009) . . . . . . . I.38 Vertical displacements after reservoir filling (©EDF-2009) . . . . . . . . . . . . I.39 Horizontal displacements after reservoir filling (©EDF-2009) . . . . . . . . . . . I.40 Experimental loading of dry-stone retaining walls: a) hydrostatic pressure done by Villemus [48]; b) backfill loading done by Colas [47] . . . . . . . . . . . . . . I.41 Numerical models developed by Oetomo [58] for the studied dry-stone retaining wall: a) Full DEM model; b) Mixed FDM-DEM model . . . . . . . . . . . . . . II.1 Numerical models of the four PEDRA cases . . . . . . . . . . . . . . . . . . . . II.2 Dimensions of the dry-stone pitching used in the numerical model . . . . . . . . II.3 Boundary conditions of the numerical model . . . . . . . . . . . . . . . . . . . II.4 Example of the modeling procedure and steps (displacement fields shown as examples with only qualitative meaning) . . . . . . . . . . . . . . . . . . . . . . . II.5 Mesh study of the FLAC3D model . . . . . . . . . . . . . . . . . . . . . . . . . II.6 Computation scheme involving equilibrium and failure criteria . . . . . . . . . . II.7 The rockfill used in PEDRA experimental campaign in addition to its particle size distribution curve . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . II.8 LEP'S diagram: Effect of the confining pressure on the friction angle . . . . . . II.9 Views of the two yield surfaces (LK Enroch model): in the meridian (left) and deviatoric (right) planes. S II is the 2nd Invariant of the deviatoric stress tensor . . II.10 Experimental triaxial tests and simulations using MC and LKE models: deviatoric stress (right) and volumetric strain (left) at 100kPa and 200kPa . . . . . . . . . . II.11 Failure envelopes of calibrated MC and LKE models at low pressures (5 kPa) in the deviatoric plane . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . II.12 Slide test involving stones to identify friction angle . . . . . . . . . . . . . . . . II.13 Wall zone coupling: a) stone pitching using infinitely rigid blocks, b) coupled walls at the interface, c) dam body . . . . . . . . . . . . . . . . . . . . . . . . . II.14 Modeling of the anchorage in case 3 . . . . . . . . . . . . . . . . . . . . . . . . II.15 Failure angles for the four cases using MC and LKE model . . . . . . . . . . . . II.16 Displacements at the 6 sensors on the downstream face; case 4 . . . . . . . . . . II.17 Displacements contour of the backfill and stone pitching for case 4 at tilting angle of 27°(just before failure): a) using MC, b) using LKE and c) experimentally . . II.18 Normal relative displacements evolution of the pitching throughout the tilting tests; simulation of case 4 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . x II.19 Simulations of PEDRA tests: normal normalized relative displacements for a tilting angle just before failure; case 1: 24°; case 2: 22°; case 3: 29°; case 4: 27°. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . II.20 Displacement increments between the last stage before failure in the static safety factor calculation by parameters degradation and the initial static equilibrium, for Case 1 with and without water pressure . . . . . . . . . . . . . . . . . . . . . . . II.21 Irregular blocks used with imperfections on the:a) bottom side, b) upper side and c) upper+bottom side of the block . . . . . . . . . . . . . . . . . . . . . . . . . II.22 Numerical model view of the irregular blocks with imperfection on the: a) upper side, b) bottom side and c) upper+bottom side of the block . . . . . . . . . . . . II.23 Normal displacements of the downstream pitching at rotation of 20 degrees for the 3 studied irregular blocks . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . II.24 Normal displacements of the downstream pitching at different tilting angles for different contact surface percentages (irregular blocks with bottom imperfection) II.25 Displacement increments of the downstream pitching (composed of irregular blocks with bottom imperfection (20% contact)) at different tilting angles . . . . . . . . II.26 View of the modeled vertical imperfections in the vertical contact planes between the blocks: (a) gap over the entire block height (no contact), (b) gap with partial contact . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . II.27 Normal displacement increments of the downstream pitching at 24 degrees of rotation (just before failure) in the cases with gaps in the vertical plane contacts between the blocks . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . III.1 FLAC3D-PFC3D coupled numerical model of Escoubous dam and representation of the various elements . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . III.2 Dimensions of the numerical model and its different elements . . . . . . . . . . . III.3 Results of a triaxial compression test for 100 kPa confinement (LKE model) -Deviatoric stress (left) / Volumetric strain (right) . . . . . . . . . . . . . . . . . . III.4 The 4 stages of the dam construction . . . . . . . . . . . . . . . . . . . . . . . . III.5 Displacement field at the end of the construction phase: a) total displacement; b) horizontal displacement; c) vertical displacement . . . . . . . . . . . . . . . . . III.6 Principal stress directions and shear stress field in the dam body after construction, and in particular in the hand-rearranged rockfill sub-layer under the pitching . . . III.7 Normal forces, shear forces and mobilized friction ratio between the blocks, within the pitching and after construction . . . . . . . . . . . . . . . . . . . . . . . . . III.8 Normal and shear forces at the pitching-backfill interface contacts . . . . . . . . III.9 Force system applied to upstream pitching stones at the end of the reservoir water filling . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . III.10Horizontal (left) and vertical (right) displacements of the dam after reservoir filling III.11Shear strain inside the dam body throughout the water filling process . . . . . . . III.12Normal forces, shear forces and mobilized friction ratio between the blocks, within the pitching after water reservoir filling . . . . . . . . . . . . . . . . . . . . . . . III.13Normal and shear forces and contacts at the pitching-backfill interface. . . . . . . III.14Intensity and direction of displacements in the body of the dam and the pitching, after construction (left) and after water filling (right) . . . . . . . . . . . . . . . . III.15Displacement magnitude at the end of the safety factor evaluation using the parameter reduction method with water pressure (F φ = 1.4) . . . . . . . . . . . . . III.39Shear strains filed after the water filling phase using: (a) parallelepipedic blocks; (b) trapezoidal blocks (with imperfections) . . . . . . . . . . . . . . . . . . . . . III.40"Active" contacts at the backfill-pitching interface, at the end of construction and after water filling using: (a) parallelepipedic blocks; (b) trapezoidal blocks (with imperfections) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . III.41a) Bulges found on the lower part of the downstream pitching around the valve chamber; b) Numerical modeling of the effect of the valve chamber. . . . . . . . III.42Normal displacement of the downstream pitching (top) and contacts at the interface (bottom) after water filling: (a) with a rigid structure (valve chamber); (b) without a rigid structure . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . III.43Normal forces (top), shear forces (middle) and mobilized friction (bottom) between pitching blocks at the end of construction (left) and after water filling (right). III.44Two sections of the pitching after construction and at tilting angle = 6°: 1) at exactly the middle of the model and 2) just 1 m far from the dam lateral boundary III.45Normal displacement of the downstream pitching with the rigid structure throughout the rotation process . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . IV.1 Seismic zones distribution in France . . . . . . . . . . . . . . . . . . . . . . . . IV.2 Dynamic boundary conditions of the numerical model . . . . . . . . . . . . . . . IV.3 Shear wave speed (V s ) and shear modulus profiles along the dam height: numerical model and Gazetas estimation . . . . . . . . . . . . . . . . . . . . . . . . . IV.4 Gazetas approach for the estimation of the shear modulus profile along the dam height . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . IV.5 Theoretical and numerical (Hysteretic) shear modulus degradation (left) and damping (right) curves for the dam rockfill . . . . . . . . . . . . . . . . . . . . . . . . IV.6 Theoretical and numerical (Hysteretic) shear modulus degradation (left) and damping (right) curves for the moraine and alluvium . . . . . . . . . . . . . . . . . . IV.7 Spectra of the five different generated signals with the target spectrum . . . . . . IV.8 Five different signals provided by EDF and with a PGA of 0.16g: acceleration (left), velocity (middle), and displacement (right) versus the dynamic time . . . . IV.9 Arias intensity of signal 1749 . . . . . . . . . . . . . . . . . . . . . . . . . . . . IV.10Choice of the cut-off frequency and its effect on the 1749 signal and its acceleration amplitude spectrum . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . IV.11Maximum shear strain in the dam body throughout the dynamic loading using an elastic model without damping . . . . . . . . . . . . . . . . . . . . . . . . . . . IV.12Shear degradation and damping values for the effective cyclic shear strain throughout the dynamic loading using an elastic model without damping . . . . . . . . . IV.13Target (1749 signal) vs numerically measured outcrop signal spectrum . . . . . . IV.14Horizontal and vertical displacements at dam crest at the end of the horizontal shaking; (a): case A ; (b) case B. . . . . . . . . . . . . . . . . . . . . . . . . . . IV.15Horizontal and vertical displacements at dam crest throughout the horizontal shaking for homogeneous (Case A) and heterogeneous (Case B) foundation cases . . IV.16Cumulative shear strain at the end of the horizontal shaking: (a): case A ; (b) case B.114 IV.17Traced blocks throughout the horizontal shaking; (a) upstream pitching ; (b) downstream pitching. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . IV.18Horizontal and vertical displacements of the traced pitching blocks throughout the horizontal shaking; (a) case A; (b) case B. . . . . . . . . . . . . . . . . . . . . . xiii IV.19Horizontal acceleration at the dam base and crest throughout the horizontal shaking; (a) case A; (b) case B . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . IV.20Horizontal and vertical displacements of the studied dam at the end of the shaking having: (a) a horizontal component (case B); (b) both a horizontal and vertical component (case C). . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . IV.21Horizontal and vertical displacements at the dam crest throughout the shaking having: (a) a horizontal component (case B, green curves); (b) both a horizontal and vertical component (case C, red curves). . . . . . . . . . . . . . . . . . . . . IV.22Shear strain of the studied dam at the end of the shaking having: (a) a horizontal component (case B); (b) both a horizontal and vertical component (case C). . . . IV.23Horizontal and vertical displacements of the traced pitching blocks throughout the shaking having: (a): a horizontal component (case B); (b): both a horizontal and vertical component (case C). . . . . . . . . . . . . . . . . . . . . . . . . . . . . IV.24Horizontal and vertical displacements of the studied dam at the end of the horizontal shaking; (a) without hydrostatic pressure (case C); (b)with hydrostatic pressure (case D) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . IV.25Influence of the hydrostatic pressure on both the horizontal and vertical displacements at the dam crest throughout the dynamic loading. . . . . . . . . . . . . . . IV.26Shear strain of the studied dam at the end of the horizontal shaking; (a): without hydrostatic pressure (case C); (b): with hydrostatic pressure (case D). . . . . . . . IV.27Influence of the hydrostatic pressure on the displacements at monitored points of the pitching throughout the dynamic loading;(a): without hydrostatic pressure (case C); (b): with hydrostatic pressure (case D). . . . . . . . . . . . . . . . . . . IV.28Horizontal acceleration at four points in the downstream pitching throughout the shaking. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . IV.29Horizontal acceleration comparison between the dam body, the sub-layer and the stone pitching: (a) rblock 1 (upper downstream pitching) and (b) rblock 3 (lower downstream pitching) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . IV.30Normal and shear forces and mobilized friction between the pitching blocks before and after the dynamic loading . . . . . . . . . . . . . . . . . . . . . . . . . . . . IV.31Normal and shear forces and mobilized friction at the interface (between the pitching and the dam body) before and after the dynamic loading . . . . . . . . . . . . IV.32Historical percentage of settlements of different types of embankment dams as a function of the peak ground acceleration . . . . . . . . . . . . . . . . . . . . . . IV.33Failure of the dam for a reduction factor equal to 1.3 at 8s and 9s of dynamic loading126 IV.34Failure of the dam at a scaling factor of 3.3 (PGA = 3.3 x 0.16g) at 12s and 13s of dynamic loading . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . A.1 Dry stone pitching deformations of case 1 at bucket rotation angle of 23°(just before failure) and sensor time displacements as a function of the bucket rotation angle . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . A.2 Dry stone pitching deformations of case 2 at bucket rotation angle of 20°(just before failure) and sensor time displacements as a function of the bucket rotation angle . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . A.3 Dry stone pitching deformations of case 3 at bucket rotation angle of 29°(just before failure) and sensor time displacements as a function of the bucket rotation angle . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . B.1 Parallelepiped (a) and tetrahedron (b) elements in FLAC3D . . . . . . . . . . . . xiv B.2 Calculation sequence performed at each time step ∆t in FLAC3D . . . . . . . . . B.3 Sequence of DEM primary operations that occur during each cycle in PFC3D . . B.4 Behavior and rheological components of the linear model . . . . . . . . . . . . . B.5 Force-displacement law for the linear component of the unbonded linear-based models: (a) normal force versus surface gap, (b) shear force versus relative shear displacement, and (c) slip envelope . . . . . . . . . . . . . . . . . . . . . . . . B.6 Behavior and rheological components of the parallel bond model . . . . . . . . . B.7 Force-displacement law for the parallel bond force and moment: (a) normal force versus parallel-bond surface gap; (b) shear force versus relative shear displacement; (c) twisting moment versus relative twist rotation; and (d) bending moment versus relative bend rotation. . . . . . . . . . . . . . . . . . . . . . . . . . . . . B.8 Failure envelope for the parallel bond . . . . . . . . . . . . . . . . . . . . . . . . C.1 MC failure envelop in the deviatoric plane . . . . . . . . . . . . . . . . . . . . . C.2 Simulation of a compression triaxial test using MC model[98] . . . . . . . . . . C.3 Deviatoric triaxial behavior using LKE model . . . . . . . . . . . . . . . . . . . C.4 LKE deviatoric failure envelop in the deviatoric plane . . . . . . . . . . . . . . . D.1 (log(p c ), D r ) plane . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . E.1 Displacements at the 6 sensors of the downstream face (case 1) . . . . . . . . . . E.2 Displacements at the 6 sensors of the downstream face (case 2) . . . . . . . . . . E.3 Displacements at the 6 sensors of the downstream face (case 3) . . . . . . . . . . E.4 Normal relative displacements evolution of the pitching throughout the tilting tests (case 1) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . E.5 Normal relative displacements evolution of the pitching throughout the tilting tests (case 2) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . E.6 Normal relative displacements evolution of the pitching throughout the tilting tests (case 3) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .List of TablesI.1 Rockfill dams built in France from 1940 to 1960 . . . . . . . . . . . . . . . . . I.2 Four cases studied in PEDRA experiments . . . . . . . . . . . . . . . . . . . . I.3 Failure tilting angle for the four studied experimental cases . . . . . . . . . . . . I.4 Critical backfill height obtained from numerical simulations (using mixed DEM-FDM and Full DEM approaches) and full-scale experiments; S: failure by sliding/ T: failure by toppling; The percentage represent the difference from the experiments II.1 Mohr-Coulomb parameters calibration at different confining pressures . . . . . . II.2 LKE peak threshold parameters calibration . . . . . . . . . . . . . . . . . . . . II.3 Calibrated sixteen parameters of LKE model . . . . . . . . . . . . . . . . . . . . II.4 Failure angles for the four cases usingMC and LKE model and their percentages of error (shown between parenthesis) with respect to the experimental results . . II.5 Safety factors for the reference case with no pitching and for the four studied cases II.6 Failure angles and positions of the dam using irregular pitching blocks with different imperfection positions . . . . . . . . . . . . . . . . . . . . . . . . . . . . II.7 Failure angle of the dam using irregular pitching blocks with different contact coefficients (bottom imperfection) . . . . . . . . . . . . . . . . . . . . . . . . . III.1 LKE model parameters for the dam body . . . . . . . . . . . . . . . . . . . . . III.2 MC parameters for the re-arranged rockfill layer . . . . . . . . . . . . . . . . . III.3 MC parameters for the different parts of the foundation . . . . . . . . . . . . . . III.4 Summary of previously studied cases: horizontal (x) and vertical (z) displacements at crest and mid-height (berm) in addition to safety factors . . . . . . . . . . . . IV.1 Main geographical zones and the reference acceleration of the seismic zones in France . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . IV.2 MC parameters for the different parts of the system for dynamic simulations (parameters different from the static model in bold) . . . . . . . . . . . . . . . . IV.3 The different studied cases in the dynamic study using MC + hysteretic model (The condition that changes in each case with respect to the previous one is put in bold) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . C.1 LKE model with its sixteen parameters definition . . . . . . . . . . . . . . . . .

  -stone structures in the Civil engineering sector are shown in Figure I.1 which are the following:

Fig. I. 1 .

 1 Fig. I.1. Different dry-stone structures: a) Rockfill dam with dry-stone pitching on its faces ©EDF-CIH, b) Dry-stone walls retaining roads ©C. Cornu, c) Dry-stone arch bridge in Canada [9], d) Dry-stone walls retaining slopes in Douro valley in Portugal ©C. Cornu

Fig. I. 2 .

 2 Fig. I.2. Sadd-El-Kafara dam: cross section on the left and its dry-stone pitching on the right[START_REF] Garbrecht | Sadd-el-kaffara, the world oldest large dam[END_REF]

Fig. I. 4 .

 4 Fig. I.4. Development of the rockfill dam height throughout the 20th century: the red color represents the period when such dams were built in France [2]

Fig. I. 5 .

 5 Fig. I.5. Examples of typical rockfill dams built with: a) hand-placed rockfill (Italian design) and b) dumped rockfill (USA design), ©EDF

  Fig. I.7. Types of downstream stone pitching (©F. Laigle)

Fig. I. 8 .

 8 Fig. I.8. Deformations observed on the stone pitching: a) Stone cracks, b) Settlements caused by the displacement of cracked stones, c) Dismantling of following extruded stones, d) Bulging of the pitching (©F. Laigle)

Fig. I. 9 .

 9 Fig. I.9. Cross-section of Escoubous dam (©EDF)

  Fig. I.12. Evolution of Escoubous dam crest settlements between 1983 and 2008 (©EDF)

Fig. I. 14 .

 14 Fig. I.14. Granular particles classification as a function of their form[START_REF] Barton | Shear strength of rockfill[END_REF]

Fig. I. 15 .

 15 Fig. I.15. Triaxial tests and plane strain tests (Monterey sand, σ 3 = 68.9kPa)[START_REF] Chen | Modélisation du comportement mécanique des grands CFRD: Identification des caractéristiques des enrochements et comportement du masque d'étanchéité amont[END_REF]

3 )Fig. I. 16 .

 316 Fig. I.16. Evolution of the internal friction angle with the normal stress[START_REF] Sarac | Shear strength of rockfill and slope stability[END_REF]

Fig. I. 17 .

 17 Fig. I.17. Shear strength of rockfill from large triaxial tests[START_REF] Thomas | Review of shearing strength of rockfill[END_REF]

  Figure I.18); • S: Equivalent particle strength as a function of simple compressive strength of a rockfill block (σ c ) and d 50 particle size (Figure I.18); • σ ′ n : normal effective stress. It can be estimated from tests such as the Schmidt hammer or the simple compression test;

Fig. I. 18 .

 18 Fig. I.18. Estimate of the rockfill particles equivalent roughness and strength[START_REF] Barton | Shear strength of rockfill[END_REF]

  Rockfill peak friction angle ( ϕ max in Figure I.19);

•

  Fig. I.19. Rockfill dilatancy friction (ϕ-ψ) versus effective mean stress[START_REF] Md Bolton | The strength and dilatancy of sands[END_REF]

Fig. I. 21 .

 21 Fig. I.21. Rockfill triaxial behavior as a function of the particles angularity[START_REF] Laigle | Proposition de modèles de comportement pour les enrochements et justification sur le cas du barrage de cogoti[END_REF]

Fig. I. 22 .

 22 Fig. I.22. Minimum and maximum void ratios vs the uniformity coefficient and particle form

Fig. I. 23 .

 23 Fig. I.23. Cyclic densification of a loose sand[START_REF] Leslie | Compaction of sands by repeated shear straining[END_REF]

Fig. I. 24 .

 24 Fig. I.24. Nonlinear cyclic stress-strain relationship of soils[START_REF] Yener | A review of considerations on seismic safety of embankments and earth and rock-fill dams[END_REF]

Fig. I. 26 .

 26 Fig. I.26. Example of the PEDRA scaled-down physical model and its dimensions

Fig. I. 28 .

 28 Fig. I.28. Experimental triaxial tests of the rockfill at a confining pressure of 100 kPa and 200 kPa

Fig. I. 29 .

 29 Fig. I.29. Dry-stone pitching of the four physical models in PEDRA experimental campaign

Table I. 3 :

 3 Failure tilting angle for the four studied experimental cases Case 1 Case 2 Case 3 Case 4 Failure angle 24°21°29.5°28°P itching angle (Failure angle + dam slope (45°)) 69°66°74.5°73°F ig. I.30. Dry stone pitching deformations of case 4 at bucket rotation angle of 27°(just before failure) and sensor time displacements as a function of the bucket rotation angle

5. 1 .

 1 Escoubous dam 5.1.1. Fully DEM modelIn 2004, Rémi Deluzarche modeled during his PhD thesis the Escoubous dam, using a fully DEM approach with PFC2D software[START_REF] Deluzarche | Modélisation discrète des enrochements: Application aux barrages[END_REF][START_REF] Deluzarche | Discrete modelling of rock-fill: Application to dams; modelisation discrete des enrochements: Application aux barrages[END_REF]. He created and tested four breakable blocks by assembling rigid elementary particles obtained by overlapping disks of different sizes, presented in Figure I.31 with breakable bonds represented in black. In this panel of possible shapes for the blocks considered in the simulations, the choice finally fell on block n°3, closer to the shape and mechanical behavior of the actual rockfill. The blocks had dimensions between 7 and 20 cm. The normal K n and tangential K s stiffnesses were taken equal to 10 8 N/m in order to meet the assumption of rigid grains (sufficiently small overlaps compared to the dimension of the blocks). The coefficient of friction at contact was fixed at 0.7 by reproducing friction tests between polished surfaces. The fracture resistance of the blocks was defined from crushing tests of individual blocks according to their orientation, then validated on biaxial compression tests by obtaining the same probability of failure in the 2D model as in the real 3D case.

Fig. I. 31 .

 31 Fig. I.31. Definition of four breakable block shapes (type from 1 to 4, from left to right).

Fig. I. 32 .

 32 Fig. I.32. (a) Actual dam; and (b) fully 2D DEM model of the dam with the different zones

Fig. I. 33 .

 33 Fig. I.33. Failure of the dam during its construction, without pitching

  The filling was simulated by gradually increasing the water level of the reservoir and applying normal forces on the blocks forming the upstream face, depending on their position below the water level of the reservoir. The displacement fields are shown in Figure I.34 with non-breakable blocks and Figure I.35 with breakable blocks.

Fig. I. 34 .

 34 Fig. I.34. Displacements after filling the reservoir (unbreakable blocks)

Fig. I. 36 .

 36 Fig. I.36. Horizontal displacement of the crest during filling and emptying of the reservoir

  Figure I.38 and Figure I.38. The behavior and displacements obtained are close to those expected (see section 2.3.2), which allowed the validation of the LK-Enroch model parameters used. Moreover, Chen[4] modeled another rockfill dam (CFRD, concrete face rockfill dams) using LKE model during his PhD thesis in 2012. He concluded that this elasto-plastic model was able to reproduce most of the rockfill characteristics (positive and negative hardening, dilatancy, plastic volume deformation under hydrostatic load, etc.

Fig. I. 37 .

 37 Fig. I.37. Elements and geometry of the Escoubous dam model (©EDF-2009)

  in 1853. The design of DSRWs has recently received new recommendations thanks to new full-scale tests conducted in France by Colas et al. (2010-2013) [46] [47] (Figure I.40a) and Villemus et al. (2006)[48] (Figure I.40b). However, Mundel et al. (2009, 2010)

  [START_REF] Colas | Yield design of dry-stone masonry retaining structures-comparisons with analytical, numerical, and experimental data[END_REF][START_REF] Colas | Yield design modelling of dry joint retaining structures[END_REF]) and other analytical tools (Alejano et al. 2012[START_REF] Lr Alejano | Stability of granite drystone masonry retaining walls: Ii. relevant parameters and analytical and numerical studies of real walls[END_REF][54], Kim et al. 2011[START_REF] Kim | Behavior of dry-stone segmental retaining wall using physical modeling and numerical simulation[END_REF], Terrade et al.

Fig. I. 40 .

 40 Fig. I.40. Experimental loading of dry-stone retaining walls: a) hydrostatic pressure done by Villemus [48]; b) backfill loading done by Colas[START_REF] Colas | Assessing the two-dimensional behaviour of drystone retaining walls by full-scale experiments and yield design simulation[END_REF]

Fig. I. 41 .

 41 Fig. I.41. Numerical models developed by Oetomo [58] for the studied dry-stone retaining wall: a) Full DEM model; b) Mixed FDM-DEM model

Figure II. 1

 1 Figure II.1 presents the four numerical models associated with the four PEDRA cases. They have the same dimensions as in the experiments (2 meters in height, 4.2 meters in length and 2.25 meters wide). Regular flawless parallelepipedic blocks with dimensions that match the average stone dimensions in each case are used to model the dry stone pitching (Figure II.2).

Fig. II. 1 .

 1 Fig. II.1. Numerical models of the four PEDRA cases

Fig. II. 2 .

 2 Fig. II.2. Dimensions of the dry-stone pitching used in the numerical model

Fig. II. 4 .

 4 Fig. II.4. Example of the modeling procedure and steps (displacement fields shown as examples with only qualitative meaning)

Fig. II. 5 .

 5 Fig. II.5. Mesh study of the FLAC3D model

Fig. II. 6 .

 6 Fig. II.6. Computation scheme involving equilibrium and failure criteria

Fig. II. 7 .

 7 Fig. II.7. The rockfill used in PEDRA experimental campaign in addition to its particle size distribution curve

Fig. II. 8 .

 8 Fig. II.8. LEP'S diagram: Effect of the confining pressure on the friction angle

3. 1 . 2 .

 12 Fig. II.9. Views of the two yield surfaces (LK Enroch model): in the meridian (left) and deviatoric (right) planes. S II is the 2nd Invariant of the deviatoric stress tensor

Fig

  Fig. II.10. Experimental triaxial tests and simulations using MC and LKE models: deviatoric stress (right) and volumetric strain (left) at 100kPa and 200kPa

Fig. II. 12 .

 12 Fig. II.12. Slide test involving stones to identify friction angle

Fig. II. 14 .

 14 Fig. II.14. Modeling of the anchorage in case 3

Fig

  Fig. II.15. Failure angles for the four cases using MC and LKE model

  . The same scale was used between simulated and experimental results (crest is at the top of Figure II.17c). LKE model gave a qualitatively better prediction with respect to experiments with a more concentrated displacement field in the middle part of the pitching than MC model. It also holds true quantitatively. On the contrary, MC gave much lower displacements at the last stage before failure which shows its disability to predict the displacement fields for such structures that are characterized by large deformations before failure.For case 4, the evolution of the relative normal displacements along the downstream stone pitching height as a function of the tilting angle was tracked (Figure II.18) to better understand the mechanical behavior of the stone pitching. The relative normal displacements correspond to the difference between the stone pitching normal displacements and the rigid body displacement of the pitching. One can note the global deformation of the downstream pitching increases and the peak displacement point moves up with the increase of the tilting angle. Similar trends are visible in the remaining cases (Figures E.4, E.5, E.6 in Appendix E) but with different peak values and positions. For a smaller tilting angle, if the zone confronting the highest displacement is situated at around 0.3*H (H represents dam height), it stabilizes at approximately 40 % of the height at failure, which is also the location of the backfill failure zone.

Fig. II. 16 .

 16 Fig. II.16. Displacements at the 6 sensors on the downstream face; case 4

Fig. II. 20 .

 20 Fig. II.20. Displacement increments between the last stage before failure in the static safety factor calculation by parameters degradation and the initial static equilibrium, for Case 1 with and without water pressure

Fig. II. 21 .

 21 Fig. II.21. Irregular blocks used with imperfections on the: a) bottom side, b) upper side and c) upper+bottom side of the block

Fig. II. 22 .

 22 Fig. II.22. Numerical model view of the irregular blocks with imperfection on the: a) upper side, b) bottom side and c) upper+bottom side of the block

Figure II. 23 Fig. II. 23 .

 2323 Fig. II.23. Normal displacements of the downstream pitching at rotation of 20 degrees for the 3 studied irregular blocks

Fig. II. 24 .

 24 Fig. II.24. Normal displacements of the downstream pitching at different tilting angles for different contact surface percentages (irregular blocks with bottom imperfection)

Fig. II. 26 .

 26 Fig. II.26. View of the modeled vertical imperfections in the vertical contact planes between the blocks: (a) gap over the entire block height (no contact), (b) gap with partial contact

Fig. III. 1 .

 1 Fig. III.1. FLAC3D-PFC3D coupled numerical model of Escoubous dam and representation of the various elements

Fig. III. 2 .

 2 Fig. III.2. Dimensions of the numerical model and its different elements

Fig. III. 3 .

 3 Fig. III.3. Results of a triaxial compression test for 100 kPa confinement (LKE model) -Deviatoric stress (left) / Volumetric strain (right)

  Fig. III.4. The 4 stages of the dam construction

Fig. III. 5 .

 5 Fig. III.5. Displacement field at the end of the construction phase: a) total displacement; b) horizontal displacement; c) vertical displacement

Figure III. 6

 6 Figure III.6 shows the shear stress field and the main stress directions in the backfill. Shear stresses are the highest in the denser rockfill sub-layer, which is made up of rock placed more tightly than the dam body. It thus concentrates extra load transfer phenomena from the pitching in the direction of the foundation (yellow tainted arrows in Figure III.6 close to the lower downstream pitching).

Fig. III. 6 .

 6 Fig. III.6. Principal stress directions and shear stress field in the dam body after construction, and in particular in the hand-rearranged rockfill sub-layer under the pitching

Fig. III. 7 .

 7 Fig. III.7. Normal forces, shear forces and mobilized friction ratio between the blocks, within the pitching and after construction

Fig. III. 9 .

 9 Fig. III.9. Force system applied to upstream pitching stones at the end of the reservoir water filling

Fig. III. 10 .

 10 Fig. III.10. Horizontal (left) and vertical (right) displacements of the dam after reservoir filling

Fig. III. 12 .

 12 Fig. III.12. Normal forces, shear forces and mobilized friction ratio between the blocks, within the pitching after water reservoir filling

Fig. III. 14 .

 14 Fig. III.14. Intensity and direction of displacements in the body of the dam and the pitching, after construction (left) and after water filling (right)

Fig. III. 16 .

 16 Fig. III.16. Evolution of shear deformation within the dam body during the parameter reduction process in the presence of water pressure

Fig. III. 18 .

 18 Fig. III.18. Evolution of shear deformation within the dam body during the rotation process (top). Shear force and mobilized friction ratio between downstream blocks after 7°rotation (bottom)

Fig. III. 19 . 4 . 2 .

 1942 Fig. III.19. Displacements (left) and shear strain (right) of the 1 st layer of the initial dam model: observed failure

Fig. III. 21 .

 21 Fig. III.21. Horizontal (left) and vertical (right) displacements of the studied model (ba-sic+pitching) after water filling

Fig. III. 22 .

 22 Fig. III.22. Horizontal and vertical displacements (top) and shear strains (bottom) after construction: (a) without berm; (b) with berm

Fig. III. 23 .

 23 Fig. III.23. Horizontal and vertical displacements (top) and shear strains (bottom) after water filling: (a) without berm; (b) with berm

Fig. III. 25 .

 25 Fig. III.25. Horizontal and vertical displacements (top) and shear strains (bottom) after construction: (a) without hand-rearranged rockfill sub-layer under the pitching; (b) with hand-rearranged rockfill sub-layer under the pitching

Fig. III. 26 .

 26 Fig. III.26. Horizontal and vertical displacements after water filling: (a) without hand-rearranged rockfill sub-layer under the pitching; (b) with hand-rearranged rockfill sub-layer under the pitching

Fig. III. 27 .

 27 Fig. III.27. Horizontal and vertical displacements after construction: (a) homogeneous foundation; (b) non-homogeneous foundation

  Fig. III.30. Horizontal and vertical displacements after construction: (a) without upstream concrete mask; (b) with upstream concrete mask

Fig. III. 33 .

 33 Fig. III.33. Post-construction displacements: (a) with downstream pitching; (b) without downstream pitching

Fig. III. 34 .

 34 Fig. III.34. Horizontal and vertical displacements (top) and shear strain (bottom) after water filling: (a) with downstream pitching; (b) without downstream pitching

4. 7 . 3 .

 73 Hand-rearranged rockfill sub-layer: friction of 55°and porosity of 0.3 (pessimistic case)

Fig. III. 35 .

 35 Fig. III.35. Displacements after construction: (a) initial case; (b) modified case (pessimistic case)

Fig. III. 36 .

 36 Fig. III.36. Horizontal and vertical displacements after water filling: (a) initial case; (b) modified case (pessimistic case)

Fig. III. 38 .

 38 Fig. III.38. Displacements at the end of construction and after water filling using: (a) parallelepipedic blocks; (b) trapezoidal blocks (with imperfections)

Fig. III. 39 .

 39 Fig. III.39. Shear strains filed after the water filling phase using: (a) parallelepipedic blocks; (b) trapezoidal blocks (with imperfections)

Fig. III. 41

 41 Fig. III.41. a) Bulges found on the lower part of the downstream pitching around the valve chamber; b) Numerical modeling of the effect of the valve chamber.

Fig. III. 43 .

 43 Fig. III.43. Normal forces (top), shear forces (middle) and mobilized friction (bottom) between pitching blocks at the end of construction (left) and after water filling (right).

Figure III. 45

 45 Figure III.[START_REF] Burgoyne | Revetments of retaining walls[END_REF] shows the normal displacement of the downstream pitching blocks throughout the tilting process. It is noticed that the intensity of the local displacements (forming the bulges) around the rigid structure increases remarkably with the increase of rotation. It reaches a value of almost 15 cm at a rotation of 6°(last step before failure) which is half of the pitching thickness (30 cm).

Fig. III. 44 .

 44 Fig. III.44. Two sections of the pitching after construction and at tilting angle = 6°: 1) at exactly the middle of the model and 2) just 1 m far from the dam lateral boundary

Fig. IV. 1 .

 1 Fig. IV.1. Seismic zones distribution in France

Fig. IV. 2 .

 2 Fig. IV.2. Dynamic boundary conditions of the numerical model

Fig. IV. 3 .Fig. IV. 4 .

 34 Fig. IV.3. Shear wave speed (V s ) and shear modulus profiles along the dam height: numerical model and Gazetas estimation

Fig. IV. 5 .

 5 Fig. IV.5. Theoretical and numerical (Hysteretic) shear modulus degradation (left) and damping (right) curves for the dam rockfill

  Fig. IV.7. Spectra of the five different generated signals with the target spectrum

Five

  Fig. IV.8. Five different signals provided by EDF and with a PGA of 0.16g: acceleration (left), velocity (middle), and displacement (right) versus the dynamic time

Fig. IV. 9 .

 9 Fig. IV.9. Arias intensity of signal 1749

Fig. IV. 11 .

 11 Fig. IV.11. Maximum shear strain in the dam body throughout the dynamic loading using an elastic model without damping

Fig. IV. 13 . 1 .

 131 Fig. IV.13. Target (1749 signal) vs numerically measured outcrop signal spectrum

Fig. IV. 14 .

 14 Fig. IV.14. Horizontal and vertical displacements at dam crest at the end of the horizontal shaking; (a): case A ; (b) case B.

Fig. IV. 17 .

 17 Fig. IV.17. Traced blocks throughout the horizontal shaking; (a) upstream pitching ; (b) downstream pitching.

Fig. IV. 18 .

 18 Fig. IV.18. Horizontal and vertical displacements of the traced pitching blocks throughout the horizontal shaking; (a) case A; (b) case B.

Fig. IV. 19 .

 19 Fig. IV.19. Horizontal acceleration at the dam base and crest throughout the horizontal shaking; (a) case A; (b) case B

Fig. IV. 21 .

 21 Fig. IV.21. Horizontal and vertical displacements at the dam crest throughout the shaking having: (a) a horizontal component (case B, green curves); (b) both a horizontal and vertical component (case C, red curves).

Fig. IV. 22 .

 22 Fig. IV.22. Shear strain of the studied dam at the end of the shaking having: (a) a horizontal component (case B); (b) both a horizontal and vertical component (case C).

Fig. IV. 23 .

 23 Fig. IV.23. Horizontal and vertical displacements of the traced pitching blocks throughout the shaking having: (a): a horizontal component (case B); (b): both a horizontal and vertical component (case C).

  Figure IV.27(a) and (b)), the horizontal evolution curves of all the points in addition to their final residual value shifted toward the downstream side. The vertical displacements increased slightly in most of the studied pitching points.

Fig. IV. 24 .

 24 Fig. IV.24. Horizontal and vertical displacements of the studied dam at the end of the horizontal shaking; (a) without hydrostatic pressure (case C); (b)with hydrostatic pressure (case D)

Fig. IV. 26 .

 26 Fig. IV.26. Shear strain of the studied dam at the end of the horizontal shaking; (a): without hydrostatic pressure (case C); (b): with hydrostatic pressure (case D).

Fig. IV. 27 .

 27 Fig. IV.27. Influence of the hydrostatic pressure on the displacements at monitored points of the pitching throughout the dynamic loading;(a): without hydrostatic pressure (case C); (b): with hydrostatic pressure (case D).

Fig. IV. 28 .

 28 Fig. IV.28. Horizontal acceleration at four points in the downstream pitching throughout the shaking.

Fig. IV. 30 .

 30 Fig. IV.30. Normal and shear forces and mobilized friction between the pitching blocks before and after the dynamic loading

Fig. IV. 31 .

 31 Fig. IV.31. Normal and shear forces and mobilized friction at the interface (between the pitching and the dam body) before and after the dynamic loading

Fig. IV. 32 .

 32 Fig. IV.32. Historical percentage of settlements of different types of embankment dams as a function of the peak ground acceleration

3 (

 3 Figure IV.33) which represents the dynamic safety factor of the dam.

Fig. IV. 33 .

 33 Fig. IV.33. Failure of the dam for a reduction factor equal to 1.3 at 8s and 9s of dynamic loading

Fig. IV. 34 .

 34 Fig. IV.34. Failure of the dam at a scaling factor of 3.3 (PGA = 3.3 x 0.16g) at 12s and 13s of dynamic loading

Fig. A. 2 .

 2 Fig. A.2. Dry stone pitching deformations of case 2 at bucket rotation angle of 20°(just before failure) and sensor time displacements as a function of the bucket rotation angle

Fig. B. 2 .

 2 Fig. B.2. Calculation sequence performed at each time step ∆t in FLAC3D

  Figure B.3), including: a) Determining the time step: the distinct element method (DEM) requires a time step small enough to guarantee model stability and to allow contacts to be created between 2 interacting bodies. In simplified terms, the critical time step depends on the ratio between the smallest mass m and the highest contact stiffness k, raised to the power (1/2) (∆t crit = m k ); b) Law of motion: the position and velocity of each body are updated according to Newton's laws of motion using the calculated time step and the forces/moments calculated in the previous cycle; c) Advance time: the modeled time advances by adding the current timestep; d) Contact detection: contacts are created/deleted dynamically, based on the current position of all model elements (particles, walls);e) Force-discplacement law: forces/moments are recalculated for each contact using the contact model applied to it.

Fig. B. 3 .

 3 Fig. B.3. Sequence of DEM primary operations that occur during each cycle in PFC3D

  Fig. B.4. Behavior and rheological components of the linear model

Fig. B. 5 .

 5 Fig. B.5. Force-displacement law for the linear component of the unbonded linear-based models: (a) normal force versus surface gap, (b) shear force versus relative shear displacement, and (c) slip envelope

  Figure B.6 shows schematically the various components of the Parallel Bond model. This element acts in parallel with the linear contact model (described above) represented by (k n , k s , µ and g s ).

Fig. B. 6 .

 6 Fig. B.6. Behavior and rheological components of the parallel bond model

Figure B. 7

 7 Figure B.7 shows the force-displacement laws applied to each force and moment component at the level of the "Parallel Bond" contact: for each component of the bond, a force/moment level leads to bond failure. These different local behaviors generate a behavior schematized in Figure B.8: we find the different numerical parameters that define the behavior of the "Parallel Bond", namely: tensile strength ( σc ), local cohesion ( c) and angle of friction ( φ ).

Fig. B. 8 .

 8 Fig. B.8. Failure envelope for the parallel bond

Fig. C. 1 .

 1 Fig. C.1. MC failure envelop in the deviatoric plane

p c = p c0 e β ε p ν (C. 7 )σ i j = λ i δ i j 3 (C. 8 )•

 ν738 The flow rule is associated (Equation C.8) in which λ i is the plastic multiplier of the isotropic mechanism.ε˙i i j = λ i ∂ f i ∂ Plastic deviatoric mechanismConcerning the plastic deviatoric mechanism, LKE was developed on the base of Hoek and Brown rule[START_REF] Hoek | Hoek-brown failure criterion-2002 edition[END_REF], HB-CJS rule[START_REF] Cambou | An elastoplastic model for granular material using three yielding mechanisms[END_REF] and LK model. The deviatoric yield surface is mainly derived from the one of the original LK model. It is represented by Equation C.9. It is adapted so that compression is at Lode angle θ = 0.f d (σ i j ) = s II h(θ )σ c h c 0 [As II h(θ ) + BI 1 ] s II = 2J 2 (C.10) h(θ ) = (1 -γcos3θ )1/6 ; h c 0 = h(0); h e 0 = h(π/3) (C.11)

3 e ḋ i j .e ḋ i j ) 1 / 2 dt, where e ḋ i j = ε ḋ i j - ε ḋ v 3 δ

 123 i j and ε ḋ v = ε ḋ kk (C.12)

η

  a(γ d ) = a 0 + ln(1 + γ d X ams γ peak )( a peaka 0 ln(1 + 1/X ams ) ) (C.13) m(γ d ) = m 0 + ln(1 + γ d X ams γ peak )( m peakm 0 ln(1 + 1/X ams ) ) (C.14)Throughout the post-peak phase, the deviatoric yield surface evolves from the peak threshold to the residual threshold, for γ peak <γ d <γ res , according to:a(γ d ) = a peak + ln(1 + 1 η γ dγ peak γ resγ peak )( 1a peak ln(1 + 1/η)) is a model parameter while σ pt corresponds to the abscissa of the intersection point of the peak and residual thresholds, in the plane (σ min ; σ max ). It is determined analytically by:

17

 17

  )

2sinψ 1 -

 1 sinψ (C.20) λ d represents the plastic multiplier of the deviatoric mechanism. More details associated to LKE model can be found in[4] A triaxial test representation using LKE model (stress-strain curve) is given in Figure C.3. One can notice that the curve is split into four different domains. The first one is the elastic domain (0) in which reversible deformations are observed. Then, plastic hardening domain (1) starts with plastic deformations until it reaches the peak deviatoric stress value: until this point, the volumetric behavior is mostly contractive. At peak, the softening domain (2) starts with dilative volumetric behavior until it reaches the residual domain (3) which corresponds to the critical state.

Fig. C. 3 .

 3 Fig. C.3. Deviatoric triaxial behavior using LKE model

Fig. D. 1 .

 1 Fig. D.1. (log(p c ), D r ) plane

Fig. E. 2 .

 2 Fig. E.2. Displacements at the 6 sensors of the downstream face (case 2)

Table I .

 I 1: Rockfill dams built in France from 1940 to 1960

	Date Dam	Height Date Dam	Height
	1942 Araing	25 m	1951 Chammet	19.3 m
	1943 Laurenti	15 m	1951 Gréziolles	30 m
	1949 Grandes Pâtures 20 m	1952 Saint-Nicolas 6 m
	1950 Auchaize	6 m	1953 Escoubous	20 m
	1950 Portillon	22.5 m 1953 Les Laquets	13 m
	1950 Vieilles Forges	10 m	1959 La Sassière	30 m

Table I . 2 :

 I2 Four cases studied in PEDRA experiments

		Case 1	Case 2	Case 3	Case 4
	Stone pitching type	Granite	Ordered	Disordered Schist Granite
		(1 layer)	Schist	with anchorage	(2 layers)
	Surface density	115 kg/m 2	114 kg/m 2 177 kg/m 2	230 kg/m 2
	Stone thickness	5.6 cm	5 cm	5 to 12 cm	11.2 cm
	Stone height	5.6 cm	5 cm	3 to 25 cm	5.6 cm
	Stone width	5 and 10 cm 3 to 40 cm 3 to 40 cm	5 and 10 cm
	Stone friction angle (Slip tests) 29°26.5°26.5°29°T
	hey were placed on the top, middle, and upper parts of it

Table I .

 I

	Wall Label	C2s	C3s	C41
	Experiments	2.3 S/T	2.78 T	2.72 T

4: Critical backfill height obtained from numerical simulations (using mixed DEM-FDM and Full DEM approaches) and full-scale experiments; S: failure by sliding/ T: failure by toppling; The percentage represent the difference from the experiments 2.46 T (7 %) 2.68 T (4 %) 2.53 T (7 %) Full DEM (2 trials) 2.38 T (3 %) 2.72 T (2 %) 2.63 T (3 %) Mixed DEM-FDM 2.41 T (5 %) 2.62 T ( 6 %) 2.52 T (7 %)

Table II

 II

	200 kPa (triaxial tests)	39.5	2	0	52	0.2
	100 kPa (triaxial tests)	42	6	0	50	0.2
	1-10 kPa (extrapolated parameters)	52	18	0	10	0.2

.1: Mohr-Coulomb parameters calibration at different confining pressures Confining Pressure φ peak (°) ψ(°) C E(MPa) ν

Table II

 II

	Table II.3: Calibrated sixteen parameters of LKE model
		General Parameters
	E	50E6 Pa	ν		0.2
	n elas	0.5	σ c		100E6 Pa
	Elasticity threshold parameters
	a 0	0.87	m 0		0.1
		Peak parameters
	a peak at 100/200 kPa	0.79	m peak at 100/200 kPa	0.92
	a peak adapted to 1-10 kPa	0.815	m peak adapted to 1-10 kPa	0.89
		Hardening parameters
	γ˙p eak	0.125	γ˙r es		0.8
	X ams	0.005	µ		2
	Volumetric deformations parameters
	φ res	38°ψ 0		38.5°I
	sotropic mechanism parameters
	p c0 (scaled-down dam model) 15kPa	β		35
	p c0 (triaxial tests)	1MPa
	.2: LKE peak threshold parameters calibration
		p= 200 kPa	p= 100 kPa
		Experimental Numerical Experimental Numerical
	Peak deviatoric stress q (kPa)	680	690	410	400

Table II

 II

	Failure angle	Case 1	Case 2	Case 3	Case 4
	Experiments	24°21°29.5°28°L
	KE	25°(4.2 %) 23°(9.5%) 30°(1.7%) 28°(0%)
	MC	23°(-4.2 %) 22°(4.8%) 30°(1.7%) 28°(0%)

.4: Failure angles for the four cases using MC and LKE model and their percentages of error (shown between parenthesis) with respect to the experimental results

Table II .

 II 5: Safety factors for the reference case with no pitching and for the four studied cases

		No pitching Case 1 Case 2 Case 3 Case 4
	F φ (wow) 1.33	2.22	2.13	3.02	3.17
	F φ (ww) 1.33	2.13	2.05	3.02	3.02
	F θ	1.28	2.7	2.5	3.7	3.3

Table III

 III The various elements of the foundation were modeled using a Mohr-Coulomb constitutive model. The parameters of each part have been deduced from previous studies on such soil types[START_REF] Dayre | Mechanical properties of granodiorite from laboratory tests[END_REF][START_REF] Lr Alejano | Considerations of the dilatancy angle in rocks and rock masses[END_REF][START_REF] Hubbard | Impact of a rock avalanche on a moraine-dammed proglacial lake: Laguna safuna alta, cordillera blanca, peru. Earth Surface Processes and Landforms[END_REF][START_REF] Laloui | Experimental and numerical investigations of the behaviour of a heat exchanger pile[END_REF][START_REF] Bourdeau | Le comportement des alluvions du rhône dans une grande boîte de cisaillement direct[END_REF] and are presented in TableIII.3.

	3.3.3. Foundation
	Table III.3: MC parameters for the different parts of the foundation
	Density (kg/m 3 ) E (MPa)	ν	φ (°) ψ (°) C (kPa)
	Granodiorite	2,700	72,000 0.30 43	18	74,000
	Moraine	1,750	90	0.30 38	0	5
	Alluvions	1,750	190	0.22 35	0	9
	Concrete	2,500	40,000 0.30 36	10	3,900
	3.3.4. Dry-stone pitching
	.2: MC parameters for the re-arranged rockfill layer
		E	ν	φ	ψ C
	75 x 11 = 825 MPa 0.25 58°18°0

Table III .

 III 4: Summary of previously studied cases: horizontal (x) and vertical (z) displacements at crest and mid-height (berm) in addition to safety factors

		Crest displacement, construction after	Crest after filling displacement,	Displacement after filling above berm	Rotational angle (°) failure	F θ	F φ
			Horiz. Verti. Horiz. Verti.
	Initial model:
	backfill + homogeneous	-	-	-	-	-	-	<1	<1
	foundation
	Initial model + pitching	24	3	9	2	7	8	1.33 1.4
	Previous model + berm	38	4	12	4	8	6	1.23 1.3
	Previous model + rockfill sub-layer	24	4	9	4	7	8	1.33 1.4
	Previous model
	+ non-homogeneous	28	6.5	11.5	5.5	10	8	1.33 1.4
	foundation
	Previous model
	+ concrete mask	32	6	10	5	8	8	1.33 1.4
	(Complete model)

  LEP'S diagram with a cohesion equal to zero. The dilatancy angle was taken equal to zero due to the loose state of the dam body rockfill. The MC parameters of all the model elements are summarized in TableIV.2. Table IV.2: MC parameters for the different parts of the system for dynamic simulations (parameters different from the static model in bold)

		G (MPa) K (MPa) φ (°) ψ (°) C (kPa)
	Dam body (at 100 kPa)	165	275	42.5	0	0
	Sub-layer	330	550	58	18	0
	Moraine	300	650	38	0	5
	Alluvions	197	286	35	0	9
	Granodiorite	72000	0.30	43	18	74000

Mots clés: Enrochement, pierres, structure hydraulique, sismique, FDM, DEM, discret-continu

The area just below the pitching (gray) is generated by assembling samples with more complex shapes. Each domain was established individually using the same process as for the central zone. All of these domains have the same void index and are made up of blocks in equilibrium.

It is defined as the ratio between the average of unbalanced mechanical forces for all nodes in the model (or at contacts for rigid bodies modeled in PFC3D) and the average of mechanical forces applied to all nodes (or to an element in PFC3D).

Shear strain is computed from the reference state which corresponds to the last initialization of displacements. It was performed at the end of the construction phase, once the system was found at equilibrium under gravity

Acknowledgements

I want to say a big thank you to all the people and organizations who played a vital role in helping me complete my Ph.D.

Displacements and shear strains after water filling are also compared (Figure III.31): horizontal and vertical displacements after water filling are slightly lower with the addition of the concrete mask. This is due to the greater stiffness of the upstream face of the pitching, which reduces the concentration of deformations induced by water pressure. Shear deformations are almost the same with and without the concrete mask.

The dam safety factors were calculated for the case with a concrete mask, but no change is observed. The static safety factor remains unchanged at 1.4 (with and without water pressure). Similarly, the angle of failure by rotation is still 8°(seismic safety factor of 1.33), as in the case without mask.

From this study, we can conclude that the upstream concrete mask has little influence on the mechanical behavior of the dam. Its main role is to create an impermeable layer on the upstream face, blocking the flow of water into the dam body. In this study, the complete dam was modeled but without the upstream and downstream dry-stone pitching in order to assess the contribution of the dry-stone pitching to the stability of the dam.

In this case, without the pitching, failure occurred at the last stage of construction: 4 th stage out of 4 (Figure III.32). A failure surface developed from the bottom of the upstream face towards the crest. This failure is provoked by the loss of the added resistance due to the pitching stones. Failure develops more rapidly on the upstream side, due to the lower thickness of the hand-placed rockfill layer: 20 cm versus 90 cm on the downstream side.

In conclusion, the dam is not stable without the dry-stone pitching, and a safety factor of less than 1.0 is thus obtained, even in the presence of the denser surface hand-placed sub-layer which once again justifies the important stabilizing role of the dry-stone pitching. c) Dynamic-solution approach: in Newton's equation of motion, the inertial terms are used as numerical means to reach the equilibrium state; the advantage of this dynamic procedure is that it allows the numerical scheme to be stable even if the physical system is unstable.

Appendix C Constitutive modeling

Two constitutive models associated to the dam body are used. Their descriptions are presented in the upcoming subsections.

Mohr-Coulomb

Mohr-Coulomb denoted MC is an elastic perfectly plastic constitutive model which is widely used in geotechnical computations for its simplicity. The associated failure surface is a straight cone. Its envelope in the deviatoric plane is given in The yield surface F which is also the failure surface reads:

where σ 1 and σ 3 are the major and minor principle stresses respectively.

The plastic potential G reads:

The flow rule is non-associated for granular soils which means friction angle φ and dilatancy angle ψ are different. Cohesion is equal to zero for purely frictional soils.

Although MC is widely used in soil modeling, it has several limits. For example, it does not take into account the effect of the pressure on the Young's modulus, friction angle and dilatancy angle. Moreover, the elastic domain is very large where no irreversible strains can be generated which do not reflect phenomena at stake in actual granular soils. It leads to the underestimation of deformations in the modeled systems. More sophisticated constitutive models may be required in • e max and e min are the maximum and minimum void ratios respectively of the studied granular material. They can be obtained from previous studies [START_REF] Biarez | Influence de la granulométrie et de son évolution par ruptures de grains sur le comportement mécanique de matériaux granulaires[END_REF].

•