The environmental conditions may play a relevant role in the stability of rock slopes. In fact, weathering can contribute either to the reduction of strength of the material or to the increasing of internal stresses in the rock massif. This Ph.D. thesis deals with the effect of atmospheric thermal cycles in rocks, from both experimental and numerical point of view.

This research is focused on the rock obtained from the cliff from La Roque Gageac, a small town located in the south-west of France, placed above a rock cliff that has a cavern at middle high. This village has been affected by several rock falls, which, based on instrumentation data, are linked to the in situ thermal variations.

With the aim of isolating the effect of thermal cycles in the La Roque Gageac limestone, an experimental study is performed. Samples obtained from two different locations in the cliff are considered (blocks felt inside the cavern and samples drilled from the face of the cliff). These samples are submitted to thermal cycles between 10ºC and 50ºC, with the aim of imitating the range of temperatures registered on the site.

In first place, I would like to express my gratitude to my supervisors for their guidance, ideas, contributions and implication in the development of this thesis. To Jean Vaunat, for his help with the implementations and his academic but also personal support. To Muriel Gasc for giving me the opportunity to develop this project and to Jose Darrozes for his help specially in the geological and mineralogical characterization of the material.

I would like to extend my gratitude to Veronique Merrien-Soukatchoff, Guilleme Remien and Maria to have accepted to be the reviewers of this work. This thesis involves the effort and work of different laboratories. I would like to thank the rock mechanical laboratory staff of the Cerema Toulouse, for their help and availability that allowed me to perform most of the experimental program. To the professors and staff from the laboratory GET in the Université Paul Sabatier and Staff from the Cerema Aix en provence and IFSTTAR for their help with the additional tests for the mineralogical composition analysis and porosimetry evaluation.

Moreover, I would like to thank the administrative personal of the Universitat Politècnica de Cataluña and Université Paul Sabatier for being so helpful to carry out with all the administrative procedures.

This research has been developed thanks to the economic support of the Centre d'études et d'expertise sur les risques, l'environnement, la mobilité et l'aménagement (Cerema) Je souhaite aussi remercier tous mes collègues du Cerema Toulouse, notamment la famille REGG ainsi que Gerard, Sébastien, Pierre, Julian et Florence pour avoir été si accueillants et pour m'avoir fait me sentir toujours chez moi, ainsi que pour votre aide constante avec mon français et m'avoir fait connaître tous les délices de la cuisine française. Merci particulièrement à Stephan, Philipe et Sébastien pour avoir toujours été prêts à m'aider, à Virginie pour son travail dans le laboratoire, à Jean Claude pour tous les « petites services » sans lesquels cette thèse n'aurait pas pu être menée à bien et à Didier pour tous ses conseils et son aide, pour avoir été le directeur dans l'ombre Quisiera agradecer a mis amigos, los que están cerca y los que no, por su apoyo y esfuerzo por hacerme sentir siempre cerca y en casa, a Mauro por siempre estar dispuesto a darme una mano y a Aurore, pour toutes les conversations autour de un café et un désert.

Por último, quisiera agradecer a mi familia por su apoyo incondicional, especialmente a mis padres y hermanas que siempre han sabido estar ahí con la palabra indicada para motivarme a seguir adelante, gracias por todo su cariño y apoyo.

The damage induced in the samples is evaluated through measurements of strains, elastic wave propagation velocities and uniaxial compressive strength. It is observed that samples experience a reduction in the elastic wave propagation velocity and material strength through the imposition of thermal cycles. In the same way, a strains accumulation is also registered. Moreover, the response obtained depends on the mineralogical composition of the rock, which is linked to the extraction point of the sample, as the cliff presents a large heterogeneity.

Based on the results obtained in the experimental program the principal characteristics of the thermal damage phenomenon are defined, proposing a constitutive model capable to reproduce the effect of thermal cycles in the La Roque Gageac limestone from a macroscopical point of view that is the effect in the mechanical response of the rock.

With this purpose, a modification on the composed material constitutive model defined by [START_REF] Vaunat | Bond degradation and irreversible strains in soft argillaceous rock[END_REF] is performed, in order to include the characteristics of the material under study. This model considers the material as composed by two different components with an established behavior for each; the response of the composite material will be defined by the interaction between those components. This constitutive model is implemented in the finite elements code CODE_BRIGHT and evaluated with the simulation of the response of the samples observed during the experimental program.

Résumé

Les conditions environnementales jouent un rôle important dans la stabilité des massifs rocheux. De fait, les variations climatiques peuvent affecter la résistance du matériau et également augmenter les contraintes internes dans le massif. Cette thèse étudie les effets des cycles thermiques atmosphériques sur les roches, dans une approche expérimentale puis numérique.

Cette recherche est centrée sur le cas de la falaise de La Roque Gageac, un village du sudouest de la France, situé au pied d'une falaise qui présente une caverne à mi-hauteur. Cette commune a été affectée par plusieurs éboulements. Les données issues de l'instrumentation du massif montrent que ces instabilités sont liées aux cycles thermiques.

Afin d'isoler l'effet des cycles thermiques dans la roche calcaire de La Roque Gageac, un programme expérimental a été mis en place. Des échantillons ont été prélevés sur deux sites de la falaise (blocs éboulés à l'intérieur de la caverne et carottés dans la face de la falaise). Ces échantillons ont été soumis à des cycles thermiques entre 10ºC et 50ºC afin de simuler les conditions enregistrées sur le site.

L'endommagement de la roche est évalué en laboratoire par le suivi de mesures de déformations, de vitesses de propagation d'ondes élastiques et par l'évolution de la résistance à la compression uniaxial. Une réduction de la vitesse de propagation des ondes élastiques de compression et de cisaillement ainsi que de la résistance à la compression du matériau sont observées. De la même façon, les échantillons enregistrent une accumulation de déformations. De plus, la réponse des échantillons est influencée par la composition minéralogique de la roche. Celle-ci est liée au lieu de prélèvement, le calcaire de la falaise montrant une importante hétérogénéité.

Les résultats du programme expérimental ont permis la détermination des caractéristiques principales du phénomène d'endommagement thermique. Une loi de comportement modélisant l'effet des cycles thermiques sur la roche calcaire de La Roque Gageac d'un point de vue macroscopique a été proposée.

À cette fin, une modification de la loi de comportement pour matériau composite présentée par [START_REF] Vaunat | Bond degradation and irreversible strains in soft argillaceous rock[END_REF] est suggérée, afin d'inclure les caractéristiques propres de la roche étudiée. Ce modèle prend en compte un matériau formé de deux constituants au fonctionnement distinct dont l'interaction décrit le comportement de la roche. Le modèle a été introduit dans le programme aux éléments finis « CODE_BRIGHT » et est évalué en simulant la réponse observée sur les échantillons pendant le programme expérimental.

Resumen

Las condiciones ambientales juegan un papel relevante en la estabilidad de taludes de roca, de hecho, las variaciones climáticas pueden afectar tanto la resistencia del material como el nivel de esfuerzos internos a los que el macizo está sometido. Esta tesis doctoral se enfoca en el efecto de los ciclos térmicos atmosféricos en las rocas, desde un punto de vista experimental y numérico.

Esta investigación se centra en el material del acantilado de La Roque Gageac, un pequeño pueblo localizado en el Sur de Francia ubicado a los pies de un acantilado de roca, que presenta una caverna a mitad de altura. Esta población se ha visto afectado por importantes caídas de rocas, las cuales, de acuerdo con registros de instrumentación parecen estar ligadas a las variaciones de temperatura registradas en el sitio.

Con el fin de aislar el efecto de los ciclos térmicos en la caliza de la Roque Gageac un estudio experimental es llevado a cabo, para esto se consideran muestras obtenidas de dos sitios diferentes del macizo (bloques caídos al interior de la caverna y muestras extraídas de la cara del macizo). Dichas muestras son sometidas a ciclos de temperatura entre 10ºC y 50ºC con el fin de mantener los rangos de temperatura registrados en el sitio.

El daño inducido en las muestras es evaluado a través de medidas de deformación, velocidad de propagación de ondas elásticas y resistencia a la compresión uniaxial. Se observa una reducción en la velocidad de propagación de onda y la resistencia del material con la imposición de los ciclos térmicos, así mismo las muestras registran una acumulación de deformación, adicionalmente, la respuesta obtenida depende de la composición mineralógica de la roca que a su vez está ligada al lugar de extracción de la muestra, ya que el acantilado presenta una importante heterogeneidad.

Con base en los resultados obtenidos en el programa experimental se definen las principales características del fenómeno de daño térmico, con lo cual es posible proponer un modelo numérico capaz de reproducir el efecto de los ciclos térmicos en la caliza de la Roque Gageac desde un punto de vista macroscópico, es decir, en el comportamiento mecánico de dicha roca.

Para esto se realiza una modificación del modelo constitutivo de material compuesto definido por [START_REF] Vaunat | Bond degradation and irreversible strains in soft argillaceous rock[END_REF], con el fin de incluir las características propias del material bajo estudio. Dicho modelo considera que el material bajo estudio está compuesto por dos componentes diferentes con comportamientos establecidos para cada uno y la respuesta del material compuesto estará determinada por la interacción entre estos componentes. Este modelo constitutivo es implementado en el código de elementos finitos CODE_BRIGHT, y evaluado a través de la simulación de la respuesta observada en las muestras durante el programa experimental. Table 3- Rock slopes instabilities can be triggered by different mechanisms, these instabilities are controlled by several factors that regulate the moment and size of the rockfall. Those factors can be divided into three types:

List of figures

List of tables

1. Predispositions, these factors represent all the characteristics of the rock slope that may help to the instability 2. Preparation, it refers to all those phenomena that favor the reduction of the rock strength, like is the case of weathering. 3. Triggering factors, those that induce the failure.

The predisposition and triggering factors presents an evolution in the time, while the predisposition factors are slope characteristics thus they don't change.

The landslide-sequencing model presented by [START_REF] Julian | Aspects of landslide activity in the Mercantour Massif and the French Riviera, southeastern France[END_REF] shows a schematic model (Figure 1-1) of a slope time evolution, represented by the evolution of the material strength, that reduces with time as a response to the preparation factors and the evolution of the disturbing forces. The failure of a slope occurs when the disturbing forces (triggering factors) are more efficient than those of slope failure resistance, points 1 to 5 in Figure 1-1 correspond to those instants where the resistance and disturbing forces are equal. These points indicate the characteristics of different triggering factors for landslides or rockfalls. Points 5 and 4 represent brief and violent phenomena that can cause the failure, like earthquakes (point 5) that highly increase the disturbing forces, or heavy rainfalls (point 4) that momentarily decrease the resistance of the material.

Furthermore, points 1 to 3 indicate conditions developed at medium to long term where the strength of the material is affected, for instance, the increasing of the slope angle due to the rock erosion that induce some stress increment (point 3); strength and stresses oscillation linked to seasonal variations (point 2) or the gradually reduction in resistances due to geochemical rock weathering or progressive damage of the rock (point 1). These three conditions are related to rockfalls that lack a recognized triggering factor [START_REF] Julian | Aspects of landslide activity in the Mercantour Massif and the French Riviera, southeastern France[END_REF][START_REF] Gunzburger | Influence of daily surface temperature fluctuations on rock slope stability: Case study of the Rochers de Valabres slope (France)[END_REF].

Rock slopes instabilities linked to the progressive degradation of the material due to atmospheric thermal cycles can be represented by point 1 in Figure 1-1. In this type of failure, the predisposition factors are also the triggering factors. Each thermal cycle induces some internal stresses that can lead to the damage of the material, even if, the effect observed for one cycle may not be significant, the accumulated consequences of several cycles can be important, an may reach a point where the failure is induced without any identifiable triggering factor [START_REF] Viles | Linking weathering and rock slope instability: Non-linear perspectives[END_REF].

Typically, the investigations related to rockfalls triggered by the effect of atmospheric thermal variations on the rock mass are mainly centered in the evaluation of freeze-thaw phenomena, typical of cold and humid areas. This phenomenon is related to the generation or propagation of fissures due to stresses generated by the freezing of water present in the pores or fractures; in this process the moisture of the rock and the duration of low temperatures are very important factors, as they control the depth where the freeze-thawing will have effect, therefore the size of the possible blocks generated. [START_REF] Gunzburger | Influence of daily surface temperature fluctuations on rock slope stability: Case study of the Rochers de Valabres slope (France)[END_REF][START_REF] Mufundirwa | A new practical method for prediction of geomechanical failure-time[END_REF][START_REF] Bakun-Mazor | Measuring thermally-induced rock block displacement inside a controlled climate laboratory[END_REF][START_REF] Hall | Thermal fatigue and thermal shock in bedrock: An attempt to unravel the geomorphic processes and products[END_REF]Hall 2004b;[START_REF] Matsuoka | Rockfall activity from an alpine cliff during thawing periods[END_REF][START_REF] Matsuoka | Microgelivation versus macrogelivation: Towards bridging the gap between laboratory and field frost weathering[END_REF]Hall 2004a).

Nevertheless, it has been observed that daily thermal variations with temperature ranges far from extreme values can induce some important thermal stresses, related to the unequal expansion and contraction of the minerals composing the rock. The cyclic nature of atmospheric thermal variations induces a degradation in the rock mass that eventually can lead to failure. [START_REF] Collins | Rockfall triggering by cyclic thermal stressing of exfoliation fractures[END_REF][START_REF] Siegesmund | Physical weathering of marbles caused by anisotropic thermal expansion[END_REF] Iñigo & Vicente-Tavera (2002) found that for different granitic rocks, the thermal variations observed under natural conditions are particularly important for the first 10 cm. Nevertheless, according to [START_REF] Collins | Rockfall triggering by cyclic thermal stressing of exfoliation fractures[END_REF], the effect of the thermal stresses associates to atmospheric temperature variations could be important for depths up to 1m.

One of the principal effects of these thermal stresses is the propagation of existing fissures, as registered by [START_REF] Bakun-Mazor | Measuring thermally-induced rock block displacement inside a controlled climate laboratory[END_REF] and [START_REF] Bakun-Mazor | Thermally vs. seismically induced block displacements in Masada rock slopes[END_REF] that used an instrumentation system to evaluate the response of a group of fissures observed in the dolomite of the Masada rock slope, showing that the seasonal atmospheric variations are enough to induce permanent plastic displacement.

In the same way [START_REF] Collins | Rockfall triggering by cyclic thermal stressing of exfoliation fractures[END_REF] survey the evolution of an exfoliation sheet in the granite of the Yosemite national park in California, with an instrumentation system. Showing that the sheet deforms synchronously with the temperature even when the light intensity is minimal. Additionally, it is observed that humidity has almost no, thus, the temperature is the dominant factor in the deformation on the rock sheets observed in the site.

Furthermore, some studies have shown that these thermal stresses may be responsible of large volume rockfalls, due to the propagation of macro-fissures and micro-fissures that reduces the strength of the rock mass [START_REF] Hall | Thermal fatigue and thermal shock in bedrock: An attempt to unravel the geomorphic processes and products[END_REF][START_REF] Hall | The role of thermal stress fatigue in the breakdown of rock in cold regions[END_REF].

For instance, [START_REF] Gunzburger | Influence of daily surface temperature fluctuations on rock slope stability: Case study of the Rochers de Valabres slope (France)[END_REF] and [START_REF] Merrien-Soukatchoff | Thermal effects on rock slopes: Case study of the "rochers de Valabres" slope (France)[END_REF] present the case of an important rockfall (almost 2000 m 3 ) observed in the southern Alpes region in France, in 2000. In absence of an evident triggering event, the study focuses on the hypothesis that this rockfall is caused by the repetition of small events, like the gradual weathering of fractures. With this purpose, an instrumentation system was installed in the site, focused on the measurement of displacements and climatological conditions. Based on these measurements and with the help of a simple numerical model it is proved that a potential triggering factor of such rockfall is the atmospheric thermal cycles. It is observed that these cycles generate some plastic strains on the internal fissures of the rock mass helping to the instability of the rock slope. [START_REF] Vargas | On mechanisms for failures of some rock slopes in Rio de Janeiro, Brazil: thermal fatigue?[END_REF] and Do Amaral [START_REF] Vargas | On the effect of thermally induced stresses in failures of some rock slopes in Rio de Janeiro, Brazil[END_REF] evaluate a number of failures in rock slopes in the region of Rio de Janeiro (Brazil). Most of those instabilities occurred during the winter (dry season), in absence of important rains or other evident triggering factor, it is evaluated the possible relationship between daily thermal fluctuations and those events, as during the winter the maximal range of daily temperatures is also registered. Based on the analysis of the characteristics of all the rockfalls registered, an experimental study is proposed in order to evaluate the thermal evolution of a block with a fissure. The experimental results are complemented with a numerical finite element code analysis, showing that thermal stresses induced by the atmospheric typical fluctuations are able to propagate the existing non-persistence fractures observed in the rock slopes, generating rockfalls like those registered in the zone.

Similar results are obtained by [START_REF] Vlcko | Rock displacement and thermal expansion study at historic heritage sites in Slovakia[END_REF] in the analysis of the stability of a rock cliff in Slovakia that jeopardize the ruins of an important medieval castle. The stability condition of this cliff is evaluated with the help of a monitoring system, it has been observed that temperature variations can be the triggering factor of the slope movements.

All the studies presented above are focused in the evaluation of temperature as the main climatological condition that affects the rock strength. The slopes under study are located in dry regions and/or the weathering effects are evaluated during dry seasons. Nevertheless, it is important to mention that when this type of analysis is performed in humid regions, where the total humidity of the air may play an important role, it is observed that the cracks preferred orientations is very similar to those registered in dry regions. However, it can be expected that in deserts regions the diurnal temperatures are larger than in more humid regions, thus higher thermal stresses may be induced, affecting more the rock strength. Still, in both cases, fracture propagation can be expected associated to diurnal temperature changes. Therefore the thermal damage cannot be only related to dry climates. [START_REF] Aldred | The influence of solar-induced thermal stresses on the mechanical weathering of rocks in humid mid-latitudes[END_REF][START_REF] Eppes | Cracks in desert pavement rocks: Further insights into mechanical weathering by directional insolation[END_REF].

La Roque Gageac is a small town located in the south-west of France, this town is placed above a cliff mainly composed of limestone, that had suffered important rockfalls. the last event was registered in 2010. The main characteristic of the instabilities registered in the zone is that they can't be attributed to any evident triggering factor, like earthquakes of rainfalls. Therefore, after 2010 event a monitoring system was installed on the site. These measurements have shown that deformation registered in the cliff has an important relationship with the atmospheric temperature [START_REF] Gasc-Babier | Termal fatigue in rocks-La Roque-Gageac's case study[END_REF].

As can be observed from the cases presented before, the thermal damage under atmospheric thermal cycles can be the triggering factor for important rockfalls. Moreover, this phenomenon can be developed in a broad type of rocks and climatological conditions, for this reason, is significant to develop a specific study focus on the main characteristics of it.

Aims of the research and thesis layout.

The overall objective of this research is to evaluate experimentally the phenomenon of thermal fatigue in rocks and propose a constitutive law able to reproduce the effect of thermal cycles on the rock.

To achieve this purpose the research is focused on the case of La Roque Gageac, mentioned in the former section. Instrumentation measurements have shown that there is a strong relationship between the temperature and the cliff deformations. For this reason, an experimental program is proposed and performed with the aim to isolate the effect of atmospheric thermal cycles in the limestone from La Roque Gageac, and evaluate the possible thermal damage induced in this rock. This laboratory study considers samples obtained from different locations on the cliff to cover the heterogeneity observed in the site, the response of those samples to thermal cycles is evaluated with different type of measurements.

Based on the experimental program observations the phenomenological response and the framework of the numerical modeling are defined. One of the principal objectives of this thesis is to propose a constitutive law capable to reproduce the behavior observed at macroscopical level, this is to say the effect of thermal cycles on the mechanical response of La Roque Gageac limestone, with this purpose, a modification of the composite material model proposed by [START_REF] Vaunat | Bond degradation and irreversible strains in soft argillaceous rock[END_REF] is presented and implemented in a finite element method code CODE_BRIGHT. This model is used to reproduce and evaluate the experimental program results.

This thesis is divided into five chapters and two appendices. Each chapter presents its own introduction and concluding remarks. A common list of cited references is presented after chapter 5. From the five chapters, three main chapters are distinguished, the main content of each chapter are as follows:

Chapter 2, presents a description of La Roque gageac case, its geology and climatological characteristics, as a description of the most important events registered on the site, as well as, the description and results of the instrumentation system installed to evaluate the response of the cliff under consideration.

Chapter 3, deals with the experimental program proposed and developed in order to evaluate the effect of atmospheric thermal cycles on the limestone from LaRG, in first place a complete characterization of the rock under study is presented, follow by the description of the complete experimental program performed, finally the results obtained are presented and analyzed, showing that thermal cycles may induce a damage on the limestone from La Roque Gageac.

In chapter 4 the numerical model proposed is presented, the macroscopic phenomena under considerations is detailed, followed by the description of the constitutive law proposed, then, the principal equations controlling the numerical model and the implementation in a finite element scheme is presented. Furthermore, the constitutive model is analyzed through an analytical solution and the general conditions on performance are analyzed. Finally, the numerical model proposed is employed to reproduce the results obtained with the experimental program.

Finally, chapter 5 gives a summary of the principal observations obtained in this thesis, the concluding remarks, and future work.

Additionally, an appendixes section is presented, where the three appendixes mentioned are assembled, the two first appendixes deal with the experimental program described in chapter 3, they are the results of some additional experimental tests performed. Appendix A presents the results of mercury intrusion porosimetry performed on samples submitted to a different number of thermal cycles, and Appendix B gives the results obtained with digital image correlation analysis for two cubical samples (sample C1 and C3) were the objective was to evaluate the generation of new cracks due to thermal cycling. Finally, Appendix C presents the numerical implementation of the model formulation described in chapter 4.

Introduction.

Rock slope instabilities may take several forms, from small rock falls involving blocks of small volume, through rock slides and block glides, to rock avalanches and rotational slumps. The weathering degree may be linked to the rock strength delimiting the size of the rockfalls, moreover, this weathering can also contributes to stress acting in the rock slope being in some cases the triggering factor of those movements [START_REF] Julian | Aspects of landslide activity in the Mercantour Massif and the French Riviera, southeastern France[END_REF][START_REF] Viles | Linking weathering and rock slope instability: Non-linear perspectives[END_REF].

Repeated weathering oscillations can generate rock breakdown enhancing the likelihood of mass movements. For instance, the atmospheric thermal cycling stresses can induce a strength reduction that help to the fracture propagation and eventually may generates the detachment of rock blocks as registered by [START_REF] Gunzburger | Influence of daily surface temperature fluctuations on rock slope stability: Case study of the Rochers de Valabres slope (France)[END_REF][START_REF] Vargas | On mechanisms for failures of some rock slopes in Rio de Janeiro, Brazil: thermal fatigue?[END_REF]. La Roque Gageac is a small town located in the south-west of France that has been submitted to several rockfalls events. These events cannot be related to evident triggering factors, like rain seasons or earthquakes. Therefore, other possible starting events or causes were evaluated, and natural thermal cycles seem to be the most related factor.

In this chapter, a complete description of La Roque Gageac case is presented, including general aspects like localization, geology and climate, as well as, the characteristics of the principal events reported, the stability conditions after the last important rockfall event and the instrumentation system installed. This information is extracted from the technical report presented by the Centre d´Etudes Techniques de l´Equipement du Sud-Ouest in 2010 (Virely andGuittard, 2010a, 2010b;[START_REF] Virely | Project DOSM Falaises de La Roque Gageac Évaluation de l ' aléa chutes de masses rocheuses[END_REF], the master thesis of D. F. [START_REF] Ruiz | Thermo-mechanical analysis of the stability of a rock-cliff under climatic actions[END_REF] and a couple of technical visits realized to the site.

Site location.

La Roque Gageac is a small town located in the south-west of France. It is placed in the northwest of the Aquitaine basin in the department of La Dordogne (Figure 2-1). Although, this village has a small population, more or less 450 people, the tourism is the principal economic activity, with more than 2 million visitors each year in the region.

Geology.

The Périgord Noir landscape has developed along the Dordogne river curse deep valleys. These valleys are not present all the river curse. In fact, the geomorphology depends on the geological material that the river cross.

This fluvial system goes through different types of formations like: massive limestone of the middle and late Jurassic where the landscape is tight (deep valleys) with small meanders. And zones of argillite rocks and marbles, where the river shows large and important meanders with gentle slopes.

When the river reaches the upper cretaceous sandy limestone, where our zone of interest is located, (Figure 2-3) the morphology presents a succession of large rectilinear zones, and tight sections with meanders. These morphological differences are related to the geological heterogeneity in the Aquitaine basin. Which is composed by soft limestones with a high content of sand, zones where the river shows a rectilinear course, and, by strong sandy limestones favorable to meanders formation. The Aquitaine basin, is composed by coniacian limestones. This geological formation can be divided in two types of materials, those belonging to the lower coniacian, that are marine impermeable limestones, and those belonging to the middle and upper coniacian, a layer of sandy limestones, porous, weak and prone to karstification phenomena. La Roque Gageac is placed over the limestone of the lower coniacian, however, the cliff is composed by an upper coniacian rock.

In fact, where La Roque Gageac is located, the Dordogne river has eroded the sandy limestone of the upper and middle coniacian, to reach the lower coniacian layer. This geological process has generated the cliff and induced a modification in the stress balance of the massif, producing a phenomenon of decompression in the cliff.

The geological formations of the middle and upper coniacian shows a succession of marine and lacustrine deposits. These deposits are characterized by oblique lithologies with different inclinations delimited in some cases by erosion superficies and showing also different types and size of material. Figure 2-4 shows a contact between different types of deposits, where the green line delimits an oblique layer and the red line an erosional surface. This geological history is responsible of the strong material heterogeneity and the poor mechanical properties of this rock. Additionally, this material is prone to karstification processes. Due to the heterogeneity of this formation, the cliff is composed by a succession of permeable and impermeable layers. Therefore, the karstification networks cannot be completely developed in vertical sense. In fact, the troglodyte fort observed in the cliff of La Roque Gageac belongs to an inactive paleo-karstification network, that was not completely developed.

Because the calcareous nature of the rock from the cliff under analysis the water might play an important role in the behavior of this material. Nevertheless, the dissolution rates registered in the region are inferior to 1 mm per year [START_REF] Trauth | Géodynamique des altérations ferralitiques sur roches sédimentaires, en bordure sud-ouest crétacée du Massif Central: paysages sidérolitiques en Quercy Blanc, Haut-Agenais, Bouriane et Périgord Noir[END_REF][START_REF] Platel | Le crétacé supérieur de la plate-forme septentrionale du basin d'Aquitaine Stratigraphie et évolution géodynamique[END_REF][START_REF] Platel | Stratigraphie, sédimentologiie et evolution géodynamique de la plateforme carbonatée di Crétacé supérieur du nord du bassin d'Aquitaine[END_REF][START_REF] Platel | Carte géologique de la France à 1/50 000, feuille Terrasson (784), terrains crétacés. Orléans : BRGM[END_REF][START_REF] Astruc | Le paléokarst quecynois au paléogéne altérations et sédimentations associées -Documents du BRGM N o 133[END_REF]Virely and Guittard, 2010a).

Climate.

The climate of the department of La Dordogne is temperate with an important oceanic influence. During the winter, the weather is influenced by the mountains due to the proximity to the Massif Central, and in the summer, it suffers the presence of warm air masses coming from the Mediterranean Sea.

At regional scale, the annual rain is around 860 mm, with an historical maximum of 100 mm in one day, registered in May 1971. The maximal temperatures range from 27ºC to 32ºC and the minimal between -10ºC and -15ºC.

This large temperature variations may be the responsible of the flaking observed in the face of the cliff, especially if we consider that the cliff is oriented to the south, and the high albedo value that this rock presents (Virely and Guittard, 2010a;[START_REF] Gasc-Babier | Termal fatigue in rocks-La Roque-Gageac's case study[END_REF].

La Roque Gageac rockfalls.

Since the beginning of the 20th century, at least four important events have been registered in this town, the latest occurred in the winter of 2010, when a part of the roof of the troglodyte cavern collapsed (Virely andGuittard, 2010a, 2010b;[START_REF] Gasc-Babier | Termal fatigue in rocks-La Roque-Gageac's case study[END_REF]. Since, there is no external load over the cliff (no buildings or any type of constructions) and no hydrological changes all this rockfalls have been attributed to the geological conditions of the rock and the possible effects of climate over this material. After 2010's event, we may consider 3 levels of risk remaining in La Roque Gageac's cliff, with different associated volume:

• Fall of small blocks detached from the cliff, this level of risk also includes the material felt from the roof of the cavern and remaining in it, after the latest important event.

• The residual rock beam inside the cavern that presents a poor stability condition.

• The general stability above the cavern.

The event of 2010 left a high volume of small blocks remaining inside the fort, generating an extra risk for the population of the town. The situation just after the collapse can be observed in Figure 2-8a. In addition, the cliff presents a high frequency of fall of small blocks (less than 0.01 m3), more than 40 events per year. This phenomenon can be associated to the flaking, observed in the face of this rock massif and shown in Figure 23456789.

This response is associated to a physico-chemical degradation of the rock related to the climate actions. The Centre d'Études Techniques de l'Équipement du Sud-Ouest attributes this type of flaking specially to the thermal stress induced by the natural thermal cycles.

In order to reduce the first level of risk, the blocks were removed from the cavern, and a wire mesh was installed in order to avoid all the small and medium (0.01 -2 m 3 ) size blocks detached reach the village, see Finally, the third level of risk considers a much higher volume of material. It considers the stability of the cliff above the cavern, which is controlled by the two principal discontinuities joints and the stratification planes (Figure 2 -11). These joints (family 1 and 2 in Figure 2-11) were developed with the decompression generated by the valley formation geological process.

Figure 2-11 Principal discontinuities families observed [START_REF] Ruiz | Thermo-mechanical analysis of the stability of a rock-cliff under climatic actions[END_REF] As can be observed in Figure 2-11 a large block is delimited by the discontinuities mentioned before. The stability of this block may be affected if the beam inside the fort fails. All the fractures that defined this block are open, and they emerge on the top of the massif.

This suggest the existence of rock-bridges inside the fracture that controls de stability, as was assumed when we consider the stability of the beam. If this large block is supported by rockbridges the climate cycle induced degradation may play an important role for the generation of important rockfalls event. Besides that, the presence of ancient karstification networks increase the instability conditions on the cliff. (Virely and Guittard, 2010b;[START_REF] Virely | Project DOSM Falaises de La Roque Gageac Évaluation de l ' aléa chutes de masses rocheuses[END_REF] To reduce the second level of instability risk and probably the third level too, several portico structures were constructed inside the troglodyte fort (Figure 2-12), with the aim of providing more stability to the roof of the cavern, by stabilizing the instable beam and for instance all the above part of the cliff. After the event registered in January 2010, with the aim to follow the stability conditions of the cliff, an instrumentation system was installed inside and around the troglodyte cavern. This instrumentation system recorded two principal features: weather conditions and deformations (Gasc-Babier, [START_REF] Gasc-Babier | Termal fatigue in rocks-La Roque-Gageac's case study[END_REF].

The closest meteorological station is the Sarlat-La Canéda, is located 10 km north-west from La Roque Gageac. Figure 2-13 presents the temperatures and rains observed in this station from June 2010 to July 2011.

The data recorded by the station indicate that the major rain period occurs during the autumn season, when any rockfall event has been registered. This condition coupled with absence of karstification evidence in the face of the cliff or the aerial surface around the cavern, suggests that water is not the principal degradation factor in this area.

Accordingly, even if the rain effects cannot be completely excluded, the weather variable recorded with the instrumentation was the temperature. The temperature registered in the cliff correspond properly with the regional temperatures recorded by the meteorological station. The deformation in the cliff was follow with two extensometers type RockTest/Télémac distofor @ . These borehole extensometers were located over the cavern (D1 and D2 in Figure 2-15) and were prepared with three displacements measuring points at 2, 4 and 6 m and two temperature measuring points at 2 and 6m.

Additionally, five jointmeters were installed in cracks inside the fort (see . These jointmeters consist in LVDT strain gauges that measure the displacement between their two extremes with an accuracy of 1/100. The extensometer registered the evolution of the rock mass and the jointmeters survey the evolution of the relative displacement of the fissures that delimit the beam:

• Relative opening of the joints located at both supports west and east (jointmeters 1 and 2)

• Vertical aperture of the horizontal discontinuity between the roof of the cavern and the beam (Jointmeter 3)

• Horizontal aperture of the vertical crack that cuts the beam close to the west support (Jointmeter 4)

• Horizontal aperture of the vertical crack in the roof of the cavern, close to the beam (Jointmeter 5) 2.8. Instrumentation results.

The instrumentation system presented in section 2.7 , has been recording the cliff displacements and external temperature for more than 6 years including several months after the construction of the stability structure inside the cavern, presented in Figure 2-12

The extensometers located in the upper part of the cliff, (Figure 2-14) measure the displacement and temperature at 2 different points: 2, and 6 meters from the face of the cliff. Figure 2-16 shows the evolution of temperature inside the cliff for extensometer D1, compared to the exterior temperature. The response observed shows a considerable phase delay, related to the distance between the points of measure and the cliff's face. In fact, between 2 and 6 meters depth, the delay observed is about 40 days.

Besides this delay, the temperature evolution shows an important attenuation in maximal and minimum values as the measure is deeper. But also, the hourly and daily thermal effects decrease. In consequence, the response observed at 6m is only related to seasonal changes. Figure 2-17 presents the displacements response recorded during 6 years with the extensometers installed in La Roque Gageac. In both devices, the response observed is directly related to the external temperature. This condition reinforces the hypothesis of natural thermal conditions as the main cause of deformations and events experienced in the cliff.

Is important to mention that as the instrumentation system was installed during the summer (july 2010) therefore the measures are taken in relation to the cliff situation at this precise moment, therefore, the amplitude for negative displacements (compression) are higher than the amplitude for positive displacements (expansion). In general, the displacements recorded with the extensometer D1 are higher than those of extensometer D2. It is interesting to note, that measures at 2 and 4 meters in extensometer D1 have the same behavior, while for extensometer D2 the displacements at 4 meters have a smaller amplitude than those registered at 2 m, as it is expected when considering the thermal conductivity and deformations inside the cliff.

The response recorded with the device D1 presents a considerable accumulation of plastic displacements at 2 and 4 m. After 6 years, the cumulated irreversible displacement is about 0.5 mm. While the displacements recorded at 6 meters don't show any accumulation, they are almost imperceptibles.

For the case of extensometer D2, displacements taken at 2 meters present a plastic accumulation of 0.2 mm after 6 years, that is less than the half of what is experienced with the device D1. As well, strains recorded at 4 m and at 6m have an elastic response, but also, displacements variation at 6 m are almost imperceptibles, just like what is observed with device D1.

The difference in the behavior observed in both devices D1 and D2 may be explained by the presence of a discontinuity, that go through the extensometer D1 in a point between 4 and 6 meters, as illustrated in Figure 2-18, inducing that displacements at 2 and 4 meters in this device present a block behavior, while in the extensometer D2 a regular thermal induced response is observed. These blocks suffer contraction during cold seasons and consequently an increasing in the joint aperture, as illustrated in Figure 2-20. The opposite behavior will be observed during hot seasons, Therefore, the negative values in joint displacements presented in Figure 2-19 represents the opening of the joint.

As was shown in Figure 2-15, the jointmeters J-1 and J-2 follow the behavior of the cracks located at the supports of the beam, and jointmeters J-3, J-4 and J-5 survey the response of the roof. It is important to mention that jointmeter J-1 was initially installed in horizontal position to measure the aperture of the vertical crack situated at the east support, but in December 2010 this device was re-localized in a horizontal segment of the same crack, which is its current position. The jointmeters located in the roof (J-3, J-4 and J-5) present a large displacements amplitude (2.5 mm). Nevertheless, the response observed is elastic; no irreversible deformation is observed. As example, Figure 2-21 presents the relation of displacements with temperature for device J-5. It can be observed that the displacements recorded with this device are completely reversible despite its important amplitude. Conversely, the devices J-1 and J-2 that follow the discontinuities located at the supports of the beam, show an important accumulation of displacements, even if the displacements amplitude register is smaller than the one in the roof jointmeters. Figure 2-22 presents the evolution of jointmeters J-1 and J-2 for the whole instrumentation period. The observed deformation can be considered as stable after July 2015, which is when the portico structures showed in Figure 2-12 where constructed. However, the response of the roof jointmeters and the extensometers don't seem to be affected by these structural elements. As mention before during the cold seasons the joints opens, these apertures reduce the area of support between the blocks that compose the joint. This can induce a local stress concentration and increase the probability of a failure. For this reason, some irreversible deformation can be associated to a low temperatures threshold. In fact, the considerable displacements drop (augmentation in joint apertures) observed in February 2012 can be explained by an unusual extreme cold period.

Actually, 2012 was the year with the more irreversible deformations recorded by the jointmeters located at the beam supports. 2.9. In situ concluding remarks.

The cliff from La Roque Gageac has been submitted to recurrent rockfalls, in the last century four important events were registered, that can't be attributed to any evident triggering factor. The last event was registered in 2010, leading a poor stability condition in the cliff, especially the zone where the troglodyte cavern is placed.

For this reason, an instrumentation system was installed in the site and have been following the cliff response for the last 6 years, suggesting that atmospheric thermal cycles have a strong relationship with the deformations of the cliff.

Other specific observations are mentioned bellow:

• Because the orientation of the cliff (south) it can be considered that sun radiation may have an important role in the response of the cliff, in fact, all the flacking observed at the face of it is considered to be related to daily thermal variations suffer by the cliff.

• Despite the calcareous nature of the rock that composed this cliff (limestone) the water doesn't seem to have a predominant effect on the response of the cliff. Actually, the dissolution rate observed in the region is less of 1 mm per year.

• The measurements of temperature measured inside the rock mass shows a considerable phase delay, related to the distance to the face cliff, moreover, as the measurement is deeper the values of maximal and minimal temperature are attenuated, and the daily effect is lost, in fact for the temperature evolution observed at 6m is only related to seasonal changes.

• The extensometers placed above the troglodyte cavern shows an important strains accumulation with the increment of seasonal thermal cycles. Furthermore, these measurements suggest the presence of an important discontinuity.

• The jointmeters measurements indicate that during cold seasons the fractures tends to open while during hot season they close. This may explain why the event registered in 2010 occurred during the winter.

• The jointmeters located at the fissures that defined the residual beam, shows an important relationship between the crack opening and closing with the temperatures. Moreover, an extraordinary cold period registered in the winter of 2012 generated an irreversible opening of these fissures.

• A portico structure was within the troglodyte cavern with the purpose of give more stability to the residual beam and possibly to the upper part of the cliff, based on the instrumentation data registered, it can be observed that the fissures around the beam present almost any additional evolution, however the response of the extensometer placed above the cavern are not affected by this structure.

Experimental Work

Introduction.

The behavior of a rock massif is the result of the combined response of the rock matrix and the discontinuities that divide the rock mass, delimiting rock blocks. The instrumentation system installed in La Roque Gageac's cliff and presented in chapter 2 is mainly focused on the response of the discontinuities observed in the cliff, showing an important relationship between of deformations and thermal atmospheric cycles.

The cliff position (oriented to the south) and the important thermal range experienced on the site, suggest that the rock matrix may also be affected by those atmospheric thermal cycles.

In fact, the flaking phenomenon observed in the cliff's face was initially related to the thermal stress induced by the weather conditions.

Therefore, an experimental program was developed to evaluate the possible damage induced in the rock matrix due to climatic thermal conditions. In order to evaluate the effect of thermal degradation, different rock properties are considered to examine its evolution through thermal cycles with an atmospheric range.

Rocks weathering due to thermal stress can appear under two different phenomena: the thermal fatigue; that is related to temperature variations that lead to repeated stress, below the rock strength, but its cycling effect induces permanent damage in the material. In the other hand, the thermal shock is related to a single event of a sudden change of temperature that generates an excess of stress that rock is incapable of absorbing other than through instantaneous failure.

The thermal degradation of rocks is mostly influenced by the different thermal properties of the minerals that compose the rock. Differences in thermal expansion, thermal conductivity and principal axis direction, controls the mineral to mineral thermal stresses induced by the thermal variations. Hence, additional stresses can be induced leading to failures across the minerals boundaries, and lines of weakness, generating some intergranular disintegration and propagation of existent fissures. [START_REF] Hall | The role of thermal stress fatigue in the breakdown of rock in cold regions[END_REF][START_REF] Siegesmund | Physical weathering of marbles caused by anisotropic thermal expansion[END_REF][START_REF] Walsh | Micromechanical modeling of thermal spallation in granitic rock[END_REF][START_REF] Hall | Thermal fatigue and thermal shock in bedrock: An attempt to unravel the geomorphic processes and products[END_REF] Under laboratory conditions, the thermal deterioration of rocks has been mostly studied for high temperatures (more than 100ºC) like those experienced during fires [START_REF] Brotóns | Temperature influence on the physical and mechanical properties of a porous rock: San Julian's calcarenite[END_REF] or in nuclear waste disposal [START_REF] Chen | Experimental study on mechanical properties of granite after freezethaw cycling[END_REF]. Nevertheless, some studies are focused on weather induced phenomena like freeze-thaw [START_REF] Malaga-Starzec | Microscopic and macroscopic characterization of the porosity of marble as a function of temperature and impregnation[END_REF], thermal shock [START_REF] Yavuz | Effect of freeze-thaw and thermal shock weathering on the physical and mechanical properties of an andesite stone[END_REF] or cold -heat cycles [START_REF] Goudie | The thermal degradation of marble[END_REF] The elastic wave propagation velocity is a common parameter used to evaluate damage in rocks, since, this parameter is extremely sensitive to the presence of defects in the sample, (like the internal fissures) due to the considerable difference between the air and the solid materials wave propagation velocity [START_REF] Martínez-Martínez | Spatial attenuation: The most sensitive ultrasonic parameter for detecting petrographic features and decay processes in carbonate rocks[END_REF][START_REF] Brotóns | Temperature influence on the physical and mechanical properties of a porous rock: San Julian's calcarenite[END_REF][START_REF] Hudson | Influence of thermal damage on linear and nonlinear acoustic properties of granite[END_REF].

There is a relation between the elastic wave propagation velocity and the dynamic Young modulus and Poisson ratio, allowing the comparison with the results obtained with mechanical tests [START_REF] Kendrick | The influence of thermal and cyclic stressing on the strength of rocks from Mount St. Helens, Washington[END_REF]. Moreover, this parameter is obtained with a non-destructive test, hence, is possible to evaluate its evolution in the same sample through the whole laboratory study.

This parameter has been widely used to evaluate the effect of thermal stresses over different lithologies, like: marbles (Luque et al. 2011;[START_REF] Yavuz | Estimating the index properties of deteriorated carbonate rocks due to freeze-thaw and thermal shock weathering[END_REF][START_REF] Malaga-Starzec | Microscopic and macroscopic characterization of the porosity of marble as a function of temperature and impregnation[END_REF][START_REF] Yavuz | Thermal effect on the physical properties of carbonate rocks[END_REF], limestones [START_REF] Yavuz | Estimating the index properties of deteriorated carbonate rocks due to freeze-thaw and thermal shock weathering[END_REF][START_REF] Martínez | Influencia De La Alteración Sobre Las Propiedades Mecánicas De Calizas, Dolomías Y Mármoles[END_REF][START_REF] Yavuz | Thermal effect on the physical properties of carbonate rocks[END_REF]( or Andesites [START_REF] Yavuz | Effect of freeze-thaw and thermal shock weathering on the physical and mechanical properties of an andesite stone[END_REF]. Even if these studies don't consider similar thermal ranges, in all of them a decay in the elastic wave propagation velocities with the imposition of thermal cycles is observed, suggesting an internal micro-fissuration of the samples, even if there is no evidence of external fissures [START_REF] Yavuz | Estimating the index properties of deteriorated carbonate rocks due to freeze-thaw and thermal shock weathering[END_REF].

Another interesting physical property to consider is the effective porosity. In fact, the initial effective porosity of the material may influence the thermal weathering of the rocks, because if the material presents a high porosity, the expansion of the minerals may be absorbed by the porosity, reducing the internal stresses induced [START_REF] Sousa | Influence of microfractures and porosity on the physicomechanical properties and weathering of ornamental granites[END_REF]. Additionally, an augmentation in the values of effective porosity and variations in the pore size distribution can also be related to the fissuration of the rock due to thermal weathering. For example, [START_REF] Koch | The combined effect of moisture and temperature on the anomalous expansion behaviour of marble[END_REF] observed an augmentation in the pore size radius and porosity percentage of marbles after imposing 40 cycles between 20ºC to 80ºC. In the same way, Luque et al. (2011) and [START_REF] Franzoni | Artificial weathering of stone by heating[END_REF] obtained a variation in the pore size distribution of marbles and granites respectively after inducing thermal degradation.

As mentioned before, the thermal weathering is a process that begins in the mineral to mineral boundaries. To evaluate the phenomena at a mineral scale, environmental scanning electron microscopy (ESEM) and scanning electron microscopy (SEM) images are generally used. For instance, Luque et al. (2011) present the evolution of a specific fissure over a marble while it is heated between 20ºC to 90ºC, showing the widening of the fissures while the temperature increase. In the same way, [START_REF] Goudie | The thermal degradation of marble[END_REF] use SEM observations to show that mineral boundary wide increase as the samples are subjected to a higher number of thermal cycles.

At macroscopical scale the mechanical properties of the rocks may also be affected by the thermal weathering. This phenomenon has been evaluated using different mechanical tests types such as: Brazilian test [START_REF] Mahmutoglu | Mechanical Behaviour of Cyclically Heated Fine Grained Rock[END_REF], Schmidt hardness [START_REF] Yavuz | Estimating the index properties of deteriorated carbonate rocks due to freeze-thaw and thermal shock weathering[END_REF][START_REF] Yavuz | Effect of freeze-thaw and thermal shock weathering on the physical and mechanical properties of an andesite stone[END_REF] or uniaxial compression test [START_REF] Mahmutoglu | Mechanical Behaviour of Cyclically Heated Fine Grained Rock[END_REF][START_REF] Franzoni | Artificial weathering of stone by heating[END_REF][START_REF] Chen | Experimental study on mechanical properties of granite after freezethaw cycling[END_REF]. The laboratory study developed by [START_REF] Mahmutoglu | Mechanical Behaviour of Cyclically Heated Fine Grained Rock[END_REF] over samples of marble and limestones, heated between 20ºC to 600ºC shows that uniaxial compression and tensile strength reduce when samples are exposed to thermal cycles. In fact, marble samples show a reduction of 50% of their Uniaxial compressive strength and 80% of their tensile strength after 16 cycles. This response was also observed when smaller thermal ranges were considered, Chen et al., 2014 presented the results obtained for granite samples subjected to thermal cycles between -20ºC and 65ºC, obtaining a reduction of 30% in the Uniaxial compressive strength value after 40 thermal cycles.

In the same way, [START_REF] Yavuz | Estimating the index properties of deteriorated carbonate rocks due to freeze-thaw and thermal shock weathering[END_REF] obtain the response of a marble and a limestone subjected to thermal cycles between -20ºC and 20ºC. This experimental program evaluates the mechanical response of the material but also some physical properties such as porosity and P-wave propagation velocity (Vp). Showing that when this material experiences 20 cycles the strength and Vp value of the rock decrease and the porosity increases.

In order to evaluate the thermal weathering over the limestone of La Roque Gageac (LaRG), two groups of samples were considered: the first group correspond to samples obtained from remaining blocks from inside the troglodyte cavern after the 2010's event, and mainly conformed by calcite. The second group is composed by samples directly drilled from the face of the cliff nearby the troglodyte cavern; these samples are composed by quartz and calcite.

Both groups of samples were subjected to thermal cycles with a thermal range of 40ºC (10ºC to 50ºC). The effect of these thermal cycles over this limestone was evaluated through different physical and mechanical properties, such as deformations, elastic wave propagation velocities, fissuration and propagation of fissures, and Uniaxial compressive strength. In this chapter, a description of the complete experimental program and its results are presented.

La Roque Gageac's limestone general characterization.

After the event of 2010 some physical and mechanical tests were performed, to have a complete characterization of La Roque Gageac (LaRG) limestone. These tests were realized on samples obtained from 2 blocks dropped from the roof of the troglodyte cavern during the 2010's event. It is important to mention that samples obtained from the block defined as A were obtained in the perpendicular direction to the stratification planes, however, block B doesn't show any evidence of the stratification planes direction.

The results obtained are presented in the report of the DO-SMS project report realized by Centre d'Études Techniques de l'Equipement du Sud-Ouest (now Cerema) (Virely and Guittard, 2010a).

Physical properties.

The density and elastic waves propagation velocity of the rock were obtained from 13 different samples. The mean density obtained is 2382 kg/m3 with a standard deviation of 45 kg/m3. The mean P-wave propagation velocity (Vp) is 4651 m/s, and the shear-wave propagation velocity, obtained from measures realized in two perpendicular directions, has a mean value of 2671 m/s. Table 3-1 gives the results obtained for the 13 samples. To complete the physical characterization of La Roque Gageac limestone, the porosity and the saturated water content were obtained. These results are summarized in Table 3-2. The mechanical properties of LaRG limestone were determined by 3 types of tests: Uniaxial compressive strength, Brazilian tensile strength and triaxial test. Table 3-3 summarizes the results obtained.

Uniaxial compressive strength was obtained from 10 different samples: 7 from block A and 3 from block B, the mean value obtained is 47.39 MPa with a standard deviation of 9.63. This last value is quite high, this characteristic may be related to the heterogeneity of the material, associated to its geological composition. From Table 3-3 can be appreciated that samples 1017 and 1018 present the lowest UCS values obtained from the 10 samples. This can be related to the fact that those two samples presented a vertical discontinuity, as shown in Figure 3-1. In fact, the mean UCS value for this two samples, is 33.30 MPa. Whereas, the mean value obtained for all the samples is 47.39 MPa. If these two samples are not considered, the mean UCS value increase to 50.87 MPa and the standard deviation becomes 7.32. Nevertheless, such a standard deviation value is still associated to an heterogenous rock. Based on the mechanical properties obtained and the in-situ observation, the Laboratoire Régional des Ponts et Chaussées de Toulouse performed a RMR evaluation. Table 34summarize the parameters determined. [START_REF] Hoek | Practical estimates of rock mass strength[END_REF] proposed a modification of the original [START_REF] Hoek | Empirical Strength Criterion for Rock Masses[END_REF] failure criterium defined for intact rocks, to include the real conditions of the massif, like joins and discontinuities. This modification considers the geological strength index (GSI) as a parameter to reduce the strength of the material based on the geological condition of the massif. The modified Hoek-Brown failure criterion is defined by Eq. 3-1:

σ′ 1 = σ′ 3 + σ ci (m b σ′ 3 σ ci + s) a Eq. 3-1
Where, ci is the Uniaxial compressive strength of the rock and mb s and a are parameters that depend upon characteristics of the rock mass, based on the GSI parameter. For LaRG limestone the GSI determined is 76 associated to the blocky structure of the massif. Under this consideration the parameters consider for modified Hoek-Brown failure criterion are: mb = 4.24, s = 0.0695 and a= 0.501.

Thermal properties.

Thermal expansion.

The thermal expansion of the material was determined with the procedure proposed by the French standard NF-EN 14581 method B (AFNOR, 2005). For this, strain gages were glued in vertical (vert) and horizontal (hor) direction on two prismatic and two cylindrical samples.

This standard considers than samples should be placed in an oven at 20ºC for two hours, then the temperature increases until 80ºC with a rate of 0.5ºc/min and remain constant at 80ºC for two more hours. Strains should be measured at each change of temperature. Finally, thermal expansion is defined as:

α i = ε i ΔT Eq. 3-2
Where, i is the thermal expansion in direction i,i represents the strains induced by the temperature variation, measured for direction i and T is the thermal variation (60ºC). Table 3-5 presents the values of thermal expansion obtained with the four samples considered. Note that for prismatic samples two different pairs of gauges were used to measure horizontal deformation. It is important to mention that Ruiz (2013) determined a valued of thermal expansion for the rock massif equal to 6x10 -6 ºC -1 . This value was obtained with a back-analysis from the instrumentation measurements in the cliff.

Samples description.

To evaluate the damage of LaRG limestone due to thermal cycles, two geometries of samples are considered: cylindrical and prismatic. The use of each geometry depends on the type of test to be performed on the sample.

Furthermore, for this experimental program, two groups of samples are considered. The first set corresponds to samples drilled from 2 blocks (see Figure 34) fell during the 2010 event and removed in December 2014 (block 1) and February 2015 (block2). Each sample obtained in vertical direction was obtained from one borehole, with cores of maximum 20 cm long, as illustrated in Figure 3456. Therefore, all vertical samples belong to the first 20 cm of the cliff.

Longer cores were obtained from the horizontal boreholes, almost 40 cm long. Thus, several horizontal samples were obtained from the same rock core. For this reason, even if these samples still be part of the most superficial zone of the cliff, some of them belong to a deeper zone.

In total 35 samples are evaluated, 30 cylindrical and 4 prismatic. Porosity of the samples is obtained based on the density of the samples and the specific weight obtained in laboratory: 2.72 g/cm 3 for samples from blocks and 2.76 g/cm 3 for samples directly drilled from the face of the cliff.

Mineralogy.

The mineralogical composition of each group of samples was determined, to evaluate if there are any mineralogical differences between them. For this, three types of tests were performed: spectrometry, X-Ray diffraction and thin sections analysis. X-ray was only performed over samples belonging to the second group. See Table 345678The samples used for the mineralogical analysis came from cylindrical samples prepared for this experimental program but rejected because they didn't follow the geometry conditions for the mechanical tests (parallelism between faces and/or dimensions). In total 5 samples were used.

M1 Group 1 -block 1 X X M2 Group 1 -block 2 X X M3 Group 2 X X M4 Group 2 X M5 Group 2 X
Infrared spectrometry.

Infrared spectrometry (IRA) analysis were performed over three samples: M1, M2 and M3. This technique is based on the evaluation of the behavior of electromagnetic waves that are reflected or absorbed by a solid body. The percentage of absorption and reflectance are controlled by physical and chemical characteristics of the body.

The spectrum obtained for each sample can be compared to typical spectra of different minerals in order to determine the principal minerals that compose the rock. This technique doesn't give a mineral content percentage but allows to identify the principal composing mineral if it is predominant.

Figure 3-7 presents the spectrometry results obtained for the three analyzed samples, compared with the typical characteristics bands observed for calcite [START_REF] Huang | The Infrared Study of the Carbonate Minerals[END_REF][START_REF] Vazquez-Moreno | Tabla de frecuencias y relacionados con la química del cemento[END_REF][START_REF] Gunasekaran | Raman and infrared spectra of carbonates of calcite structure[END_REF]. Note that spectra obtained for sample M1 and M2 follows precisely the calcite bands. However, spectra of sample M3 doesn't match either with typical calcite or quartz bands. Therefore, it can be concluded that samples from the face of the cliff (sample M3) are not composed by just one principal mineral and samples M1 and M2 are mainly composed by calcite. Thin sections.

In order to complement the information obtained from the IRA and XDR tests and to evaluate mineralogical characteristics of the two groups of samples, thin sections were created from samples of each group. Those thin sections were scanned under natural and polarized light, the mineralogical analysis is performed over these images.

The analysis of the thin section images allows to obtain a mineral composition percentage, but also, the shape and size distribution of minerals composing the rock. For this, an imagebased analysis is performed using the software SPO [START_REF] Launeau | Shape preferred orientation of object populations: automatic analysis of digitized images[END_REF][START_REF] Launeau | Fabric analysis using the intercept method[END_REF]. This program works with processed images, in which each color represents a mineral. Therefore, it allows to determine the area percentage composed by each mineral. Additionally, based on the intercepts method [START_REF] Launeau | Magmatic fabric acquisition mechanisms in a syenite: Results of a combined anisotropy of magnetic susceptibility and image analysis study[END_REF] and the shape preferred orientation (SPO), the mean shape of the minerals and their shape ratio distribution can be determined.

As the minerals have an irregular shape, the SPO analysis transforms this shapes in ellipsoids.

In this way is possible to determine the mean shape of the minerals and the shape ratio distribution, considered as the relation length -width of the ellipsoids. The principal results obtained with the SPO analysis are: the percentage area cover by each mineral and the shape ratio distribution. From these characteristics is possible to determine the predominant minerals and the relative size of the minerals that compose the rock.

Two thin sections were created from both samples belonging to group 1 (sample M1 and M2) and only one from group 2 (M3). From thin sections of samples M1 and M2, it is observed that these two samples are mainly composed by calcite, just like was obtained from the spectrometry test. Calcite is present in two types of crystallinity: micrite, darker and without crystalline shape and sparite, that is the lighter part with coarse crystalline shape. In fact, almost any other type of mineral is observed, just a little portion of quartz (colorless minerals with low positive relief in natural light).

For instance, Figure 3-11 presents the typical appearance of the calcite and quartz observed in sample M2. Some fossils are identified in samples M1 and M2. these fossils are mainly composed by calcite in micrite crystallinity. In the same way, sample M3 also shows minerals of calcite, in both types of textures, sparite and micrite. However, an important percentage of voids are observed in the thin section, as illustrated in Figure 3-12 (light grey/light pink zones). These voids can represent either the sample porosity of the presence or minerals that were lost during the preparation of the samples. This last cause is the one assumed to be the explanation for the voids observed in sample M3 because the porosity of LaRG limestone is not that high. In fact, for this sample, the porosity measured is only 9%.

Considering the DXR results it was observed that the sample M3 presents an important quartz content. For this reason, it is assumed that the voids observed in the thin sections represent the quartz content. Unfortunately, from this sample is not possible to evaluate the quartz crystals shape ratio, as it is not possible to know if each void was composed by just one crystal or several crystals. Taking into account the considerations presented before, the modal composition of the samples is presented in Figure 3-13. It can be observed that calcite is the predominant mineral for the three samples. Nevertheless, for samples obtained from inside the cavern (M1 and M2) calcite is mainly observed in sparite texture and for the sample obtained from the face of the cavern (M3) in micrite. In the other hand, sample M3 presents a considerable quartz content in comparison with samples M1 and M2. From the shape ratio distribution, it is possible to determine the predominant shape ratio for a specific mineral. Figure 3-14 presents the shape ratio distribution for the calcite in micrite shape, it can be observed that samples obtained from inside the cavern have a predominant shape ratio of 1.4, while this value for the sample drilled in the face is 1.6, which means that micrite crystals for sample M3 are longer than for samples M1 and M2. It is interesting to point out that the large number of objects observed in sample M1(b) may be related to the fact that this sample also experienced the highest percentage of micrite content. In the same way, the number of micrite objects observed for sample M3 is similar to the values observed for samples M1(a) and M2, nevertheless, its micrite percentage is higher. Suggesting that the micrite crystals observed in the face of the cliff are bigger than those of samples drilled from blocks within the cavern.

In the other hand, the shape ratio distribution for sparite is presented in Figure 3-15, where can be observed that for all samples evaluated the predominant shape ratio is 1.5. In other words, most of sparite crystals are 50% longer than width. Note that the number of objects observed for sample M2 are the lowest of all the samples analyzed. If we consider that this sample also present the major sparite content (Figure 3-13), can be concluded that the sparite crystals observed in this sample are bigger than in the other samples. Sample M3 is the sample with the lowest sparite content, only 20%, nevertheless, its number of objects is quite important, more than sample M2 and almost the same that sample M1(b), suggesting that its crystals are smaller than those observed in samples obtained from within the cavern.

Finally, Figure 3-16 presents the shape distribution for quartz, it is important to remember that because of the thin sections characteristics from sample M3 is not possible to obtain the quartz shape ratio distribution for this sample. It can be observed that the predominant shape ratio in samples M1 and M2 varies from 1.3 to 1.5, some differences are observed between thin sections a and b of samples M1 and M2, nevertheless, as the number of objects (Figure 3-16) and quartz content (Figure 3-13) is very low, it can be considered that these differences are negligible.

Besides the mineralogical composition, another important condition that can be analyzed from the thin sections is the minerals fissures density. This characteristic may be important for our analysis, because the presence of pre-existing fractures may indicate a more fragile behavior, hence, a large damage level after the thermal weathering. Not to mention that the existing crack density level may control the mineralogical thermal response of the rock [START_REF] Hall | Thermal fatigue and thermal shock in bedrock: An attempt to unravel the geomorphic processes and products[END_REF]. To evaluate this fissuration density an equal grid should be draw in each thin section image.

Then, the number of cracks is count in each cell. Finally, the density of pre-existing cracks is calculated, considering the area of each cell. Figure 3-17 presents the cracks density observed for the thin sections evaluated, where any particular difference is observed between them, except by the fact that sample M1 shows a slightly inferior quantity of fissures by mm2.

From visuals observations can be determined that sample M1 also shows a major quantity of fossils. This condition may explain why the cracks density observed is slight inferior, as the fossils don't have fissuration.

Mineralogical analysis remarks.

The three types of mineralogical analysis performed shows similar results. Spectrometry and DRX results corroborate what is observed with the most visual technique, thin sections images. Samples obtained from blocks inside the cavern (group 1) are principally constituted by calcite, and samples from the cliff face (group 2) are composed by calcite and quartz in a relationship that can reach the 50-50. This mineralogical characteristic may lead to a different mechanical response. More precisely samples from group 2 can present more thermal degradation associated with the strong thermal interaction between the calcite and quartz minerals, as mentioned by [START_REF] Hall | Thermal fatigue and thermal shock in bedrock: An attempt to unravel the geomorphic processes and products[END_REF].

On Table 3-6 it can be observed that samples from group 1 present a mean elastic wave propagation value of 4400 m/s, a similar value to the one obtained by the Laboratoire Régional des Ponts et Chaussées de Toulouse (chapter 3.1.1) over samples obtained from blocks from inside the troglodyte caver. Nevertheless, the mean P-wave propagation velocity obtained for samples from group 2 is 2900 m/s. This variation in P-wave velocity can be associated with the mineralogical composition of the rock, but also, to the fact that samples from group 2 were obtained near the face of the cliff and thus present a higher initial deterioration.

Experimental layout.

In order to evaluate the thermal weathering of LaRG limestone, an experimental program is developed. This program is based on the laboratory study performed by Gasc-Barbier et al. (2014), in which four samples of limestone were submitted to 700 thermal cycles and measurements of deformation and elastic wave propagation velocity were performed to evaluate the effects of thermal fatigue in the samples.

The experimental program consists in submitting the samples described in Table 3-7 andTable 3-6 to thermal cycles. As the aim of this analysis is to evaluate the weathering of the rock matrix, the thermal rang employed is close to the one experienced in atmospheric conditions.

The effect of thermal cycles was evaluated through the measures of deformation (Def), elastic wave propagation velocity (Vp-Vs) Uniaxial compressive strength (UCS) and the evaluation of generation and propagation of fissures by digital image correlation technique (DIC). The number of cycles applied for each sample depends on the type of measurement that would be performed. Table 3-9 summarizes the experimental layout of each sample.

X X B6 828 X X B7 828 X X B8 624 X B9 624 X B10 624 X B11 420 X B12 420 X B13 420 X B14 210 X B15 210 X B16 210 X B17 0 X B18 0 X B19 0 X B20 396 X CV1 910 X X CV2 1222 X CV3 1222 X CV4 910 X X Sample Cycles applied Def Vp-Vs UCS DIC CV5 1222 X CH1 1222 X CH2 910 X X CH3 1222 X CH4 1222 X CH5 910 X X C1 948 X C2 948 X C3 948 X C4
1080 X

Heating -cooling treatment.

The thermal cycles were imposed with the climatic chamber presented in Figure 3-18. With this climatic chamber, thermal cycles can be program for the whole experimental program, it allows temperatures between -45ºC to 180ºC. Additionally, it is equipped with a humidor that allows to impose wet thermal cycles. However, the relative humidity cannot be controlled. Our experimental program is developed under dry conditions. During the experimental process, it was notice that the temperature felt by the samples at the surface was not exactly the one imposed, because of the thermal regularization system of the climatic chamber. To assure the temperature imposed on the samples, the variations of temperature were measured with a thermal probe as shown in Figure 3-20. The probe was fixed to the face of sample B4, to measure the real imposed temperature, then the probe was placed in a small borehole performed in the center of sample B4, to measure the temperature in the specimen. Therefore, all the possible damage that may be induced in the samples because the heating/cooling cycles will be related to a process of thermal fatigue. This means, that in this experimental program it is only considered the rock deterioration related to typical thermal variations, and not to exceptional changes of temperature, like the one registered in February of 2012 that generated a significant plastic deformation, observed in the joint-meters installed in the cavern. (see chapter 2.8).

Figure 3-21 Temperature variation experienced by the samples Thermal conductivity could not have been measured in the laboratory, as no layout was available. Nevertheless, as this parameter may present some valuable information, especially to the modeling phase of this work, a back-analysis was performed to estimate it, based on the temperature measurements on sample B4. This back-analysis was performed with the finite element method software CODE_BRIGHT [START_REF] Olivella | Numerical formulation for a simulator (CODE_BRIGHT) for the coupled analysis of saline media[END_REF].

The value of thermal conductivity obtained was: 1 W.m.K-1 considering a specific heat of 800 J.Kg-1.K-1, because the sample is dry. Figure 3-22 presents the good agreement obtained from the modeling with these values of thermal conductivity and specific heat. Likewise, Ruiz (2013) obtained from back-analysis from the instrumentation data, a value of thermal conductivity for the rock cliff, a value of 2 W.m.K-1 with a specific heat of 1200 J.Kg-1.K-1. Deformation measurements are performed with strain gauges. Two pairs of strain gauges, with a nominal resistance of 120 Ω, were glued on each sample: a pair in vertical position, to measure the axial deformation, and a pair in horizontal position, to measure the radial deformation, as illustrated in Figure 3-24. As the strain gauges are submitted to thermal variations, they are installed in half-bridge circuit with a thermal compensation.

Elastic wave propagation velocity.

The elastic wave propagation velocity is performed at least once a week. It is also a preliminary measure that it is performed for all cylindrical samples to have an initial characterization and comparison parameter (see Table 3456).

The measure is performed by the ultrasonic pulse technique described by the French Standard NF P 94-411 [START_REF] Afnor | Détermination de la vitesse de propagation des ondes ultrasonores en laboratoire[END_REF], that requires the equipment presented in Figure 3-25. For this, an electric pulse generator that produces waves higher than 250 volts for a maximum rise duration of 0.2 ms is used to induce compression waves (P-waves) and shear waves (S-waves) in two perpendicular directions (S1 and S2).

This generator is connected to the transmitter and receptor transducers with a resonance frequency from 25 kHz to 1 MHz. The sample is fixed between the transmitter and receptor transducers in a frame provided with a uniaxial pump, to assure good contact between the sample and the transducers, without applying a mechanical load to the sample. A force transducer is used to control that the applied load is 2 kN. To have a proper wave transmission between the sample and the transducers, very thin sheets of lead are fixed with honey on the transducers face. Finally, the signals are registered by an oscilloscope with a time resolution of 10 ns connected to a personal computer equipped with a signal registration software. From the signal, registered by the equipment presented before, is possible to obtain the time of arrival of the wave, as illustrated in Figure 3-26.Then the elastic wave velocity is calculated as:

V i = L i t i Eq. 3-3
Where, Vi is the elastic wave propagation velocity for samples i, Li is the sample large and ti the wave arrival time. From the measures of compression (Vp) and shear (Vs) elastic wave propagation velocity it is possible to determine the dynamic Young modulus (Ed) and dynamic Poisson ratio (d) following the next equations, where is the bulk density:

E d = ρ V s 2 (3V p 2 -4V s 2 ) (V p 2 -V s 2 )
Eq. 3-4

υ d = (V p 2 -2V s 2 ) 2(V p 2 -V s 2 )
Eq. 3-5

Uniaxial compressive strength.

Uniaxial compressive strength test is performed only on samples from group 1, following the ASTM D7012-10 standard (ASTM, 2010). An MTS system (Figure 3-27) is used, with a stiffness of 109 N/mm, able to perform compression and shearing tests. To evaluate the static Young modulus loading/unloading cycles were realized at 50% and 80% of the Uniaxial compressive strength previously estimated (47.3 MPa) and presented in chapter 3.1.2 using a loading/unloading rate of 0.15 MPa/s.

Fifteen cylindrical samples are considered to evaluate the effect of thermal cycles in the rock strength. They are divided into groups of tree samples submitted to a different number of cycles (0, 210, 420, 624 and 828) as detailed in Table 3-9. During the laboratory study presented by Gasc-Barbier et al. ( 2014) also performed on LaRG samples, some fissures were observed in the samples, as shown in Figure 3-28. In order to evaluate if the thermal cycles are capable to induce that kind of fissures, a specific test is considered, to determine the apparition of fissures and the propagation of existing ones. For this, the samples (see Table 3-9) are photographed each day. These photos are analyzed with the Digital image correlation software Ncorr [START_REF] Blaber | Ncorr : Open-Source 2D Digital Image Correlation Matlab Software Ncorr : Open-Source 2D Digital Image Correlation Matlab Software[END_REF], an open source code that works in a Matlab environment.

The digital image correlation (DIC) technique has been widely used to evaluate the behavior of rocks under compressive loads [START_REF] Dautriat | Localized deformation induced by heterogeneities in porous carbonate analysed by multi-scale digital image correlation[END_REF][START_REF] Yang | Experimental investigation of the delayed behavior of unsaturated argillaceous rocks by means of Digital Image Correlation techniques[END_REF]. It allows to determine the apparition of punctual defects, but also, the general strains field of the sample even when they are low. DIC analysis is mostly used to evaluate the apparition and evolution of fissures. For instance, it is commonly used to analyze the responses observed with the Brazilian disk test, from photos taken with a regular digital camera: more than 50 m/pixel [START_REF] Hild | Digital image correlation: from displacement measurement to identification of elastic properties-a review[END_REF][START_REF] Lagneau | Crack propagation in reservoir rocks in the CO 2 storage context[END_REF]. But also from macro-scale photographs as in the study presented by [START_REF] Casperson | Investigation of thermal effects on fatigue crack closure using multiscale digital image correlation experiments[END_REF] in which DIC analysis on macro photographs (2 m/pixel) were used to observe the closure of a crack under high temperatures (600ºC).

The principal goal of this technique is to obtain displacements and strains fields for a material undergoing deformation, through images processing, as the images are taken while the material deforms.

In brief, the DIC method consists in obtaining one to one correspondence between points in the image of reference and in the current condition. For this, the images are divided in small sections called subsets and deformation of each subset is obtained by the transformation used to match the position of the subset in the current image [START_REF] Hild | Digital image correlation: from displacement measurement to identification of elastic properties-a review[END_REF][START_REF] Bornert | Assessment of digital image correlation measurement errors: Methodology and results[END_REF][START_REF] Blaber | Ncorr : Open-Source 2D Digital Image Correlation Matlab Software Ncorr : Open-Source 2D Digital Image Correlation Matlab Software[END_REF], as follows:

x ̃cur i = x refi + u rc + Where x refi and y refi are x and y coordinates of an initial reference subset point, x refc and y refc are the X and Y coordinates for the center of the initial reference subset point x ̃curi and y ̃curi are the x and Y coordinates of a final current subset point, and, P is the deformation vector.

Then, the deformation of each point in the subset will be defined by the optimal vector P. That is to say, a deformation vector that gives the best match between xref i and xcur i and yref i and ycur i.

To realize properly this analysis, the photos should be taken with the samples exactly in the same position. For this purpose, a specific device to fix the samples was created (Figure 3

Mechanical response.

To evaluate the damage induced in LaRG limestone due to atmospheric thermal cycles, the mechanical response of samples of this rock was evaluated through the temperature variations. The response obtained for deformations, elastic wave propagation velocity and Uniaxial compressive strength is presented below.

Deformations.

Three samples obtained from blocks fell inside the cavern (group 1 -samples B1, B2 and B3) and four samples drilled from the face of the cliff (group 2), two horizontally drilled: samples CV1 and CV4 and two vertically drilled: samples CH2 and CH5 were considered for deformation measurements (see Table 3456789).

As mentioned before, these measures are performed with strain gauges, allowing to have the response even when the thermal cycles are imposed. Is important to consider that positive values of deformation represent expansion and negative values contraction.

General tendency.

Axial strains.

Axial deformations evolutions are presented in Figure 3-30 for samples from group 1 and in Figure 3-31 for group 2 samples.

In first place, it can be observed that all the samples don't evolve in the same way. Some of them accumulate positive deformations (expansion) due to thermal cycles while some of them contract. For instance, samples B1, B3 and CV4 are similar: an augmentation of deformation is observed through the imposition of thermal cycles. Whereas, samples B2, CV1, CH2 and CH5 presents a global decay response, even if sample B2 shows a positive tendency until cycle 250. All samples considered shows an erratic response in the first cycles, until cycle 70 in group 1 and until cycle 50 in group 2. This is evidenced mostly on samples B1 and B2, where a big accumulation of deformations is observed and in samples CV1, CH5 and CH2 where an important drop of vertical strains is observed. By eliminating the data obtained for this firsts cycles (70 for group 1 and 50 for group 2) table

3-10 presents the final accumulated deformations. It can be observed that samples drilled in horizontal direction (samples CH2 and CH5) show a higher amount of accumulated axial deformations than samples obtained in vertical direction (samples CV1 and CV4). Another remarkable observation is that, between 250 and 400 cycles an important value of irreversible deformations is observed, especially in samples B1, B2 and CH2. Moreover, after this interval the strains remain constant. Nonetheless, this can't be observed on sample B3 and CV4 response, in this case the deformation keeps a constant positive tendency since cycle 200 until the end.

In addition, sample CV4 presents a variation in the axial deformations amplitude, it decreases until 300 cycles where it becomes stable until the end of the experimental program.

Radial strains.

The evolution of radial deformations through the imposition of thermal cycles is presented in Figure 3-32 and Figure 3-33 for samples from group 1 and group 2 respectively. For this last group, the evolution of sample CH5 was not considered in the analysis because it presents a peculiar response that it will be describe afterwards. From a global point of view, it can be observed that except for sample CH2, (the only sample drilled in horizontal direction) all of them behave in the same way, radial deformations increase. Nevertheless, sample B2 presents a much higher value of final deformation than the other 4 samples (sample B1, B3, CV1 and CV4).

By considering the samples from group 2, where the direction of extraction in relation with the cliff's face is known, it can be remark that samples vertically drilled (sample CV1 and CV4) radial deformation presents a positive tendency with similar final strains accumulation value: 2.7e-4 m/m for sample CV1 and 3.2e-4 m/m for sample CV4. In the other hand, sample CH2 (horizontally drilled sample) radial strains decrease, just as its axial strains, but with a lower accumulated strain, only 4.6e-4 m/m in comparison to the 27e-4 m/m obtained in axial direction.

As it was observed on the axial strains the radial strains shown an erratic response until cycle 70 for samples obtained from blocks and until cycle 50 for the ones that were obtained from the cliff face. In this case it is mostly evidenced in samples B2 and CV1 by a considerable increment and in sample B3, by a strain decrease.

As mention earlier, the experimental programs for samples from blocks and drilled from the face of the cliff were not performed at the same time. Which means that the strain gages glued on samples of group 1 and samples of group 2 were not connected to the acquisition system at the same moment. In fact, it is only possible to connect 4 pairs of strain gages in the equipment used. Therefore, the erratic response observed with both groups of samples for the first cycles cannot be related to electrical punctual variations, but may be related to a stabilization time of the gauges.

After excluding the response observed in the first cycles the obtained accumulated radial strains are presented in Table 3-11, where, it can be observed that sample B1 and CV1 presents a radial accumulated strain very small in relation with the other values obtained, it can be assumed that this samples don't accumulate strains in this direction. Moreover, sample CH2 presenters a variation in the strains amplitude between 200 to 300 cycles (see Figure 3-33), after the 300th cycle, the response stabilizes and returns to its original amplitude value.

As mention before, sample CH5 was not considered for the analysis because its evolution shows a peculiar response, see Figure 3-34. The strains amplitude observed this sample is very high in comparison with the other radial deformation measurements observed and its daily response is very different from the behavior observed for the other radial measurements; it is more similar to a vertical direction behavior.

In addition, its strains amplitude increases radically in cycle 270, and then decreases, stabilizing only after the cycle 400, with a higher value than the one observed before cycle 270. These variations in amplitude are also observed in other samples response (see daily response), but in these other samples, the strains amplitude remains stable after cycle 300. These characteristics may be related to electrical issues, or to gauges installation differences. Despite different manipulations, it was not possible to purpose a reliable reason. Altogether, Samples drilled in vertical direction presents general positive strains tendencies in both directions. Conversely, samples obtained in horizontal direction presents a negative tendency. Note that even if results of radial strains for sample CH5 were rejected, their tendency is negative. Regarding the accumulated strains, samples CV1 and CV4 (vertical) shows similar values for axial and radial directions, otherwise, sample CH2 (horizontal) presents a major strains accumulation in axial direction.

By considering the response observed in axial and radial deformation for all the samples evaluated it can be observed that in general terms, the highest strain accumulation is observed in vertical direction. Except for sample B2 and CV4 where the absolute values of the cumulated deformations were very similar in both directions: 8.0e-4 mm/mm and 9.0e-4 mm/mm for vertical strains and radial strains respectively for sample B2, and 2.5e-4 mm/mm in both directions for sample CV4.

However, from results observed in axial and radial strains (Figure 3-30, Figure 3-31, Figure 3-32 and Figure 3-33) it is not possible to determine a specific behavior regarding the direction of the sample to the cliff's face, that may indicate the global direction for the samples from group 1 regarding to the face of the cliff. The strains evolution in time don't follow any identifiable response.

Daily response.

In order to, evaluate the response of the gauges during the heating-cooling process, as it is not possible to observe it from the general tendency because of the scale, the daily strains response is evaluated for some punctual days.

Here it is important to point out that the heating treatment used for samples from group 2 is slightly different from the one performed over samples of group 1: the interval time used to perform measures pass from 3 hours to 2 hours, this allows to perform 6.5 instead of 6 cycles per day.

For this reason, the experimental program for samples obtained from blocks (group 1) with a total of 948 thermal cycles, lasts 158 days. While samples drilled from the cliff face (group 2) were subjected to 900 cycles in 138 days. An important consequence of this variation is that, for example, thermal cycle number 100 belongs to testing day 16 for group 1, and to testing day 15 for group 2.

Deformations measured during 5 different days (0, 2, 10, 100 and last day) are compared, in the first place to evaluate if the heating/cooling cycles may induce a variation in the deformation response and in the second place to compare the direction of strains during the heating/cooling cycles.

Results blocks samples are presented in Figure 3 From these figures, it can be observed that strains augments during heating (sample expand) and decrease during cooling, just as expected. In fact, all the vertical responses show an evolution that follows closely the temperature variations. Nevertheless, in some radial evolution, a slight delay is observed.

It is interesting to note that even if the general tendency of strains (like radial deformation for sample CH2 and axial strains for sample B2), the daily response continues to be expansion for heating and contraction for cooling.

Samples obtained from blocks shows an interesting behavior, note that for radial deformation evolution (Figure 3-35b, Figure 3-36b, and Figure 3-37b) a slight strains decrease is registered in all the samples, when the temperature starts to increase after the 3 hours interval at 10ºC, that is used to perform the measures. In fact, this drop is registered between 12:00 to 13:00 and the begging of the heating ramp is registered in the face at 12:00.

It is also important is to mention, that deformations recorded for radial deformation in sample B2 (Figure 3-36b) present a peculiar behavior for days 2 and 10. In first place, the slightly drop registered in all the samples after the 3 hours interval is not observed. Moreover, it seems that the sample expands during cooling and contract during heating. This phenomenon is only observed for this sample. Nevertheless, remember that for the first 70 cycles (12 days) an erratic response is observed.

From Figure 3-38 it can be noted that sample CV1 shows a delay between radial stains and the temperature variations, additionally, the vertical strains amplitude presents a slight decrease through the chosen days.

In the other hand, from Figure 3-40b it can be observed that radial deformations for sample CH2 follow the temperature variations: there is no delay between thermal changes and deformations, but a slight increase in the strains amplitude is observed through the test duration. As observed from the general tendency for sample CH2 (Figure 3-33) between 200 cycles and 300 cycles the amplitude of strains decrease, in Figure 3-40b, the response obtained at day 40 (cycle 260-266. 5) is included, it can be observed that the amplitude obtained for day 40 is 5 times lower than the value obtained for day 100.

In the same way, sample CV4 shows a decrease in vertical strains amplitude, just as it was observed from the general tendency response (Figure 3-31). In fact, vertical strains amplitude observed for sample CV4, for one cycle of day 10, cycle 65 for instance, is 3 times higher than the value observed for cycle 390 (day 60). Nevertheless, the amplitude remains constant for days 0, 2 and 10; and afterward for days 60, 100, 138. Actually, in the general strains evolution (Figure 3-31) it can be observed that this change in amplitude is quite fast: it happens between cycle 200 (day 30) and cycle 300 (day 46).

Generally speaking, samples obtained from blocks shows similar strains amplitudes for radial and axial direction, and any evolution is observed through the thermal cycling. Nevertheless, samples drilled from the cliff's face presents a higher amplitude strains in the vertical direction. Furthermore, sample CV1 and CV4 show a change in the vertical strains amplitude in a short time, between day 30 (cycle 200) and day 46 (cycle 300), this characteristic is also observed in radial strains for sample CH2.

Although, these variations in amplitude between 200 to 300 cycles are observed in several samples response. This phenomenon is reversible only in sample CH2 in the radial response, in other words, the strains amplitude before and after the interval between 200 (day 30) to 300 (day 46) cycles are similar.

Cyclic response.

With the aim of evaluating the deformation response for one thermal cycle (10ºC-50ºC-10ºC), temperature-deformation graphs are considering, using 5 punctual cycles to be analyzed. The temperature used for this analysis is the one measured at the face of the samples and presented in Figure 3 In fact, the evolution observed in cycle 500 is as expected: expansion during heating. However, radial strains at cycle 100 and 900 for sample B2, have almost none evolution. Deformation register at maximal and minimal temperatures are almost the same: 2.03e-4 m/m and 2.01e-4 m/m respectively for cycle 100 and 1.03e-4 m/m and 1.04e-4 m/m for cycle 900.

It is interesting to point out, that even if the strains globally decrease the temperature-strains response shows an expansion of the sample during heating. As example, the response of vertical strains for sample CV1 (Figure 3-46a) can be observed. Even if, for each cycle the deformations due to thermal variations are positive, the evolution for cycle 100 is placed under the evolution of cycle 50 and so on. In other words, the initial strain for each cycle decrease with the increase of number of cycles. Similar response is observed for vertical and radial strains for sample CH2 (Figure 3-48), vertical strains for sample CH5 (Figure 3-49a)

Conversely, the opposite response is observed on sample CV1 for its radial deformation and on sample B1, B3 and CV4 for both radial and vertical deformation. That has a general positive strains tendency, and the initial strain for each cycle increase with the number of cycles.

Moreover, sample B2 in vertical direction has a mixed response related to its general strains evolution. From Figure 3-30 it is observed that radial strains tend to increase until cycle 200, but then, they decrease. Thus, the initial strains for each cycle tend to increase from cycle 1 to cycle 200. In other words, cycle 50 is placed above cycle 1 and so on. Nevertheless, the response for cycle 500 is placed under cycle 1 response.

Another remarkable effect that can be observed is that strains don't follow the same path for heating and cooling. It exists like a delay in the strains evolution. In fact, at the same temperature, deformations are higher when the temperature decreases (cooling) than when the temperature increases (heating). Thus thermal-strains evolution has an oval shape instead of a line. Furthermore, the amplitude of this oval is not the same for all the samples and all the directions.

If we consider that the maximal amplitude is located at the mean temperature (29.5ºC), the mean radial and vertical amplitude are presented in Table 3-12. Regarding to these amplitudes is important to realize that samples mostly composed by calcite presents a higher strains amplitude in radial sense. Whereas, all samples (vertically and horizontally drilled) composed by calcite and quartz show a larger amplitude in axial direction. The thermal expansion coefficient can be defined as the slope of a linear tendency adjusted to the temperature-strains results. Table 3-13 presents the values obtained for these coefficients at the cycles considered for the analysis for samples from group 1 and Table 3-14 for samples from group 2. In general, the values of thermal expansion found are similar to those found by the Laboratoire Régional des Ponts et Chaussées de Toulouse during their characterization of LaRG limestone, between 14.2e-6 ºC-1 and 1.4e-6 ºC-1 (see chapter 3.1.4). Nevertheless, the value obtained for sample B2 in radial direction is negative and smaller than the one found for the LaRG limestone. This condition may be associated to the response observed for this sample in this specific direction.

In the same way coefficients found for sample CV1 and CH5, are quite different from the values obtained during the characterization of LaRG. In fact, thermal expansion coefficient for radial strains in sample CV1 is 2.5 times lower than the minimal value formerly presented (1.4e-6 ºC-1). Moreover, the coefficient obtained for axial direction in samples CV1 and CH5 is 2.5 times higher than the maximal value obtained from the thermal expansion coefficient test.

Nonetheless, it is important to remember that the strains amplitudes observed in vertical direction for samples CV1 and CH5 are exceptionally high, in comparison with what was found for the other samples considered. Additionally, is important to remember that the tests performed by Laboratoire Régional des Ponts et Chaussées de Toulouse were realized on samples obtained from blocks fell inside the cavern, therefore the mineralogical composition of samples used should be similar to group 1 samples.

Moreover, from Table 3-13 and Table 3-14 it can be remarked that the thermal expansion coefficient remains constant, despite the general strains evolution observed through the thermal cycles, like the initial deformation at each cycle.

Radial

(ºC -1 ) 1 4.0e-6 3.0e-6 7.0e-6 -2.0e-6 3.0e-6 1.0e-6 2 4.0e-6 4.0e-6 7.0e-6 -2.0e-6 3.0e-6 2.0e-6 10 4.0e-6 4.0e-6 7.0e-6 -2.0e-6 3.0e-6 2.0e-6 100 3.0e-6 3.0e-6 7.0e-6 -0.9e-6 3.0e-6 3.0e-6 500 3.0e-6 4.0e-6 10.0e-6 3.0e-6 2.0e-6 3.0e-6 900 6.0e-6 5.0e-6 7.0e-6 0.2e-6 2.0e-6 2.0e-6 mean 4.0e-6 3.8e-6 11.3e-6 -0.9e-6 2.7e-6 2.2e-6

Vertical

(ºC -1 ) 1 40e-6 1.0e-6 10.e-6 3.0e-6 8.0e-6 2.0e-6 30e-6 2 40e-6 0.6e-6 10e-6 3.0e-6 8.0e-6 3.0e-6 30e-6 10 40e-6 0.4e-6 10e-6 2.0e-6 7.0e-6 2.0e-6 30e-6 100 40e-6 0.2e-6 10e-6 4.0e-6 6.0e-6 2.0e-6 30e-6 500 30e-6 0.2e-6 3e-6 2.0e-6 3.0e-6 3.0e-6 30e-6 900 30e-6 1.0e-6 3e-6 2.0e-6 6.0e-6 4.0e-6 40e-6 mean 36.7e-6 0.6e-6 7.7e-6 2.7e-6 6.3e-6 2.7e-6 31.7e-6

As mentioned before, samples from group 2 were drilled form the face of the cliff, considering two different directions (parallel and perpendicular to the cliff's face). It can be observed that these samples present higher thermal expansions in axial direction, which means that samples expand the most in perpendicular direction to the face cliff, which can be related to the direction where the cliff suffered a decompression during the valley formation geological process.

It can be remarked that samples B1 and B3 from group 1 present a similar thermal expansion coefficient in both directions and have a higher temperature-strains evolution amplitude for radial deformation (oval amplitude). Conversely, sample B2 has higher thermal expansion coefficient in axial direction and its temperature-strains evolution is similar for vertical and radial direction

These conditions may be related to a preferential deformation direction in the sample (anisotropy), perhaps related to its direction in relation to the cliff. However, as mention before, because of the nature of these samples, there is no information about its original position in the cliff and from the results obtained is not possible to infer their directions.

From Figure 3-43 to Figure 3-49 it may be assumed that there is no strains accumulation after each thermal cycle, which means that the deformation at the beginning and the end of each cycle are the same. Nevertheless, because the strains evolution observed thought the imposition of thermal cycles this is not the case.

Therefore, even if small a strains accumulation should be observed in each cycle, as was observed by calculating the difference between initial and final strains at each cycle. Figure Negative values of accumulated strains represent cycles were the final deformation is smaller than the initial one. This condition is observed when the cycle belongs to an interval where the global strains decrease. For instance, sample CH2 that present a negative tendency in the strains evolution for both directions (see Figure 3-31 and Figure 3-33) has negative values of cumulative strains for all the cycles considered, except for radial strains at cycle 2 and 200. However, despite axial deformation at sample B3 have a positive tendency the irreversible deformations after each cycle presents negative values.

No tendency is observed for the evolution of irreversible deformations; nevertheless, for all the samples and directions at 900 cycles, this value is very close to zero. Moreover, the variation of radial plastic strains through thermal cycles for samples from group 2, tends to decrease with the imposition of thermal cycles.

From the general tendency and the daily and cycle strains response is not possible to determine any feature regarding to the mineralogical composition or the original sample direction to the face cliff, that may indicate for instance the original direction of group 1 samples.

Elastic waves propagation velocities.

Elastic waves propagation velocities are measured at least once a week for 15 samples: 5 obtained from blocks, 5 drilled in the face of the cliff in vertical direction and 5 in horizontal direction. The samples used to measure deformations, except sample B1 that is only used for strains measurements, are included. Furthermore, this parameter is measured for all the cylindrical samples at least before and after the heating/cooling treatment, to have an extra characterization of the rock.

For the purpose to evaluate if strain gauges have any effect in the elastic wave propagation velocities, measurements were performed on samples B2, B3, CV1, CV4, CH2 and CH5, with and without the gauges strains glued; results obtained are presented in Table 3-15. It can be observed that similar values were found for both conditions. That means that any elastic wave propagation changes are induced by the strain gauges system. La Roque Gageac cliff is characterized by its heterogeneity (see chapter 2.3-Geology). This condition is highlight by the mineralogical compositions differences found between samples from inside the cavern and outside of it. Yet, this heterogeneity is also evidenced in samples from the same mineralogical composition. For instance, from Table 3-6 it can be observed that initial Vp value obtained for samples B2 and B3 are higher than the Vp value found for samples B5, B6 and B7.

The difference observed in the initial Vp value may explained, why after 950 thermal cycles almost any evolution of this parameter is observed for samples B1 and B2, while, samples B5, B6 and B7 experience a decay of 4% of the initial Vp value, for 828 cycles, as can be observed in Figure 3-54.

In the other hand, Figure 3-55 presents the evolution of compressive elastic wave velocity (Vp) for samples obtained from the face of the cliff, these samples have a lower initial value of Vp than block samples, as is presented in Table 3-6, they shows a decrease in the Vp value through the thermal cycles; after 910 cycles, this value decreases about 8% of its initial value. Samples oriented vertically present a mean decrease for the Vp value around 10%, while horizontal samples show a mean decay of 7%. That means that samples drilled vertically present more damage associated to thermal fatigue. From Figure 3-54 and Figure 3-55, it can be observed that the reduction observed in the compression wave velocities through thermal cycles is not constant for the whole experiment time. In fact, it can be observed that the Vp values of seem tend to a stable value, at the end of the experimental program, after cycle 600. Moreover, it can be observed that the highest Vp evolution is registered in the firsts cycles, samples CV2, CV3 and CV5 evolves the most between the cycle 200 and 300 while samples B5, B6 and B7 presents a steeped evolution for the first 100 cycles.

It is important to remember that, as mentioned in chapter 3.2 (samples description), samples obtained in horizontal sense may belong to the same rock core. This is the case of samples CH1, CH2 and CH3 and samples CH4 and CH5. Moreover, sample CH3 shows a higher Vp reduction than the other samples that composed the rock core 10, after 900 cycles the Vp value for sample CH3 decreases a 9 while sample CH2 and CH1 decreases only a 7%. It is also remarkable that, even if samples CH2 and CH5 were originally located deeper than samples CH1 and CH4 respectively, their initial Vp values and evolution through thermal cycles are very similar. This may indicate that the effect of climatic solicitations in the cliff is mainly evidenced in the first 25 cm.

In the same way, it is interesting to observe that despite vertical samples were all obtained in the first 20 cm of the cliff, they show more heterogeneity in the values of initial Vp, with a range of 685 m/s in comparison with the range of 226 m/s observed for horizontal samples, (excluding sample CH3).

Furthermore, horizontally drilled samples show more heterogeneity in the amount of damage induced by the thermal cycles, evidenced by the decrease in the Vp value (between 5% to 14% of the initial measurement), as observed with the horizontal samples. This drop of wave propagation velocity is stronger in samples with higher initial velocity.

Figure 3-57 presents the evolution of shear elastic wave propagation velocity (Vs) for group 1 samples and Figure 3-58 for samples from group 2. The Vs value is calculated as the mean Vs value from the measures realized in two perpendicular directions. From these figures, it can be remarked that the Vs evolution is similar to the one found for compressive waves, samples from group 1 presents almost any evolution of Vs with the thermal cycles, samples B1 and B2 show a decay of only 1%, and for samples B5, B6 and B7 this parameter reduces a 3% from the initial value.

Moreover, samples from group 2 show an important heterogeneity in the value of initial Vs value in samples vertically drilled, even if all these samples are obtained very close to the cliff face. Samples with a higher decrease in the Vs value are samples CV2, CV4 and CV5 (mean decay of 12% of their initial value). Moreover, sample CV1, the one that has a higher initial value shows the lowest evolution through heating/cooling cycles, only 4% after 900 cycles.

In the other hand, as was observed for the compressive waves evolution, sample CH3 shows the higher initial value of the horizontally drilled samples. It is also the sample were the Vs value evolves the less, only 4% in comparison with the 7% observed for the other horizontal samples.

From Figure 3-58b, it can be observed that initial Vs values and their evolution for samples CH1, CH2 and CH3, follow a response that can be related to the position of the samples in the cliff. Note that Vs for sample CH3 is higher than for sample CH2 and this last value is higher than the measurement obtained in sample CH1. Nevertheless, for samples CH4 and CH5 this cannot be observed.

The measurements of elastic wave propagation can be used to directly calculate the dynamic Young modulus and Poisson ratio, through the equations 3-4 and 3-5 Finally, the dynamic bulk modulus can be calculated by the equation 3-9, where Ed and d are the dynamic Young modulus and Poisson ratio, respectively.

K d = E d 3(1 -2 υ d )
Eq. 3-9

The bulk modulus evolution through the imposition of thermal cycles can be observed in Figure 3-59 for samples obtained from blocks and in Figure 3-60 for samples directly drilled from the face of the cliff.

It can be noted the bulk modulus evolves very similar to the values of VP and Vs. Samples B5, B6 and B7 show a mean decrease of 6% from the initial value, for samples from group 2 a mean decreases of 20% of the initial stiffness value is observed. However, samples B2 and B3 experience a mean stiffness value decrease of only 2% of the initial value. The stiffness variations between the initial values and the measurements after the thermal treatment are detailed in Table 3-16. The mean UCS value found for samples submitted to 828 thermal cycles is 11% lower than the mean value obtained for samples without thermal treatment. In fact, considering a linear model the variation of Uniaxial compressive strength with the number of cycles follows the equation 3-10, σ c (MPa) = 0.0075 * Number of cycles + 45.044

Eq. 3-10

From the stress-strain evolution, the Young modulus can be calculated, as the slope of the virgin loading path. The ASTM standard stated that this modulus can be obtained by any of the next methods:

a) The tangent slope at a stress level, which is fixed as a percentage of the maximum strength. b) The tangent average slope of the strain line portion of the stress-strain curve. c) The secant slope, from zero stress to some percentage of the maximum strength.

Moreover, as mentioned before for the uniaxial compression test performed, two unloading/ loading cycles were considered at 50% and 80% of the Uniaxial compressive strength estimated in the characterization of LaRG limestone (47.5MPa). The Young modulus can be obtained also as the slope of the unloading paths.

In order to evaluate the effect of thermal cycles in the rock Young modulus, this value was calculated from the virgin stress-strain paths as the average slope of the stress-strain portion (Etan) from 1/3 c and 2/3c, and as the secant modulus (Esec) from zero to 70% c. But also, as the slope from the unloading path at 50% (E50) and at 80% (E50). The obtained young module values are listed in Table 3-18. It can be observed that the Young modulus value obtained as the slope of the linear virgin stress-strain is a 13% higher than the one obtained as the secant slope. Moreover, the Young modulus observed in the unloading cycles are higher than, that at the virgin loading, and very similar between them. Similar results are observed in 99% calcite limestones by [START_REF] Al-Shayea | Effects of testing methods and conditions on the elastic properties of limestone rock[END_REF], where is also stated that the Young modulus value has less variation when is obtained in the unloading paths than that obtained from the loading virgin curve, and it is almost constant regardless of the stress level.

Figure 3-63 presents the evolution of Young modulus through the imposition of thermal cycles, as was observed whit the Uniaxial compressive strength, this modulus shows a reduction with the increment of thermal cycles. The reduction observed is very similar for all the methods considered to determine the young modulus, in general, terms after 900 cycles this value decays in a 10%. Therefore, to be able to obtain the static young modulus from elastic wave propagation velocities, a dimensionless parameter K= Ed/Est should be included in the equation Eq. 3-4.

The new equation to calculate the static young modulus from the elastic wave propagation velocities is defined by the equation

E st = 1 K ρ V s 2 (3V p 2 -4V s 2 ) (V p 2 -V s 2 )
Eq. 3-11

The ideal K value will be 1, which means that the static young modulus is perfectly calculated from the equation Eq. 3-1. However, in the literature, different values of K are found. For instance, [START_REF] Ciccotti | Differences between static and dynamic elastic moduli of a typical seismogenic rock[END_REF] found for a mudstone-limestone with a calcite content of 50% that the K value varies from 1 to 1.3.

In the same way, [START_REF] Martínez-Martínez | Comparison of the static and dynamic elastic modulus in carbonate rocks[END_REF] found for limestones composed of 90% of calcite that K varies from 0.5 to 1.9. It is important to mention that for [START_REF] Martínez-Martínez | Comparison of the static and dynamic elastic modulus in carbonate rocks[END_REF] study the static young modulus was obtained as the secant slope to the 70% of the maximal strength. Moreover, in the study performed by [START_REF] Najibi | Empirical relations between strength and static and dynamic elastic properties of Asmari and Sarvak limestones, two main oil reservoirs in Iran[END_REF] on a similar limestone, a higher range for K value was found, between 1.4 to 3.8.

The mean K value obtained by considering 3 methods to calculate the static Young modulus is presented in Figure 3-64, the parameter K is higher when considering the secant method for the Young modulus calculation. If the static Young modulus is defined as the slope of the unloading path, which is the more consistent method (with less variations, as mentioned before) the parameter K varies from 1.9 to 2.5.

Discussion on the mechanical response.

In addition to the measurements of strains, elastic wave propagation velocity and Uniaxial compressive strength, the effect of thermal cycles over LaRG limestone porosity was evaluated through mercury intrusion porosity test performed on samples that were submitted to a different amount of cycles: 0, 216, 420 and 630. However, any remarkable response is obtained from this test, the obtained results are presented in appendix A.

After more than 1000 thermal cycles between 10ºC to 50ºC, no visible cracks were observed on the samples. Nevertheless the decrease in the Vp suggest a thermal fissuration of the samples, [START_REF] Martínez-Martínez | Spatial attenuation: The most sensitive ultrasonic parameter for detecting petrographic features and decay processes in carbonate rocks[END_REF][START_REF] Brotóns | Temperature influence on the physical and mechanical properties of a porous rock: San Julian's calcarenite[END_REF][START_REF] Hudson | Influence of thermal damage on linear and nonlinear acoustic properties of granite[END_REF]) related to the propagation or generation of internal microfissures and voids, that generally appear in the boundary of the grains [START_REF] Malaga-Starzec | Microscopic and macroscopic characterization of the porosity of marble as a function of temperature and impregnation[END_REF][START_REF] Yavuz | Estimating the index properties of deteriorated carbonate rocks due to freeze-thaw and thermal shock weathering[END_REF]. In fact, P-wave velocity has shown to be very sensitive to quantify and evaluate the thermal damage on rocks [START_REF] Yavuz | Estimating the index properties of deteriorated carbonate rocks due to freeze-thaw and thermal shock weathering[END_REF][START_REF] Martínez | Influencia De La Alteración Sobre Las Propiedades Mecánicas De Calizas, Dolomías Y Mármoles[END_REF][START_REF] Wanne | Bonded-particle modeling of thermally fractured granite[END_REF]Luque et al., 2011).

In this experimental program, the elastic waves propagation was always measured at least at the beginning and the end of the thermal treatment to have an extra characterization parameter, in addition, the porosity of the samples was also calculated. Table 3-6 presents all the characteristics of the 30 cylindrical samples considered, including the Vp and porosity.

It is interesting to point out that samples obtained from blocks (group 1) present a lower porosity and high elastic wave propagation velocity, a relation that may be expected as the measures of Vp is highly influenced by the porosity of the rock. Conversely, samples obtained from the face of the cliff (group 2) shows a higher porosity coupled with lower values of elastic waves propagation velocities. Nevertheless, samples from group 1 present a linear relation between porosity and Vp, whereas, the VP-porosity relation observed for group 2 samples has a more disperse relationship, as illustrated in Figure 3-65. It is important to note that horizontally drilled samples, present more porosity the closer they were placed to the face of the cliff, as presented in chapter 3.5.2, the original location in relation to the cliff face has an impact on the evolution of Vp through thermal cycles, samples placed deeper (lower porosity) present more important decrease of their Vp.

Another key point is that sample B7 shows a similar initial porosity than sample CV1 and CH2. Nevertheless, the reduction in Vp and Vs, and therefore in bulk modulus, is less important in sample B7 than in samples CV1 and CH2. Which may be related to the mineralogical composition of the samples, remember that sample B7 is mainly composed by calcite, while samples CV1 and CH2 are composed by calcite and quartz. The thermal response between quartz and calcite minerals may induce some additional stresses on the grains boundaries that can be the responsible for the presented response [START_REF] Hall | Thermal fatigue and thermal shock in bedrock: An attempt to unravel the geomorphic processes and products[END_REF][START_REF] Hartlieb | Thermo-physical properties of selected hard rocks and their relation to microwave-assisted comminution[END_REF]. In fact, other experimental analysis [START_REF] Malaga-Starzec | Microscopic and macroscopic characterization of the porosity of marble as a function of temperature and impregnation[END_REF][START_REF] Tiskatine | Experimental evaluation of thermo-mechanical performances of candidate rocks for use in high temperature thermal storage[END_REF] found that the presence of quartz may control the response of the rock under thermal cycling due to its large and variable thermal expansion coefficient, even below 650ºC, where the mineral suffers important changes.

As mention before, for samples used for deformation measurements (except sample B1) the evolution of the elastic wave propagation velocities is also evaluated. The deformation and Vp evolution are compared for one sample of each type: sample B3 obtained from blocks, sample CV4 drilled in the face of the cliff in vertical sense and sample CH2 drilled from the face of the cliff in horizontal sense. Furthermore, between cycles 800 to the end of the test any Vp measure is performed and the strains evolution has less dispersion. Nonetheless, from Figure 3-30 and Figure 3-32, it can be observed that this variation in the strains evolution is also observed in sample B1 where any measure of Vp is realized.

In the same way, Figure 3-68 shows the strains and Vp evolution for sample CV4. In first place, it can be remarked that vertical and radial strains have similar evolutions. An important reduction in Vp for the first 300 cycles is observed, for this period the radial and vertical strains present a higher accumulation. Additionally, vertical strains present a variation in the strains amplitude, after cycle 300 it decreases. After cycle 300, thermal cycles accumulate less deformations. From Figure 3-52 and Figure 3-53 it can be observed that plastic strains tend to zero, especially for radial strains, and, Vp tend to a value of 2690 m/s.

Finally, Figure 3-69 illustrates the strains and Vp evolution for sample CH2, drilled in the cliff in horizontal direction, for this case it can be observed that vertical and radial strains have different responses, even if in both cases their tendency is negative (compression), the final accumulated vertical deformation is 6 times the radial.

The radial strains evolution for sample CH2 is more important between 200 to 600 cycles, where the Vp begins to have a less steep response. After cycle 600 Vp tends to 2830 m/s. In the same way, from cycle 600 to the end of the test, radial and vertical strains remain nearly constant.

Generation and propagation of fissures.

From the mechanical results, it is observed that samples submitted to thermal cycles experience damage that can be identified by a drop in P-wave propagation velocity and reduction in the strength. This damage can be related to the generation or propagation of fissures developed generally in the grains boundaries [START_REF] Goudie | The thermal degradation of marble[END_REF]Luque et al., 2011;[START_REF] Hayatdavoudi | Effect of cyclic temperature, water vapor and exposure time on micro fracture propagation in shale[END_REF] Moreover, during a former experimental program developed on samples from LaRG, some important macro-fissures were observed (Gasc-Babier, Girmaand and Gendre, 2014). Therefore, to evaluate the generation and propagation of fissures (at photo-scale) due to thermal damage, photos were taken on samples (see Table 3-9) submitted to thermal cycles. These photos were analyzed with the digital image correlation technique through the software N-corr [START_REF] Blaber | Ncorr : Open-Source 2D Digital Image Correlation Matlab Software Ncorr : Open-Source 2D Digital Image Correlation Matlab Software[END_REF]. Results obtained with this analysis are presented in this section.

Brazilian test.

The digital image correlation (DIC) technique is widely used to analyzed samples submitted to tensile strength test, like the Brazilian test. Therefore, to evaluate the results obtained with the N-corr software photos were taken to samples 1343 to 1348 (see Table 3-3) during the loading path of Brazilian test.

Moreover, a common practice for DIC analysis is to generate a speckle pattern in the samples with white and black paint to enhance the natural contrast from the rock texture [START_REF] Hild | Digital image correlation: from displacement measurement to identification of elastic properties-a review[END_REF][START_REF] Tung | Application of digital-image-correlation techniques in analysing cracked cylindrical pipes[END_REF][START_REF] Nguyen | Fracture mechanisms in soft rock: Identification and quantification of evolving displacement discontinuities by extended digital image correlation[END_REF]. For this reason, half of the samples considered were painted and the other half was left with the natural contrast, to evaluate if different responses are obtained with and without the pattern. Moreover, the response obtained with this technique agree with the experimental observations: the left side of the sample moves to the left (negative displacements values) while the right side moves to the right (positive displacements value), the division of these two parts of the samples is delimited by the apparition of the crack. In addition, sample 1345 presents more displacements in the upper part of the sample, while, sample 1346 has more displacement in the lower part. In fact, the crack observed in sample 1345 propagates from the top of the sample, while in sample 1346 it propagates from the bottom, this characteristic is also observed in the horizontal strains field.

Even if the strains values obtained from the DIC analysis are quite high, almost 30%, they show a clear evolution of the crack appearance and aperture, even when the fissure can't be visually identified, see photo 3 on Figure 3-72. The strains field results are helpful to identify the appearance of cracks in the sample and to evaluate the aperture and propagation of them, as presented by [START_REF] Morgan | An experimental and numerical study on the fracturing processes in Opalinus shale[END_REF]. The imposition of a speckle pattern in the samples don't necessarily help to have better results in the digital image correlation analysis on samples from LaRG. This can be explained by the fact that this rock is composed by minerals with different colors and sizes (see section 3.3.1), therefore it has already good contrast.

Nevertheless, one face of samples C1, C2 and C3 were painted to have this pattern, because these samples were used to evaluate the generation of fissures associated with thermal cycles that may be thinner than those observed in the Brazilian tests, thus an additional contrast may help in the analysis. In the other hand, samples B20 and C4 were not painted since there were used to evaluate the propagation of fissures and it was observed that the natural contrast of the samples is enough to evaluate the propagation of a crack visually identifiable.

Fissures generation.

To evaluate the appearance of cracks due to thermal cycles, 3 prismatic samples (C1, C2 and C3, see Table 34567) were photographed at least once a week and evaluated with the DIC technique, as detailed in section 3.4.5. These samples were subjected to 948 cycles of temperature.

Any new fissure was detected on the surface of samples C1, C2 and C3, at least from a visual examination, as was observed on the cylindrical samples (presented before). Nevertheless, DIC analysis was performed to evaluate if any pixel-scale fissures can be observed. Is important to mention that photos considered present a scale of 0.08mm/pixel. Note that no identifiable crack can be observed from these analyses, similar results are observed by [START_REF] Leiss | Fabric anisotropy and its influence on physical weathering of different types of Carrara marbles[END_REF] on different types of marbles, subjected to heating until 50ºC.

The strains values observed are very small in comparison with what was observed in Brazilian tests. The small variations observed in the strains fields can be related to slight variations in the illumination imposed to the samples. Those variations may generate that some pixels are darker than others. Nevertheless, it is interesting to observe that sample with the speckle pattern shows less noise associated to changes in light. Note that the strains obtained from the DIC analysis present less dispersion and have a smaller range for the face with the speckle pattern (Figure 3

Fissures propagation.

Failure of rocks due to thermal damage associated with a process of thermal fatigue generally occurs where a pre-existing weakness is located, such as minerals boundaries, or pre-existing fractures. In fact, when the rock matrix present cracks, the thermal cycles will generate preferable the propagation of these fractures [START_REF] Hall | The role of thermal stress fatigue in the breakdown of rock in cold regions[END_REF][START_REF] Gómez-Heras | Surface temperature differences between minerals in crystalline rocks: Implications for granular disaggregation of granites through thermal fatigue[END_REF][START_REF] Walsh | Micromechanical modeling of thermal spallation in granitic rock[END_REF][START_REF] Hall | Thermal fatigue and thermal shock in bedrock: An attempt to unravel the geomorphic processes and products[END_REF] With the purpose of evaluating the propagation of cracks due to thermal cycling, two types of samples (with an artificial fissure and with a natural one) subjected to thermal treatment are photographed at least once a week. These images are also analyzed with the DIC technique.

Artificial fissured sample.

A crack was generated in cylindrical sample B20, thought a hammer impact, the fissure obtained is presented in Figure 3-76. This sample was submitted to a total of 396 cycles between 10ºC to 50ºC. Based only on visual inspection, a clear evolution of the crack can be observed during the test. Figure 3-77 presents the evolution of the crack for three specific moments: 0, 150 and 396 cycles. It can be observed that the thermal cycles induce a propagation of the crack but also a widening of it. Furthermore, Figure 3-78 presents the results obtained from the DIC analysis performed on samples B20, first it can be observed that the maximal and minimal strains values obtained are higher than those observed in samples C1, C2 and C3, but much lower than the strains obtained for the analysis with the Brazilian test. Moreover, from this analysis it is possible to identify the location and response of the fissure, vertical strains (Eyy) shows that the crack gets wider with the imposition of cycles of temperature, additionally, horizontal strains (Exx) concentrate above and below the tip of the main crack. The results of Eyy and Exx suggest that at the tip of the crack the fissure opens. From the combination of the strains, it can be concluded that the main strain presents a direction that is perpendicular to the crack. In the same way, shear strains (Exy) also concentrate at the tip of the main crack, showing the shear effect associated with the propagation of the crack. Figure 3-77 shows that between the initial condition and cycle 150 a small fissure located at the top of the sample propagates. In fact, it seems that this fissure could reach the main crack, (Figure 3-76), delimiting a block. These phenomena can also be observed with the DIC analysis (Figure 3-78). Note that from cycle 150 a concentration of strains is observed at this location, delimiting the mentioned block, from vertical and shear field strains from cycle 280 and over.

Figure 3-78 DIC analysis for sample B20

Natural fissured sample.

A small block that has a natural fissure was collected from the troglodyte cavern located in the cliff of La Roque Gageac during a technical visit to the site. This block was prepared to obtained flat surfaces, more favorable for the DIC technique, obtaining a prismatic sample with a crack (samples C4), see Figure 3-79. This sample was submitted to 1046 thermal cycles, with a rate of 6.5 cycles per day, and photos were taken once a week on two opposite faces of the sample. Even if from visual inspections was not possible to determine the evolution of the cracks observed in sample C4, a DIC analysis was performed over the two faces of the sample. Despite the precautions taken to assure that the sample was placed in the same position, between test day 45 (cycles 290 to 296.5) and 52 (cycles 322 to 328.5), the camera slightly moves, leading that photos taken before and after this days were not comparable. Therefore, the DIC analysis should be performed in two sets: from the begging of the test until test day 45 (Cycle 290) and from test day 52 (cycle 325) until the end. Results obtained from the DIC analysis are presented below.

Sample C4 -face 1.

The evolution of the main crack observed in samples C4 is followed by photos taken of face 1.

Results of the DIC analysis are presented in Figure 3-82 for cycles 0 to 290 and Figure 3-83 for cycles 325 to 1047. Note that for cycles 0 to 293 horizontal (Exx) strains are concentrated along the fissure. It is interesting to point out that the value of strains is negative, suggesting a closure of the crack associated with the thermal cycles, which is a completely opposite response from what is observed in sample B20 (Figure 3-82). Moreover, shear strains (Exy) shows a concentration with positive values along the fissure but mostly in the lower part of the fissure, where the crack is composed by two different paths. From vertical strains (Eyy) not clear evolution can be observed for the complete fissure, nevertheless, it can be observed that Eyy present a higher value at the tip of the crack. This characteristic is also observed in Exx and Exy, indicating a propagation of the fissure.

From cycles 325 to 1047, a less evident evolution of the complete crack can be observed. The strains fields obtained here are related to the condition at cycle 325. Nevertheless, it can be observed that all strains fields considered (Exx, Eyy and Exy) show higher strains values where the fissure is composed by two paths, just as observed from cycle 0 to 293.

Furthermore, from this last photo set it can be observed that, from cycle 820, a point of high strains is observed at more or less the middle of the crack length. This strains accumulation is associated to the felt of a grain, as can be observed on Figure 3-81, that focus on this point at cycle 730 and cycle 820. This observation reinforces the idea that thermal damage induces fissures along the grain boundaries. The second face considered in the analysis of sample C4, is the one that presented the detachment of the small block. As for face 1, the DIC analysis was performed in two stages: the first one goes from cycle 0 to cycle 290, and the second one goes from cycle 425 to 975, that is the last moment registered before the complete detachment of the block, the DIC analysis of this last group is related to the sample condition at cycle 325.

Figure 3-84 and Figure 3-85 present the DIC analysis obtained for the two sets of photos, it can be observed that at cycle 975 a high strains concentration that clearly defines the block detached is observed, even before it really fell. Moreover, it is interesting to observe that, even if from visual observations it is not possible to identify a crack in this zone, the results of DIC analysis, see Figure 3-84, show a strains concentration at the block level since cycle 240.

Since Cycle 390 (see Figure 3-85) in addition to the strains that delimit the block a second strains concentration is observed, that may be related to the appearance of a fissure, however, from visual inspections is not possible to identify a crack at this place.

Experimental work conclusions.

An experimental program has been performed to evaluate the possible damage induced by thermal cycles. Samples from La Roque Gageac (LaRG) cliff were submitted to cycles similar to the atmospheric range of temperature, from 10ºC to 50ºC. Their response is evaluated through different laboratory measurements.

The results obtained from this experimental work suggest that atmospheric thermal cycles can induce damage in LaRG limestone, evidenced by strains accumulation, reduction in the elastic waves propagation velocities and strength parameters, and the propagation of fissures.

Some specific observations are summarized below:

• As a consequence of the heterogeneity of the cliff, the composition of the samples is related to the location of extraction. Specimens from blocks felt within the cavern are mainly composed by calcite while samples drilled from the cliff's face by calcite and quartz.

Additionally, the average values of porosity are 14% and 21% respectively. The exposition of the material to the atmospheric conditions can be the reason for the above characteristics, as the calcite can be more easily dissolve by water generating a lower calcite content and higher porosities in the cliffs face.

• Samples obtained from the face of the cliff have an anisotropic thermal expansion. Higher values coincide with the perpendicular direction to the face cliff, which can be related to t the direction where the cliff suffered a decompression during the valley formation geological process

• The thermal expansion remains constant thought the thermal cycles; therefore, the slope and amplitude of the strain-temperature relationship is not affected by the thermal damage (chapter 3.5.1). In the other hand, it is observed that thermal cycles generate a small plastic strains.

• The initial Vp values are also affected by the heterogeneity of the rock. The samples obtained from blocks show higher Vp values than the samples obtained from the face of the cliff. This can be linked to the mineralogical composition and porosity but also to the initial alteration level of the sample due to the exposition of the material to the atmospheric conditions of the cliff. In fact, for samples taken from the same borehole (CH1, CH2 and CH3) the highest Vp value is obtained for the deepest sample (CH3). Moreover, samples CH1 and CH2 shows similar values of Vp, which suggest that the effect of atmospheric conditions affect mostly the first 25 cm of the cliff.

• The evolution of elastic waves propagation velocities through the accumulation of thermal cycles was follow for samples from both locations. All samples present a decrease of this parameter, especially for the P-wave. This response can be explained by the generation and propagation of microcracks [START_REF] Yavuz | Effect of freeze-thaw and thermal shock weathering on the physical and mechanical properties of an andesite stone[END_REF]. Samples obtained from blocks that presents higher values of initial Vp, and therefore higher stiffness, shows the lowest decrease of Vp.

• Another evidence that the natural alteration is limited to the first centimeters of the cliff face is that the higher reduction of Vp is found in the deepest sample, indicating that this material was probably more isolated from the atmospheric actions.

• By comparing the Vp values evolution of samples with the same porosity but different mineralogical composition it is found that samples compose by calcite and quartz develop more thermal damage, indicated by a larger Vp decrease.

• One of the main consequences of the thermal damage is the change of the mechanical response of the rock. After 800 cycles, the Uniaxial compressive strength (UCS) and Young`s modulus present a reduction of 10% in comparison with samples with no thermal treatment.

• A potential law can be used to describe the relation between UCS and Vp for samples mainly composed by calcite, indicating a direct relationship between the reduction in Vp (thermal damage) and the weakening of the rock (reduction in UCS).

• The Digital Image Correlation (DIC) technique has shown to be a helpful tool to evaluate the generation and propagation of cracks.

• Thermal cycles at atmospheric range are not able to generate new cracks in LaRG limestone.

Results of DIC analysis on prismatic samples without fissures (C1, C2 and C3) don't shows any new crack. In the same way, it is not possible to determine any fissure in the 30 cylindrical samples tested, from visual inspections. Nevertheless, the results obtained in the mechanical response of the rock suggest the generation or propagation of existent micro-fissures.

• DIC analysis on naturally and artificially fissured samples indicate that thermal cycles induce the propagation of existent defects, like cracks, especially when they have no suffer formerly any thermal cycle (artificial crack). Note that the evolution of the artificial crack was observed even from visual inspections, whereas the main natural crack evolves very little. Moreover, from sample C4 it can be observed that thermal cycles induce the detachment of crystals across their boundaries which are weakness zones

• The prismatic sample with a natural crack present the detachment of a small block after 970 cycles, even if any fissure was observed at the beginning of the test program at this level.

Based on the results obtained on the other prismatic samples, it can be assumed that when this sample was collected from the site, a local defect existed at this corner and the thermal cycles generated the propagation of it until the detachment of the small block.

Introduction.

From a general point of view, rocks are composed by several minerals with different mechanical and thermal properties. These local heterogeneities result in the development of differential grain expansion during heating, which, in turn, is the source for the build-up of internal stresses.

The same effect can develop in single mineral material as a result of heterogeneity in grain or crystal shape and orientation (Le [START_REF] Per | Désagrégation des parios rocheuses et climat: approche thermique et thermodynamique[END_REF][START_REF] Wanne | Bonded-particle modeling of thermally fractured granite[END_REF][START_REF] Hall | Thermal fatigue and thermal shock in bedrock: An attempt to unravel the geomorphic processes and products[END_REF]. et al., 2006) In the same perspective, Royer-Carfagni (1999) found, by microstructural analyses, that thermally induced stresses are often concentrated at grain interfaces and not inside the grains, leading primarily to grain decohesion rather than fracturation. The author indicated also that microcracking effects tends to disappear when the thermal expansion coefficient of the different minerals present in the rock are similar. Other authors [START_REF] Wong | Microcracking and grain size effect in Yuen Long marbles[END_REF][START_REF] Zimmermann | Damage evolution during microcracking of brittle solids[END_REF] pointed out the grain size distribution of the rock as an important factor controlling microcracks generation and propagation. On the other hand, several authors [START_REF] Bodin | Continuum damage approach to asphalt concrete fatigue modeling[END_REF][START_REF] Wanne | Bonded-particle modeling of thermally fractured granite[END_REF][START_REF] Malésys | A probabilistic model to predict the formation and propagation of crack networks in thermal fatigue[END_REF][START_REF] Dubois | Finite element model for crack growth process in concrete bituminous[END_REF][START_REF] Lu | Microcrack-based coupled damage and flow modeling of fracturing evolution in permeable brittle rocks[END_REF][START_REF] Walsh | Micromechanical modeling of thermal spallation in granitic rock[END_REF] emphasize the time dependent character of fissure propagation.

Microcracking is accompanied by a progressive damage of the rock that may affect significantly its macroscopic mechanical properties [START_REF] Mazars | From damage to fracture mechanics and conversely: A combined approach[END_REF][START_REF] Lu | Microcrack-based coupled damage and flow modeling of fracturing evolution in permeable brittle rocks[END_REF]. As a result, after some time, microcracks may coalesce and eventually culminate in the development of a failure surface across the rock sample, leading to failure (Hall 2004b;[START_REF] Hall | The role of thermal stress fatigue in the breakdown of rock in cold regions[END_REF][START_REF] Tiskatine | Experimental evaluation of thermo-mechanical performances of candidate rocks for use in high temperature thermal storage[END_REF].

The modeling of such process may be indistinctively tackled within the continuum, discrete or hybrid frameworks available in Rock Mechanics. The choice of the method depends mainly on the characteristics of the problem under study [START_REF] Jing | Numerical methods in rock mechanics[END_REF], and, particularly, on the scale focused for the analysis.

As far as concerns thermal damage, two types of approaches are generally considered. The first one focuses on the modeling at the microscale and represents explicitly the grain responses within idealized microstructural patterns [START_REF] Lu | Microcrack-based coupled damage and flow modeling of fracturing evolution in permeable brittle rocks[END_REF]. These models are generally developed using the distinct elements method and represent the rock as an assembly of particles tight by separable joints. When the model is fed with experimental microscopic data on material heterogeneity and mineralogical arrangement, results have proven to compare well with laboratory observations. They generally show that fissures develop mostly at grains contacts [START_REF] Royer-Carfagni | On the thermal degradation of marble[END_REF][START_REF] Zimmermann | Damage evolution during microcracking of brittle solids[END_REF][START_REF] Wanne | Bonded-particle modeling of thermally fractured granite[END_REF][START_REF] Shushakova | Microcracking in calcite and dolomite marble: Microstructural influences and effects on properties[END_REF][START_REF] Walsh | Micromechanical modeling of thermal spallation in granitic rock[END_REF]Hoek and Martin, 2014).

On the other hand, the problem can be tackled at macroscale by idealizing the rock as a continuum provided with degradation laws for its mechanical properties. These laws aim at representing the macroscopic effect of grain decohesion, fissuration and, more generally, all types of microstructure degradation [START_REF] Lu | Microcrack-based coupled damage and flow modeling of fracturing evolution in permeable brittle rocks[END_REF].

This approach has been widely used for the modeling of microcracking phenomena in quasibrittle materials, like, rocks [START_REF] Mazars | From damage to fracture mechanics and conversely: A combined approach[END_REF][START_REF] Shao | Microcrack-based continuous damage model for brittle geomaterials[END_REF][START_REF] Lu | Microcrack-based coupled damage and flow modeling of fracturing evolution in permeable brittle rocks[END_REF][START_REF] Chang | A combined continuous-discontinuous approach for failure process of quasi-brittle materials[END_REF] or concrete [START_REF] Mazars | Continuum Damage Theory-Application to Concrete[END_REF][START_REF] Jirásek | Nonlocal models for damage and fracture: Comparison of approaches[END_REF][START_REF] Bodin | Continuum damage approach to asphalt concrete fatigue modeling[END_REF][START_REF] Dubois | Finite element model for crack growth process in concrete bituminous[END_REF][START_REF] Rojas-Solano | Interaction-based non-local damage model for failure in quasi-brittle materials[END_REF]. These materials share common characteristics like, heterogeneity due to the presence of different components and reduction in stiffness and strength due to microcracks development. Damage theory is the main tool used for this type of analysis. It consists in considering an irreversible degradation of stiffness parameters according to the evolution of the density of microcracks. The variable controlling damage evolution varies from one formulation to another, and can be driven by the strain rate [START_REF] Bodin | Continuum damage approach to asphalt concrete fatigue modeling[END_REF], the strain [START_REF] Mazars | Continuum Damage Theory-Application to Concrete[END_REF][START_REF] Lu | Microcrack-based coupled damage and flow modeling of fracturing evolution in permeable brittle rocks[END_REF], the stored energy [START_REF] Carol | On the formulation of anisotropic elastic degradation. I. Theory based on a pseudo-logarithmic damage tensor rate[END_REF][START_REF] Chang | A combined continuous-discontinuous approach for failure process of quasi-brittle materials[END_REF] or microcrack propagation characteristic [START_REF] Shao | Microcrack-based continuous damage model for brittle geomaterials[END_REF].

The experimental results presented in chapter 3 have evidenced that 'La Roque Gageac' limestone is susceptible to experienced mechanical damage under thermal cycles whose amplitudes are in the range of daily atmospheric variations. This chapter is devoted to the presentation of a numerical model able to reproduce this phenomenon. The macroscale approach is the one considered for thermal damage in rocks presented in this chapter.

Model scope.

The numerical model considered is based on the constitutive law initially developed for soft argillaceous rocks by [START_REF] Vaunat | Bond degradation and irreversible strains in soft argillaceous rock[END_REF] (see also [START_REF] Gens | In situ behaviour of a stiff layered clay subject to thermal loading: observations and interpretation[END_REF] and extended during the European TIMODAZ project [START_REF] Arson | Thermal Impact on the Damaged Zone Around a Radioactive Waste Disposal in Clay Host Rocks -Deliverable D10: Constitutive models development, ASTM[END_REF], is developed.

Bases of the original model for bonded soils.

The original constitutive model was defined for soft argillaceous rocks which are characterized by a low porosity, significant bonding and a stiffness decrease upon loading. Because of the high bonding existing between the clay minerals particles, the response of this type of rocks is controlled by micro-cracks development, fracture propagation, and existing discontinuities. Thus, the behavior of soft argillaceous rocks intermediates between a structured soil and a rock [START_REF] Leroueil | The general and congruent effects of structure weak rocks[END_REF].

It considers the material as made of two components: the soil matrix endowed with a soil-like behavior and the bonds whose response is controlled by a quasi-brittle type of law typical of material. When an external load is applied to the medium, bonds and soil matrix will carry part of the stresses according to their stiffness. Each component will then experience a strain controlled by the stress it knows and its proper constitutive law. The resulting strain of the medium is finally obtained as the one resulting from the strains of each component and must be compatible with the externally applied deformation.

According to this scheme, the model needs the definition of the constitutive law of each component as well as the partition of external stress and strain into the stress and strain of each component.

The main characteristic of the original model is that it considers the bonds as responsible of all structuring effects. Hence, the matrix has the same local behavior (expressed in terms of local matrix stresses and strains) independently of the structuration level. Therefore, its response can be expressed by a law of the unbonded material. A variety of models has been employed for this purpose, such as, Critical state models (Cam-Clay or CASM) or elastoplastic model provided with linear elasticity and different shear failure criteria (Mohr-Coulomb, Hoek and Brown).

On the other hand, bonds are modelled by a quasi-brittle constitutive law based on the isotropic damage elasticity developed by [START_REF] Carol | On the formulation of anisotropic elastic degradation. I. Theory based on a pseudo-logarithmic damage tensor rate[END_REF]. This framework aims to model the degradation in material stiffness with micro-cracks development with a reduced number of parameters (4 in total).

Partition of stress and strain requires considering the conceptual scheme shown in Figure 4-2 for the arrangement of bonds and matrix within the composite material. From this figure, several ratios between voids, bond and matrix can be defined:

a) The void ratio e of the composite material. It is the ratio between void (VV) and solid volumes (Vs=Vm+Vb, where Vm is the volume of solid grains and Vb the volume of bonds). The variation of void ratio by unit volume of composite material (V=Vs+Vv) is given by the so-called internal volumetric strains increment, defined as . dε vm corresponds also to the volumetric strain increment externally applied. Based on these definitions, the strain partition is given by equation Eq. 4-1, where Cb is the bond content defined as Vb/V:

dε v int = dV V V = - de (
dε v ext = dε v int + C b dε vb Eq. 4-1
As for strains, three local stresses are defined.  b are the stresses of the bonds. They are related to  b by the constitutive law defined for the bond material (isotropic damage).  m are the stresses of the matrix and are related to  m by the constitutive law defined for the matrix. the stresses of the composite medium. By applying the principle of virtual work on the representative element presented in Figure 4-2, the stress partition appears to be given by equation Eq. 4-2:

σ = (1 + X)σ m + Xσ b Eq. 4-2
where χ is the coefficient defined by Eq. 4-3

C b dε vb = χdε v int
Eq. 4-3

The coefficient χ states that changes in porosity ((1 + 𝑒)dε v int ) are directly proportional to the strain in the bonds. According to Eq. 4-2, this constrain on internal deformation provides also the redistribution of external stresses over the bonds and the matrix in order to conserve the work. It can be considered as a descriptor of the effect of microstructural arrangement on the stress-strain relationship of the composite medium.

Bases of the updated model for two-minerals rocks.

The main characteristic of the model developed for bonded soils relies on the fact that the load externally applied on the composite material is first carried by the matrix and then transferred to the bonds through the matrix only. In other words, as drawn in Figure 4-2, bonds never receive directly the load externally applied on the composite because they do not intersect the boundary of the representative element. This assumption has strong implications in terms of stress and strain partitions, since it states that the strain suffered by the matrix is equal to the external strain (or strain of the composite material). This microstructural arrangement has been considered in order to represent essentially materials composed by a soil matrix reinforced by cements at particle contacts.

The picture is quite different in case of rocks where all the mineral components are quite evenly distributed in the material and there are no components acting as a shield for the other ones. The configuration shown in Figure 4-2 has then to be modified for these materials. This modify the stress and strain partition rules and changes significantly the behavior of the composite material.

The new representative element scheme is presented in Figure 4-3. In this case bonds play the same role as the matrix and also receive directly part of the external load. In this new configuration, components can be exchanged without modifying the response of the material. This configuration is closer to a system in parallel than the one proposed by [START_REF] Vaunat | Bond degradation and irreversible strains in soft argillaceous rock[END_REF] although differences still exist with the parallel system because the use of virtual work forbid the possibility to set the same strain in both components. Based on Figure 4-3, the following ratios can be defined:

a) The total volume (V) is equal to the sum of voids volume (Vv), matrix volume (Vm) and bonds volume (Vb): V=Vv+Vm+Vb b) The void ratio for the composite material (e) is the ratio between voids volume (VV) and solids volume (Vs=Vm+Vb). Variations in void ratio by unit volume of composite material (V=Vs+Vv) is given by the internal volumetric strains increment, defined as:

dε v int = - dV v V - de (1+e)
c) The bond ratio eb, gives the bond volume (Vb) in solid volume (Vs). Its change over itself corresponds to the bond volumetric strain increment:

dε vb = - dV b V b = - de b e b
d) The matrix ratio em, gives the matrix volume (Vm) in solid volume (Vs). The change in matrix ratio over itself corresponds to the matrix volumetric strain increment:

dε vm = - dV m Vm = - de m e m
e) The quantities Cb and Cm give the bond and matrix volumetric contents (i.e. per unit volume of composite): Cb= Vb/V and Cm= Vm/V. These two quantities are related to the porosity of the composite () by: Cb+Cm=1- These different ratios provide the bases for the stress and strain partition developed in the next section.

In this configuration, we can formally exchange the role of each component (that is of subscript m and b) without modifying the structure of neither the stress and strain partitions nor the tangent matrix.

At microscopic level, the thermal damage effect is mostly evidenced at the boundary between minerals (Figure 4-1), nevertheless in order to reproduce the behavior observed at macroscopic level (weakening of the material) this effect is only attributed to the material defined as bond. For the limestone of LaRG that is mainly composed by calcite and quartz, it is considered that the bonds represent the calcite while the matrix defined the presence of quartz.

Model formulation.

This section presents the formulation and numerical implementation of the updated model, distinguishing the constitutive laws used for the matrix and the bonds and the procedure to couple them.

Matrix constitutive law.

The matrix corresponds to the mineral component that does not include the effect of microcracks. It is thus provided with any elastoplastic constitutive law. For this purpose, they are several models: linear or non-linear isotropic and cross-anisotropic elasticity with two types of failure criterion (Mohr-Coulomb or Hoek & Brown) or critical state models of Cam Clay or CASM type [START_REF] Yu | CASM: a unified state parameter model for clay and sand[END_REF]. In all cases, the flow rule can be considered as associated or no-associated and several types of hardening and softening relations can be considered. Effect of suction and temperature [START_REF] Gens | Soil-environment interactions in geotechnical engineering[END_REF] are also considered in most of these models, see CODE_BRIGHT User's guide (Departament d'enginyeria del terreny, 2017) In this work, the linear elastic model provided with Mohr-Coulomb yield criterion (see Eq. 4-4) is considered.

F = (cos θ m + 1 √3 sin θ m sin ϕ′ m ) J m -sin ϕ ′ m (P m + P t m ) ≥ 0 Eq. 4-4
where:

P t m = 𝑐´𝑚 cot ϕ′ m is the matrix tensile strength, 𝑐´𝑚the matrix cohesion, ϕ ′ m matrix friction angle, P m is the matrix mean stress, J m the second invariant of matrix deviatoric stress tensor and  m the Lode angle.

Nevertheless, it will be see during the analysis that, because of test conditions, matrix never reaches the yield criterion and remains therefore in the elastic regime.

Bond constitutive law.

Bond corresponds to the mineral component where micro-cracks effect is included. It is consequently provided with a damage law. For the sake of the present analysis, the scalar damage proposed by [START_REF] Carol | On the formulation of anisotropic elastic degradation. I. Theory based on a pseudo-logarithmic damage tensor rate[END_REF] and enhanced with suction and temperature dependency is considered. This model has been chosen because of its simplicity and the conceptual basis on thermodynamic principles that prevent any spurious energy dissipation during load cycles.

According to damage theory, bond response is elastic, but with modulus that progressively decrease during loading as micro-cracks develop inside the material. This phenomenon is modeled by distinguishing: The relation between these quantities is defined through the parameter D, equal to the ratio of fissured bond area and the whole bond area. [START_REF] Cordebois | Endommagement anisotrope en élasticité et plasticité[END_REF] proposed to use a relation based on the equivalence of the energies computed externally (

σ ij b ε ij b 2
) and internally

( σ ̂ij b ε ̂ij b 2
). This condition imposes the following relation between the apparent and effective stress and strain:

σ ̂ij b = 1 √1-D σ ij b ε ̂ij b = √1 -D ε ij b Eq. 4-5
As mention before, the damage theory considers that the response of the material is linear elastic, therefore the intact (unfissured) part of the material follows the law:

σ ̂ij b = D ijkl eb0 ε ̂kl b Eq. 4-6
where D ijkl eb0 is the stiffness tensor for the undamaged material, which is considered constant in this formulation. From equations Eq. 4-5 and Eq. 4-6, it can be obtained the relation between apparent bond stresses and apparent bond strains, defined by the equation Eq. 4-7.

σ ̂ij b = D ijkl eb ε ̂kl b Eq. 4-7
Where:

D ijkl eb = (1 -D)D ijkl eb0
Eq. 4-8

D ijkl eb
is the stiffness tensor of damaged material, whose coefficients decrease as micro-fissuration increases. [START_REF] Carol | On the formulation of anisotropic elastic degradation. I. Theory based on a pseudo-logarithmic damage tensor rate[END_REF] proposed to use an alternative damage variable L defined as:

L = ln ( 1 1-D )
Eq. 4-9

The damage variable L has a null value for the undamaged state and infinitive value when the material is fully damaged. By considering equation 4-9, equations 4-5, 4-7 and 4-8 can be rewritten as:

𝜀̂𝑖 𝑗 𝑏 = 𝑒 -𝐿 2 ⁄ 𝜀 𝑖𝑗 𝑏 Eq. 4-10 σ ̂ij b = e L 2 ⁄ σ ij b Eq. 4-11 D ijkl eb = e -L D ijkl eb0
Eq. 4-12

From equation 4-12, L can be simply related to the degradation of bond Young modulus by:

L = ln ( E 0 E ) Eq. 4-13
where E0 and E are the undamaged and damaged bond Young's modulus, respectively.

Therefore, total relative change in volume is directly obtained by the addition of the relative changes in void, matrix and bonds volumes:

dV V = dV v V + dV m V + dV b V = dV v V + C m dV m V m + C b dV b V b
Eq. 4-20

Where

C m = V m V , C b = V b V = V b V v +V b V v +V b V = S b ϕ m and C b + C m = V b V + V m V = V-V V V = 1 -ϕ
Cm corresponds to matrix volume concentration, Cb to bonds volume concentration,  is the porosity of the composite medium, ϕ m is the porosity between matrix grains (

V v +V b V
) and Sb the degree of saturation of bonds

( V b V v +V b
). From equation 4-20 the strain partition can be defined as (soil mechanics convention is considered):

dε v ext = - dV V = - dV V V -c b dV b V b -c m dV m V m = dε v int + c b dε v b + c m dε v m Eq. 4-21
where dε v ext is the volumetric strains applied to the composite medium, dε v int the volumetric strain used to update the porosity, and dε v b and dε v M the volumetric strains inside the bonds and matrix respectively.

The cornerstone of the model is the statement of two coefficients relating bond and matrix strain increments to the strain increment used to update porosity:

X b = c b dε ij b dε ij int Eq. 4-22 X m = c m dε ij m dε ij int Eq. 4-23
Xm and Xb express the sharing existing between solid component and pore deformations during an incremental straining process. They depend generally on microfabric patterns and relative stiffness between porous medium and its solid components.

These coefficients are considered as descriptors of rock microfabric at a given state of damage.

It is moreover expected that part of the load carried by the bond is transferred to the matrix when bonds degrade, leading to a decrease in the value of Xb. To account for this effect, coefficient Xb is set as a function of damage variable L. The function selected is given by:

X b = X b0 e -L 2
Eq. 4-24

Using Xm and Xb and extending equation 4-21 to all strain components, it is possible to compute bond and matrix strain increments from the external one:

dε ij ext = (1 + X b + X m )dε ij int → dε ij int = dε ij ext 1 + X b + X m Eq. 4-25 dε ij b = X b /c b 1 + X b + X m ⏟ b ε dε ij ext = b ε dε ij ext Eq. 4-26 dε ij M = X m /c m 1 + X b + X m ⏟ U ε dε ij ext = U ε dε ij ext Eq. 4-27
Stress partition are obtained by stating that the work input into the medium under any compatible external strain increment is equal to the sum of the work input into the matrix and the bonds. Labelling σ ij ext , σ ij m and σ ij b the stresses acting on the composite medium, matrix and bond, respectively, this condition reads:

σ ij ext dε ij ext * = c m σ ij M dε ij M * + c b σ ij b dε ij b *
Eq. 4-28

where

dε ij ext * = dε 𝑖𝑗 int + c m dε ij M * + c b dε ij b *
Eq. 4-29

By dividing both sides of equation by dε ij ext * , the stress partition reads:

σ ij ext = X m 1 + X b + X m ⏟ U σ σ ij m + X b 1 + X b + X m ⏟ b σ σ ij b Eq. 4-30
In this framework, matrix strain tensor ε ij m is associated to matrix stress tensor σ ij m through the constitutive law of the soil matrix while bond strain tensor ε ij b is associated to bond stress tensor σ ij b through the constitutive law for bonds:

dσ ij m = D ijkl m dε kl m
Eq. 4-31

dσ ij b = D ijkl b dε kl b Eq. 4-32
Finally, from equations 4-26, 4-27, 4-30, 4-31 and 4-32, the constitutive law for the composite material reads:

dσ ext = U σ dσ m + b σ dσ b + σ m dU σ + σ b db σ = U σ D m dε m + b σ D b dε b + σ m dU σ + σ b db σ
or:

dσ ext = U σ U ε D m dε ext + b σ b ε D b dε ext + σ m dU σ + σ b db σ Eq. 4-33
The change in U and b are related to the variations in Xb because the change in the damage multiplier. The equations governing this process are:

dX b = - 1 2 X b0 e -L 2 dL = - X b 2 dL
Eq. 4-34

db σ = d ( χ b 1+χ b +χ m ) = - X b (1+χ m ) 2(1+χ b +χ m ) 2 dL
Eq. 4-35

dU σ = d ( χ m 1+χ b +χ m ) = X b X m 2(1+χ b +χ m ) 2 dL
Eq. 4-36

Finally, the tangent matrix of the composite medium is obtained as:

D ext = U σ dσ m dε ext + b σ dσ b dε ext + σ M dU σ dε ext + σ b db σ dε ext
Or, by reworking the equation:

D ext = U σ U ε D m + b σ b ε D b + (σ M dU σ dL + σ b db σ dL ) dL dε ext D ext = U σ U ε D m + b σ b ε (D eb -σ b 1 b ε dL dε ext ) + (σ M dU σ dχ b + σ b db σ dχ b ) dχ b dL dL dε ext D ext = U σ U ε D m + b σ b ε (D eb -[σ b ( 1 b ε - 1 b σ b ε db σ dχ b dχ b dL ) -σ M 1 b σ b ε dU σ dχ b dχ b dL ] dL dε ext ) D ext = U σ dσ m dε ext + b σ dσ b dε ext + σ M dU σ dε ext + σ b db σ dε ext D ext = U σ U ε D m + b σ b ε D b - X b 2(1 + χ b + χ m ) 2 ((1 + X m )σ b -X m σ m )dL Eq. 4-37
In absence of damage (dL = 0), the tangent matrix reduces to:

D ext = U σ U ε D m + b σ b ε D b Eq. 4-38
It is therefore the sum of the tangent matrix of both components, weighted by the product of stress and strain partition coefficients.

Model implementation.

The numerical implementation of the model within the Finite Element framework requires to develop strategies to cope with the discrete equations composing it, that are:

-Damage constitutive model -Elastoplastic constitutive mode -Strain partition -Stress partition -Evolution equation of Xb coefficient Moreover, the tangent matrix of the model must be computed to feed the tangent matrix of the global iterative scheme that solves the non-linearities of the Finite Element solution. This calculation is based on the expression given by equation 4-37 that requires the computation of:

-The tangent matrix of the damage model -The tangent matrix of the elastoplastic model -dL as a function of the external strain increment One of the problems found during the integration of the damage model relies on the use of an exponential law for the variation of the stiffness and compliance matrix for the bonds. When linearized, this type of law requires the use of very small time steps to reach acceptable precision. To solve this issue, an alternative non-linear discrete approximation of the exponential has been implemented. It is presented in Appendix C.

Bounds for Xm and Xb coefficients.

Two central parameters of the model are the coefficients Xm and Xb. They provide indeed the relative importance of the matrix and bond models for the composite. These parameters are related to microstructural aspects of the material and can be usually only indirectly derived through back-analysis of the rock stiffness.

The stiffness of the rock (composite material) can be determined through mechanical tests or calculated from elastic wave propagation velocities (as presented in section 3.4.3). Nevertheless, these parameters can also be estimated from the properties of the components of the rock and the pore fluids. In fact, this value depends on the individual stiffness of the constituents, their volume fraction and some geometrical details of how the components are arranged.

These geometrical conditions are the most difficult characteristics to consider, thus in absence of the stiffness moduli of the composite material the best option is to estimate an upper and lower bounds, for it Hashin & Shtrikman (1963) proposed a method to obtain these bounds without making any assumption about the geometrical disposition of the phases (spheres embed in a matrix or intercalations of phases bands), and defining them through the strains energy stored in the material when subjected to uniform stress or strains.

According to Hashin & Shtrikman (1963) the upper and lower bounds for bulk modulus are defined according to equations 4-39 and 4-40 respectively:

K 1 * = K 1 + A 1 1 + α 1 A 1 Eq. 4-39 K 2 * = K n + A n 1 + α n A n
Eq. Where:

α 1 = - 3 3K 1 + 4G 1 , α n = - 3 3K n + 4G n A 1 = ∑ C r 1 K r -K 1 -α 1 r=n r=2 , A n = ∑ C r 1 K r -K n -α n r=n-1 r=1
K1 is the highest bulk modulus of the components and Kn is the lowest value, G1 and Gn are the shear moduli of the correspond components and Cr the volume fraction of component r, and n represents the number of components.

In order to estimate the upper and lower bounds of bulk modulus for the limestone of LaRG three different components are considered (n=3 for Eq. 4-39 and Eq. 4-40): calcite, quartz and air, this last component represents the porosity of the material under dry conditions.

According to [START_REF] Carmichael | CRS practical handbook of physical properties of rocks and minerals[END_REF], the typical values of bulk modulus (K) and shearing (G) modulus for the calcite and quartz, obtained from characteristic values of elastic wave propagation velocities are K= 75 GPa G=30.6 GPa for the calcite and K=38 GPa and G= 44 GPa for the quartz.

Considering these values, and taking into account the two mineralogical compositions observed in samples from the site under study (see section 3.3.1) the upper and lower stiffness bounds for LaRG limestone are estimated according to the method presented by Hashin & Shtrikman (1963).

According the mineralogical analyses, samples obtained from blocks are evaluated with a quartz content of 10%, while samples obtained from the cliff face are considered with a quartz content of 40%. From equation 4-38 it can be easily proven that the undamaged, (that is when dL = 0) bulk modulus of the composed material depends on the bulk modulus and content of matrix and bond materials and on the coefficients Xm and Xb, following the equation 4-41. Hence, it can be considered that the coefficients Xm and Xb represents those geometrical characteristics of the arrangement of components, this is to say the microstructural features of the composed material [START_REF] Mitchell | Soil deposits -Their formation, structure, geotechnical properties, and stability[END_REF].

K = ( X m 1 + X b + X m ) 2 1 C m K m + ( X b 1 + X b + X m ) 2 1 C b K b
Eq. 4-41

Consequently, if Cm, Cb, Km, Kb, and K are known, it is possible to determine all possible pairs of Xm-Xb that guarantee the value of bulk modulus for the composite material (K). Moreover, even if Xm and Xb can reach theoretically any value, these coefficients are restricted by the lower and upper bounds of K, that can be calculated according to equations 4-39 and 4-40

For instance, samples from La Roque Gageac obtained from blocks with a porosity equal to 0.2 should present a bulk modulus that varies between 0.75 GPa to 42.2 GPa, thus the coefficients Xm and Xb must be defined in the range presented in Figure 4-6. It is important to mention that the lower bound is mainly control by the porosity and air parameters. Thus, if the sample was fully saturated instead to unsaturated this boundary will be higher, as the bulk modulus for the water is higher than for the air. Moreover, the coefficients Xm and Xb also controls the strains (Eq. 4-25) and stress (Eq. 4-30) partitions, thus variations in Xm and Xb may generate higher strains and stresses in the bonding inducing a higher damage in the composite material, even if the initial stiffness is kept constant.

Model verification.

Finite Element formulation.

In order to model Boundary Value problems, the constitutive model has been implemented in the Finite Element code CODE_BRIGHT developed at the Universitat Politècnica de Catalunya [START_REF] Olivella | Numerical formulation for a simulator (CODE_BRIGHT) for the coupled analysis of saline media[END_REF]. It deals generally with the solution of thermo-hydro-mechanical coupled problems in geological media. However, since there is no fluid present in the medium under study, only the thermo-mechanical coupling will be considered in the analysis. As discussed by [START_REF] Gens | In situ behaviour of a stiff layered clay subject to thermal loading: observations and interpretation[END_REF], in case of stiff materials like rocks, the thermo-mechanical coupling acts mainly in one direction (Figure 456). Under the presence of energy fluxes and temperature changes, material experiences thermal strain which affects in turn the stress state, and eventually, the rock stiffness (bold downward arrow in Figure 456). On the other hand, since rock thermal properties depends essentially on mineralogy and porosity, whose changes are low in stiff materials, there will be almost no effects of the mechanical response of the material on the thermal field (thin upward arrow). Equation Eq. 4-42 allows relating the volumetric deformations in the medium, given by the divergence of solid velocity js, with the porosity changes. Equation Eq. 4-43 is a reduced equation of the general equation of energy balance, where terms of energy fluxes convected by the fluid phases have been removed, since there are no fluids in the medium. As well, because of the assumption of small strain, the energy flux of the solid phase is neglected. As a result, the energy flux is only due to conduction (ic) in the problem considered.

Test characteristics.

For the purpose of model verification, a simple Boundary Value Problem is considered. The geometry of the model considers the cube shown in Figure 4-8. Modeling parameters considered are summarized in Table 4-1. Material is provided with porosity equal to 0.2, matrix and bond contents equal to 0.4 (Cb=Cm=0.4) and matrix and bond bulk modulus equal to 50 GPa. Under these conditions, the bulk modulus of the composite material is 25 GPa. Thermal expansion of both components is considered equal to 6. 10 -6 ºC -1 . A thermal ramp, starting from the initial temperature (0 ºC) up to 50 ºC and performed in 1h, is simultaneously applied to all the volume under null stresses. Under these conditions, the strain, stress and thermal fields are constant over the mesh and the response output by the model is a direct reproduction of the thermo-mechanical constitutive law material implemented at Gauss points. The solution obtained by the numerical model will be compared with an algebraic expression obtained by integrating explicitly the constitutive law along the imposed thermal path under null stress. The integration has been carried out considering temperature increments T equal to 1ºC. To make easier the comparison, no strain rate effect has been considered in the numerical model and algebraic expression.

Comparison with the algebraic expression.

The algebraic expression is obtained on the basis of the following equations:

1) The strains partition defined by Eq. 4-21:

δε = δε int + C b δε tot b + C m δε tot m
Eq. where:

δε tot = δε mec + δε term
Eq. 4-46 δε mec are the deformations associated to mechanical stresses and δε term the strains related to thermal changes:

δε term = α v ΔT Eq. 4-47
Where α v is the volumetric thermal expansion of the component and ΔT is the temperature variations.

It is important to mention that, as a result of the strain partition, thermally induced strains and thus stresses will develop in the two different components even if there are not mechanical loads externally applied.

2) The definitions for Xm, Xb and Ue(Eq. [4][5][6][7][8][9][10][11][12][13][14][15][16][17][18][19][20][21][22]:

δε tot m = X m C b δε tot b X b C m
Eq. 4-48

δε tot = 1 U e δε tot m
Eq. 4-49

3) The elastic law for both components:

δP b = K b (δε mec b ) = K b (δε tot b -δε term b )
Eq. 4-50

δP m = K m (δε mec m ) = K m (δε tot m -δε term m )
Eq. 4-51

4) The stress partition expression in absence of external loads (Eq. 4-30):

δp = U σ δp m + B σ δp b = 0
Eq. 4-52

5) The bond damage law (Eq. 4-24):

K b = K b0 e -λ db
Eq. 4-53

X b = X b0 e -λ db 2 ⁄
Eq. 4-54

6) The bond damage loading function:

F = u b -r = 0 Eq. 4-55
where U b is the energy in the bond defined as:

U b = 1 2 P b ε mec b Eq. 4-56
and r is the evolution of damage locus, considered linear:

r = r 0 + r 1 λ db
Eq. 4-57 r0 is the initial damage locus and r1 is the damage rate.

As the variations in temperature is known, it is possible to obtain the thermal strains for the bonding and matrix, related to thermal expansion coefficients. Considering these values and the expressions Eq. 4-48 and Eq. 4-51 the δε tot b can be obtained from equation Eq. 4-52.

Then, following the expression Eq. 4-48 the δε tot u can be obtained, is then possible to calculate the bond and matrix stresses (p m and p b ), allowing the computation of the new Xb and Kb related to the damage law of the bond. It is important to mention that here in after the strains are draw with rock mechanics convention, in other words, positive values express expansion. Results obtained with both solutions are very similar, however, some small differences are observed at the end of the heating interval, possible related to the integration scheme of the finite element model.

From Figure 4-10 it is important to point out the strong relationship between the bulk modulus of the composed material and the damage multiplier, when this last value increases a reduction in the bulk modulus is registered. Note that until 0.15 hours the damage multiplier is constant and equal to zero, therefore K doesn't present any evolution, remaining equal to its initial value (25 GPa). According to the equation 4-14, damage envelope F can be rewritten in terms of p b , q b and db by considering the linear damage evolution law (Eq. 4-15):

F = 1 2 (p b -p 0 b ) 2 K b + 1 2 (q b -q 0 b ) 2 G b -(r 0 + r 1 λ db ) = 0
Eq. 4-58

Figure 4-11 presents the damage envelope F = 0 at three stages during the first thermal ramp: 1) at t = 0 (red curve). Damage locus is given at the initial value of r0 (energy equal to 1. 10 -5 kPa); 2) at t = 8 min, when temperature reaches the value of 10ºC. At that time, damage locus has already experienced some hardening due to bond stress increase resulting from thermal expansion. It appears from the graph that the onset of damage occurs for a temperature close to 5ºC, time at which pb reaches 1 MPa and the corresponding energy reaches r0. On the other hand, when temperature reaches 30ºC, the bond mean stress increases until 5.5 MPa, bringing with it the damage locus and thus creating damage. During the experimental work (see chapter 3) it was observed that damage accumulates during cycles. Moreover, as mentioned before, the thermal damage is mostly related to the propagation of existent microcracks, which is a time dependent phenomenon.

This indicates that the capture of damage accumulation during thermal cycles requires the introduction of an additional mechanism related to the time dependency of thermal stress induced micro-cracks propagation.

For this reason, a rate dependent damage locus has been introduced in the model. The rate dependency is controlled by the bond viscosity  b , taken equal to 100 in the present model. All other test characteristics are kept equal. Figure 4-13 presents the evolution of the visco-damage loading function F during the thermal loading.

The rate-dependency allows the bond mean stress to stand outside the loading function for a time that depends on the velocity of hardening of the surface, controlled by the value of  b . In this way, a delay is introduced in the development of damage, leading to a more progressive transition between two stages of damage. As well, if temperature cycles are applied at time intervals shorter that the time required to reach the full damage, bond stresses will be incomplete during a single interval and will proceed during the next one, simulating the rate dependency of micro-cracks development. During cooling, bond mean stress comes back to lower value inside the elastic zone. The model predicts therefore no damage during this stage. The constitutive law presented in section 4.3 considers an important number of parameters, some of them are physical properties of the rock, such as porosity, matrix and bond content and stiffness parameters of each component. Others are microstructural characteristics like Xm and Xb, or modeling adjustment parameters as r0 and r1.

In general, from all these parameters, only the porosity and mineralogical content can be easily defined. For this reason, it is interesting to make a sensitivity analysis of the effect of different parameters in the damage evolution of the composite material.

In this perspective, the simplified problem considered in the last section will be analyzed with different sets of parameters. In order to make clearer the effect of each parameter, no time rate is considered.

The sensitivity analysis holds on the ratio between matrix and bond bulk modulus (Km/Kb), thermal expansion (am/ab ratio) and bond content (Cb). Each analysis will be performed for three different Xm/Xb ratios: 0.5, 1 and 2.

The effect of r1 and r0 are not affected by Xm and Xb, in general, a reduction of any of these two parameters will generate higher values of the damage multiplier, therefore a higher reduction in the composed material bulk modulus.

Bulk modulus.

Three different scenarios for the bulk modulus are considered, Km/Kb = 0.5, 1 and 2. For the sake of comparison, the value considered for the highest component stiffness is 50 GPa. As such, Km/Kb = 0.5 means that Kb= 50 GPa and Km=25 GPa whether Km/Kb = 2 means Km= 50 GPa and Kb = 25 GPa. Volumetric content and thermal expansion are taken equal to respectively 0.4 and 6. 10-6 ºC for both the matrix and the bonds.

Furthermore, it has been verified that coefficients Xm and Xb belong to the admissible zone of bulk modulus. The first case considers a material where the bond (damaging component) has a higher effect than the matrix (no-damaging component) on rock response.

Figure 4-15 shows the time evolution of total, bond and matrix strains. When Xm is lower than Xb, the total strains registered in the bonds will be higher than those obtained for the matrix, as can be observed in Figure 4-15. Moreover, the highest total strains (whatever is for the composite material, for the matrix or for the bonds) are obtained when the bulk modulus of the matrix is higher than the one of the bonds. Under null applied external stress, it will be expected that total strains (etot) were only controlled by thermal expansion of the materials. Nevertheless, because the microstructural characteristics of the composed material (related to parameters Xm and Xb) the total strains can be higher or lower to thermal strains as a result of restrictions imposed by the presence of two heterogeneous components. The first restriction comes from the stress partition (Eq. under the applied external stress (if the present case equal to 0). Figure 4-16 shows the time evolution of mean bond and matrix stresses resulting for the thermal load. It can be observed that, whatever is the ratio between the bulk moduli of the components, matrix presents compressive stresses, while bonds develop stresses in extension. This is clearly because, when one component extends, the other must compress to maintain the overall stress equal to zero. Moreover, the lowest stress magnitude is observed for the case Xm = Xb, that is when both components have a similar contribution to rock response. This indicates that internal stresses develop essentially due to the difference existing between the responses of the two components. Additional restrictions are the constitutive laws of the bonds and the matrix that relates internal stresses to internal mechanical strains. Time evolution of bond and matrix volumetric strains is shown in Figure 4-17. As the ratio Km/Kb increases, more mechanical strains develop in the bonds and less in the matrix. Moreover, for this specific case, bond mechanical strains appear to be higher when the bulk modulus of the matrix is higher than the bulk modulus of the bond. Conversely, matrix strains are lower than bond mechanical strains when Km = 2 Kb. This result cannot be directly extrapolated to stresses since they are related to the strains through the stiffness of the component. For example, matrix strains are higher than bond mechanical strains for Km = 0.5 Kb while stresses are lower because of the low stiffness of the matrix.

The final restriction is the additivity of mechanical and thermal strains. It imposes that the strain measured externally is not only due to the thermal expansion of the rock but depends also on the values of mechanical strains experienced by rock components. It is why the external strain depends on the ratio between bond and matrix stiffness.

The same explanation holds for the occurrence of rock damage under thermal load. Figure 4-18 shows the evolution of bond damage multiplier ( db ) together with the degradation of rock bulk modulus K. For ratios Km/Kb equal to 0.5 or 2, bonds are more stressed than for ratio Km/Kb = 1 and  db variable increases. Therefore, the bulk modulus of the rock experiences a higher damage. There is no proportional relationship between rock degradation and Km/Kb, but it must be investigated looking at the values of stresses in the damage-prone component of the material Consequently, no damage will be registered in the bond (damage multiplier equal to zero) and the composite material bulk modulus will be constant for the complete thermal load (see . This corroborates the fact that, in the model under consideration, development of thermally induced damage under constant stress is caused by the differential response of its components. The last case (Xm larger than Xb) considers a material where the matric (not damaging component) has a more pronounced effect than the bonds (damaging component) on rock response.

In this case, the general response obtained will be the opposite of the behavior found for Xm lower than Xb. In first place, the total strains are higher for matrix strains are higher than bond strains and close to the values of external strains Additionally, the highest total strains (either in the bond, matrix and composite material) occurred when the bonds are provided with the largest bulk modulus., see Figure 4-22. For this case, total strains in the matrix are higher than thermal strains, therefore, the mechanical strains and internal stresses generated for this component will be positive (extension), as presented in Figure 4-23 and Figure 4-24. Conversely, the bond develops contraction strains and stresses. Nevertheless, just as it was observed for Xm/Xb=0.5 the bond strains are greater (in magnitude) when the matrix bulk modulus is higher than the bond bulk modulus. Because of the different bulk moduli considered, the maximal internal stresses are not directly related to the maximal mechanical strains. In fact, the highest values of internal stresses are developed when Km is equal to Kb, which is not the case where the higher mechanical strains are developed, as can be observed in Figure 4-24 and Figure 4-23. It is interesting to point out that, even if the highest mechanical strains for the bond are developed when Km/Kb=2, the lowest P b is related to this case. Which can be attributed to the fact that for this specific case the initial Kb is equal to 25 GPa while the other two cases under analysis present a Kb equal to 50 GPa, yet, it is important to remember that Kb is affected by the damage multiplier. From Figure 4-25, it can be observed that the damage multiplier increases the most when the two components (matrix and bond) present the same bulk modulus, thus, for this case the reduction in the bond bulk modulus is higher too. Nevertheless, the reduction of bulk modulus for the composite material is higher when the matrix bulk modulus is lower than the bond bulk modulus, this characteristic is related to the values of Xm, Xb and the matrix bulk modulus considered. As for the situation under analysis Xm is larger than Xb, the stiffness of the composite material is going to be strongly related to the matrix bulk modulus, and remember that for the case km/kb=1, Km is equal to 50 GPa while, for the case km/Kb=0.5 it is equal to 25 GPa.

Moreover, it is important to remark that the increase in the damage multiplier obtained for Xm/Xb=2 is much higher than for the case where Xm/Xb=0.5 (Figure 4-18) because the magnitudes of P b obtained for this last case. Nevertheless, the stiffness of the composite material is going to be strongly related to Kb if Xm/Xb=0.5 and to Km if Xm/Xb=0.5. Consequently, the reduction in the composite material bulk modulus is very similar for both cases.

Thermal expansion coefficient.

In order to evaluate the effect of the thermal expansion coefficients, three different situations are evaluated: mb=0.5 (m = 3e-6 ºC -1 and b = 6e-6 ºC -1 ), mb=1 (m = b = 6e-6 ºC -1 ) and mb=2 (m = 6e-6 ºC -1 and b = 3e-6 ºC -1 ). For the sake of simplicity, thermal expansion is considered isotropic. The value of the volumetric thermal expansion coefficient is thus equal to three times the value of the linear one.

The sensitivity analysis on m and b will be performed for cm = cb = 0.4, Km = Kb and three different relations for Xm/Xb (0.5, 1 and 2).

Similar to the evaluation of the bulk modulus, this analysis is performed by considering three different relations for Xm/Xb, in this case they will be related to the bound defined for Km=Kb (red line in Figure 4-14)

Xm/Xb=0.5 (Xm=0.2 Xb=0.4).

When Xm is lower than Xb the total strains obtained for the bond are higher than those obtained for the matrix, just as was observed before (see section 4.7.1). See Figure 4-26. Furthermore, the highest total strains are related to higher values of m+b, obtained for the case mb=1. It is interesting to point out that for the cases mb=0.5 and mb=2 where the addition of the coefficients is equal (m+b=9e-6 ºC -1 in both cases), the highest total strains are obtained when m is lower than b. Moreover, for this case (mb=0.5) the bond and matrix total strains are equal to the strains due to thermal changes. Thus, no mechanical strains and stresses are obtained, as can be observed in Figure 4-27.

In addition, from Figure 4-27 it can be observed that while the thermal expansion coefficient ratio (mb) is larger than the ratio Xm/Xb, the bond will develop positive stresses (expansion) which mean that the total strains are greater than the strains obtained only by thermal changes. Otherwise, if mb is lower than Xm/Xb the total strains obtained for the bond are lower than the strains related to the temperature variations. Thus, the bond present contraction stress (negative values for P b ). The opposite response is observed for the matrix. Furthermore, the highest magnitude of P b are obtained when the thermal expansion coefficient is larger for the matrix than for the bond, as the initial Kb considered for all the cases is the same, this is also related to higher mechanical strains. Hence, the highest increment in damage multiplier is also developed for mb=2, as presented in Figure 4-28. Consequently, the bond bulk modulus will decrease the most for this case having as consequence the highest decrease in the composite material bulk modulus, because Km is equal for all the cases under analysis. If Xm is equal to Xb, the maximal total strains are going to be related to the case where m+b is higher (m+b=1). If the addition of thermal expansion coefficients of the components is equal, the total strains for the bond ( b tot) and the matrix( m tot) developed are the same no matters if the thermal expansion of the bond is larger than the thermal expansion for the matrix or vice versa (see mb=0. Moreover, if m/b is higher than Xm/Xb the bond will develop expansion stresses and mechanical strains, while the matrix will develop contraction stresses and mechanical strains. The opposite response is observed when m/b is lower than Xm/Xb. It is important to realize that because of the values of thermal expansion considered for the cases m/b=0.5 and m/b=2 the magnitude of internal stresses developed in the matrix and in the bond, are similar.

This condition may explain why the evolution of the damage multiplier obtained for these two cases is very similar (Figure 4-31) just like the evolution of the composed material bulk modulus.

As any mechanical strains is developed when m is equal to b no evolution for the damage multiplier or the composed material bulk modulus is obtained. Just as observed in the analysis of bulk modulus, for Xm larger than Xb, the response obtained is contrary to what it is obtained when Xm is lower than Xb. The total strains are higher for the bond than for the matrix, and the larger total strains are obtained when m/b is superior to Xm/Xb. Similar to the two scenarios evaluated before, if m/b is equal to Xm/Xb the total bond and matrix strains will be equal to the deformation related to the thermal variations. Therefore, no mechanical strains and internal stress are developed for this case, as observed in Figure 4-33.

In the other hand, if m/b is higher than Xm/Xb the bond will develop expansion and the matrix compression stresses, the opposite response is obtained if m/b is lower than Xm/Xb, nevertheless the largest (in magnitude) internal stress are obtained (for the matrix and the bond) for the lowest m/b ratio. In the same way, the thermal damage multiplier increases the most as the ratio between m and b is lower (Figure 4-34), hence the composed material bulk modulus reduces the most when m is inferior to b. It is interesting to point out that, for all the cases evaluated the maximal damage multiplier is obtained when b is higher than m (m/b=0.5) and Xm is higher than Xb (Xm/Xb=2). However, for this case the final bulk modulus for the composite material is only 95% of the initial one. Whereas, if b is lower than m (m/b=2) and Xm is lower than Xb (Xm/Xb=0.5), the final bulk modulus is 80% of the initial one. This characteristic is related to the values of Xm and Xb considered. If Xm is higher than Xb the matrix bulk modulus controls in the calculi of the composite material bulk modulus, otherwise, the bond bulk modulus (the one affected by the damage law) will be predominant.

Bond content.

In order to evaluate the effect of bonding content (Cb) and therefore matrix content (Cm) in the response obtained with the composite material constitutive law, three different scenarios of Cb are considered: Cb=0.2 (Cm=0.6), Cb=0.4 (Cm=0.4) and Cb=0.6 (Cm=0.2).

According to the equation 4-41 the initial bulk modulus of the composite material depends on Cb and Cm, hence, the choice of parameters Xm and Xb must be adjusted to the bounds associated to the Hashin & Shtrikman (1963) bulk modulus bounds. Figure 4-35, presents the ratios considered in the present analysis. As observed in the former analysis, when Xm is lower than Xb the bond developed higher total strains than the matrix, yet, this is observed only if the bond content is lower or equal to the matrix content (Cb=0.2 and Cb=0.4).

Additionally, the highest values of total strains for the composite material are obtained for Cb=0.6, this is also observed for the matrix total strains. Nevertheless, the highest bond total strains are obtained for the lowest bond content, as illustrated in Figure 4-36.

For the case under evaluation (Xm/Xb=0.5), the bond develops expansion stresses when its content is lower or equal to the matrix one, whereas, the matrix is submitted to compression stresses, that means that the matrix cannot develop its complete thermal expansion. Otherwise, if the bond content is superior to the matrix, the bond is not able to develop its complete thermal expansion and is submitted to compression stress, as presented in Figure 4-37.

The highest (in magnitude) internal stress are generated when Cb is equal to 0.2, which is the minimal bond content considered. This may explain why this is the case where the highest damage multiplier is obtained and therefore the highest reduction in the composite material bulk modulus. See Figure 4-38. When equal values of Xm and Xb are considered, the highest total strains for the composite material are obtained when the material is composed by the same content of bond and matrix (Cb=0.4 in Figure 4-39). Moreover, for this case, the bond and matrix total strains are equal to the thermal strains, therefore no internal stresses are developed, as presented in Figure 4-40.

For a bond content of 0.2 and 0.6 the total strains for the composite material are very similar, it is important to remember that for these two cases the content of bond and matrix are alike too, for Cb=0.2, Cm=0.6 and for Cb=0.6, Cm=0.2. This characteristic generates that the bond stress for Cb=0.2 be equal to the matrix stress for Cb=0.6, as the thermal expansion coefficients considered are equal for both components. Moreover, for Cb=0.6 the matrix is submitted to expansion stress while the bond is submitted to compression. The completely opposite response is observed when the bond content is equal to 0.2, see Figure 4-40. Because the magnitude of bond stress and strains are equal when the bond content is 0.2 and 0.6, the damage multiplier coefficient evolves very similarly for these two cases, as illustrated in Figure 4-41, therefore, the bond bulk modulus reduces the same too. Nevertheless, for the composite material bulk modulus, a much higher reduction is obtained for Cb=0.2, because the bulk modulus for the composite material is affected by the bond bulk modulus in a rate that is inversely proportional to Cb, as the calculi of K is directly proportional to the coefficient Uthat is itself inversely proportional to Cb and defined according to the Eq. 4-27.

In other words, if the bond content is lower than the matrix content the bulk modulus of the composite material will be strongly related to the bond bulk modulus, (affected by the damage law) As observed before, the response obtained when Xm is larger than Xb is the opposite to what is obtained for Xm lower than Xb (Xm/Xb=0.5) In first place, the total strains for the composite material are higher for the lowest bond content, while for the case Xm/Xb=0.5 this is obtained for Cb=0.6. Moreover, the matrix total strains are higher than the bond total strains if the bond content is equal or higher to the matrix content, as shown in Figure 4-42 Because of the total strains developed in the components, some mechanical strains are generated and therefore internal stresses appear too. when Cb=0.4 and 0.6 the matrix develops expansion stresses, and if Cb is equal to 0.2, the matrix will be submitted to compression stress (see Figure 4-43). 4-44 it can be observed that the maximal bond stress is obtained when the bond content is 0.6, as usual, this is related to a higher increase of damage multiplier, therefore to a higher reduction of bond bulk modulus. Nonetheless, the composite material bulk modulus is not much affected by this reduction. This can be associated to two different phenomena: first, the value of Xm considered is larger than Xb, which generates that K be principally controlled by Km, additionally K is less affected by Kb as higher is Cb. This last condition explains why even if the damage multiplier for Cb=0.2 increases the less, the reduction of the bulk modulus of the composite is the highest. • If the coefficients related to the material microstructure, Xm and Xb, are equal and the composite material is form by two materials with similar properties, each component will be able to expand or contract without generating any internal stresses, therefore no damage is observed.

• The response obtained when Xm is higher than Xb is the opposite of the one observed when Xm is lower than Xb

• According to the stiffness of the composite material-K expression, this parameter is directly proportional to the matrix and bond stiffness and inversely proportional to the matrix and bond content. An interesting feature is that K will be strongly related to the matrix stiffness if Xm is higher than Xb and to the bond stiffness if Xb is higher than Xm.

• In general, the evolution of damage multiplier is directly proportional to the stresses developed in the bonding.

• The degradation of the bond stiffness is directly related to the damage multiplier increasing.

fraction with exception of sample CV5, whose proportion is around 40%/60%. Porosity has been obtained from laboratory density tests (see section 3.2).

In order to reduce the number of parameters to estimate, it has been considered that matrix and bonding materials are provided with the same thermal expansion coefficient. Under this assumption, the thermal coefficient can be back-analyzed for each sample from the external strains measurements presented in section 3.5.1. Is important to point out that, according to test results, thermal expansion coefficient has been considered anisotropic and two different values have been defined for the axial and radial directions.

Typical values have been taken from the literature for the bulk modulus of quartz and calcite phases [START_REF] Carmichael | CRS practical handbook of physical properties of rocks and minerals[END_REF]. Once fixed Km, Kb, cm and cb for each sample, it is possible to backanalyzed Xm and Xb from the measurements of rock bulk modulus presented in section 3.5.2 and under the restrictions provided by Hashin & Shtrikman bounds. Rock thermal properties have been taken from the measurements performed on sample B4 and presented in section 3.4.1.

Finally, damage properties have been back-analyzed from the experimental. A special care has been given to the search of unique values for samples belonging to the same type of material. between block and cliff samples the principal difference is the damage rate (r1), the cliff samples present a lower damage rate than the block samples, additionally for the block samples the parameter Xb is higher than Xm while the opposite condition is obtained for the cliff samples. The numerical modeling of samples considered was performed with two different approaches, as a Gauss Point and with an axisymmetric geometry. For the Gauss Point (GP, the geometry considered is a cube, in which the thermal variations are imposed to the complete volume, like the one presented in section 4.5 (Figure 45678). In other words, for this first approximation no boundary conditions are evaluated, therefore the effect of temperature propagation is not included.

The main result obtained with this approach is the evolution of the bulk modulus thought the imposition of thermal cycles. To reproduce the real experimental conditions, each sample is modelled with an axisymmetric geometry, that reproduce the shape and size of samples (cylinders with 5 cm of diameter and 10 cm high) and the temperature variations are imposed in the boundary of the geometry, as illustrated in Figure 4-48. Therefore, the effect of boundary conditions is included. Results obtained for each sample will be evaluated at the boundary and within the sample, points A and B in Figure 4-48, respectively. Results obtained with the axisymmetric models are presented next. It can be observed that the bulk modulus evolves less within the sample (point A) than at its face, where the thermal variations are imposed. Moreover, the stiffness evolution obtained with the Gauss Point model is located between the evolution observed on the face and inside the sample, nevertheless it is closer to the face one. Similar response can be observed in Figure 4-50, where the distribution of bulk modulus inside the samples is presented for different cycles.

It can be noted that in the bulk modulus is the lowest in elements where the temperature variations are imposed at two boundaries (the external corners).

This condition can be related to the thermal propagation inside the sample, in the center of it the temperature variations appear later than at the face. Moreover, the thermal cycle registered in the center of the sample can be a little bit smaller than in the face, which leads to smaller deformations, therefore to smaller bond stresses, finally to less damage and lower reduction in bulk modulus.

Similar results are obtained in the evolution of damage multiplier thought the imposition of thermal cycles, presented in Figure 4-51. Note that damage multiplier increases the most on the boundary of the sample (point B) than inside the sample (point A), this can be linked to lower values of bulk modulus at point B than at point A.

As mention on section 4.7, the evolution of the composite material bulk modulus when Xb is higher than Xm is directly related to the evolution of the damage multiplier. This may be explained why, the bulk modulus distribution inside the sample is very similar to the evolution of damage multiplier, as can be observed by comparing Figure 4-50 and Figure 4-52. Remember that thermal cycles generate an augmentation in the damage multiplier (higher values in red) related to a reduction in the bulk modulus (lower values in blue). The numerical modeling generates a small strains accumulation., contrary to what was observed in the experimental program (see section 3.5.1), where a high strains accumulation is registered. With the numerical modeling and after 1000 thermal cycles the strains accumulation in axial direction is 8.22e-6 m/m and in radial direction 1.26e-5 m/m both in negative direction (compression). Because the small strain variation due to thermal cycles, the reduction in bulk modulus is mostly related to a stress accumulation. This can be related to the conditions of the problem under evaluation, as the only external load imposed is the thermal variations, the strains associated will be related to the thermal expansion coefficient of the material, which is invariable. Moreover, the damage generates a stress transfer between the components, because the bond stiffness reduction.

Figure 4-55 presents the evolution of the mean stress with the imposition of thermal cycles at the two points of analysis within the sample (Point A) and at the boundary (Point B). It can be observed that higher stresses are developed inside the samples than at the boundary, moreover, the stresses at the boundary developed an important increment through the imposition of thermal cycles, while the stress at point A remains more or less constant for the whole analysis time. This can be related to the damage evolution obtained at these two points, remember that sample B3 presents more stiffness reduction and damage multiplier increment at the boundary of the sample than inside of it (see Figure 4-49 and Figure 4-51).

Sample B7.

During the experimental program, the response of sample B7 thought the imposition of thermal cycles was evaluated with the evolution of the elastic wave propagation velocities, no deformation measurements were performed in this sample. Therefore, the evaluation of numerical modeling results is mainly related to its stiffness evolution. Nevertheless, it is important to point out that the maximal damage multiplier obtained for sample B7 (0.12) is much higher than the one obtained for sample B3 (0.05). Because the higher values of initial bulk modulus of sample B3 and the lower reduction of this parameter observed for this sample, even if both samples were obtained from blocks.

The axial and radial thermal expansion coefficients selected for the numerical modeling of sample B7 were obtained from the deformation measurements obtained in other block samples (sample B1, B2 and B3). The radial and axial strains obtained are presented in Figure 4-60. In this case a small strains accumulation is obtained (-3.19e-5 m/m in axial direction and -4.15e-5 m/m in radial direction), nonetheless, the reduction in bulk modulus is mostly related to stress increasing. The results obtained for the bulk modulus evolution on sample CV4 is presented in Figure 4-61. It can be observed that the bulk modulus reduces more in the center of the sample (point A) than at the boundary (point B). The evolution obtained with the numerical model presents some differences with the measurements obtained for the first 500 cycles, after the results obtained follow closely the laboratory data. The response observed with the bulk modulus evolution, higher at point B than at point A, (opposite to what was observed in samples B3 and B7) can be related to the thermal expansion coefficients considered in the modeling of this sample, from Table 4-2 it can be observed that the thermal expansion coefficient in axial direction is 3 times the one considered in radial direction. This may generate some additional internal stresses inside the sample, inducing more damage, as can be observed in Figure 4-63 (higher damage multiplier at point A than B), and as result, higher reduction in the bulk modulus for point A than point B.

Moreover, it can be observed that the bulk modulus obtained for the Gauss Point model (Figure 4-63) is higher than the values obtained with the axisymmetric model, since for the GP model a volumetric value of thermal expansion was considered. In other words, no distinction between radial and axial coefficients was taken into account, thus no additional stresses that may affect the bulk modulus evolution were generated. As mentioned before, for the numerical modeling of samples drilled from the face of the cliff the parameter of Xm is higher than Xb, which makes that the evolution of bulk modulus is not directly related to the evolution of damage multiplier (see section 4.7).

Note that the damage multiplier for points A and B (Figure 4-63) increases for the whole modeling time (1200 cycles) while the bulk modulus (Figure 4-61) remains constant after cycle 500. Nevertheless, is important to point out that the damage multiplier presents a very steep tendency for the first 400 cycles, which can be related to a higher reduction in the stiffness.

Moreover, by comparing Figure 4-62 and Figure 4-64, it can be remarked that the distribution of stiffness in the sample CV4 is different to the damage multiplier distribution (conversely to what was observed in the block samples), but the elements with lower bulk modulus coincides with higher damage multiplier values. Additionally, the maximal values of damage multiplier are much higher (10 times higher) than those observed with the modeling of block samples.

The experimental program developed included the strains measurements in sample CV4, Therefore, the stiffness reduction observed is mainly related to an increment in the internal stresses. Figure 4-68 presents the evolution of mean stress at point B an inside the sample (Point A), it can be remarked that the mean stress developed in this cycle is higher within the sample than at the boundary, similar to what is observed on block samples (see Figure 4-55), and the stresses at point A presents shows an important reduction while the stresses at point B shows an slight increment.

Moreover, it can be observed that the mean stress reduction in point A is higher for the first 300 cycles, this can be related to the damage multiplier evolution (Figure 4-63) that shows an important increment for this period. Sample CV5 presents the highest decrease (34%) in bulk modulus for all the samples evaluated in the experimental program (section 3), moreover, the initial bulk modulus value for this sample is very low. From the mineralogical composition of samples (section 3.3.1) it was observed that samples obtained from the cliff presented a quartz content varies from 10% to 50%. This is why, the bonding and matrix content is considered equal for the numerical modeling of samples obtained from the cliff. However, with this condition it is not possible to obtain the initial bulk modulus value and the evolution observed for sample CV5. Consequently, the matrix content (quartz) was reduced, for this reason the relation between the bond and matrix contents considered for this sample is different to those used for the other cliff samples.

Despite this correction in the bond and matrix content, the bulk modulus obtained with the numerical modeling after 1200 thermal cycles is slightly higher than the value registered in the experimental study (500 MPa higher), see Figure 4-69. Nevertheless, it can be considered that the results obtained with the numerical model is a good approximation for the response of this sample, after 1200 thermal cycles a reduction of 27% in the bulk modulus is obtained. The important reduction on bulk modulus observed in the axisymmetric model for the first 200 cycles is linked to the evolution of the damage multiplier, presented in Figure 4-70. It can be remarked that the first cycles present a considerable increasing in the damage multiplier and after cycle 500 (where the bulk modulus tends to a stable value) the damage multiplier continue to increase but with a less steeper response.

Analogous to sample CV4 and as can be assumed from the bulk modulus evolution, the damage multiplier is higher at point A than at point B. Moreover, the value obtained after 1200 cycles is higher than the one registered for sample CV4 but with similar magnitude order.

During the experimental program, the deformation was not measured in sample CV5, thus the thermal expansion coefficient consider for this sample was adjusted form the strains measurements performed in other vertically drilled samples (CV1 and CV4). Figure 4-71 shows the strains response obtained with the axisymmetric model of sample CV5 in axial and radial direction. It can be observed that this sample shows a strains reduction with thermal cycles imposition, especially for the first 200 cycles, the same period that presents higher damage multiplier increasing. Again, the stiffness reduction is mostly related to the accumulation of stresses. The response of bulk modulus obtained for sample CH2 is presented in Figure 4-72, it can be observed that the numerical model reproduces adequately the data measured during the experimental program. In the same way to what was observed in the other cliff samples, the bulk modulus reduces the most in the inside the sample than at its boundary, and the results obtained with the GP model shows a lower reduction in the stiffness than the one observed with the axisymmetric model.

From Figure 4-73 it can be observed that the bulk modulus is smaller as closer is the element to the symmetry axis, as is expected from the bulk modulus evolution at points A and B. Moreover, for this sample the bulk modulus continues to decrease after 500 thermal cycles Moreover, as observed in the modeling of samples CV4 and CV5, the higher increments of damage multiplier observed for the first cycles (until cycle 300) is related to a higher reduction in the stiffness for this period. Similarly, to what was observed with the former numerical model results, while the measured data presents some important accumulation of deformation, a slightly strain accumulation is obtained with the numerical modeling, -7.44e-5 m/m in axial and -5.36e-5 m/m in radial direction.

Therefore, the variation in bulk modulus due to thermal cycles is mostly related to stresses accumulation in the composite material, as presented in Figure 4-78, where it can be observed that the mean stress at point A (within the sample) shows an important reduction, especially for the first 200 cycles. Sample CH3 presents the higher reduction in bulk modulus observed for the 5 horizontally drilled samples evaluated during the experimental program, a reduction of 29% from the initial value is obtained after more than 1000 cycles of temperature. Additionally, this sample presented the higher value of bulk modulus for all the samples drilled from the cliff, see section 3.5.2.

Figure 4-80 shows the evolution of the bulk modulus obtained with the numerical models considered for this sample in comparison with the experimental data measured. It can be remarked that for the first cycles, the numerical model results present some differences with the measurements, still the total variation in stiffness is adequately reproduced by the numerical models, note that after cycle 600 the data are follows closely by the modeling results. In absence of strains measurements, the thermal expansion coefficients considered for this sample are defined according to the deformation obtained for other horizontally drilled samples (CH2 and CH4). Thus, the thermal expansion considered in axial direction is higher (3 times) than the radial one.

This anisotropy generates, as observed for all samples drilled from the face of the cliff, that the bulk modulus reduces the most in the center of the sample (point A) than on the boundary (Point B), see Figure 4-80. Moreover, as for the GP model the thermal expansion coefficient considered is equal to a third part of the volumetric thermal expansion, less internal stresses are generated, leading the stiffness to be less affected by the thermal cycles.

From Figure 4-80, it can be observed that the stiffness of the material presents an important reduction for the first 200 cycles at point A, while the evolution of this parameter at point B is softer. This characteristic is linked to the evolution of the damage multiplier at these two points, see Figure 4-81. It can be remarked that the damage multiplier evolution for point A presents a more steeper response than the one observed at point B. As mentioned before, during the experimental program no strain gauges were installed on sample CH3, thus there are no strains measurements on this sample. Nonetheless, the axial and radial strains obtained with the axisymmetric model are presented in Figure 4-82, were it can be observed that a slightly strains accumulation is obtained with the imposition of thermal cycles, for sample CH3 an accumulation of -6.57 e-5 m/m and -4.85 e-5 m/m is obtained in radial and axial direction after 1000 thermal cycles. Moreover, the axial strains are higher than the radial ones, this is related to the thermal expansion coefficient values considered for this sample.

It is important to point out that all samples considering in the modeling shows a negative strains accumulation in both directions considered, this is to say compression strains, moreover even if small, the numerical modeling of samples obtained from the cliff shows almost 10 times higher accumulation than the modeling of samples obtained from blocks, which can be related to higher damage.

Conclusions for the numerical modeling

A composite material constitutive law has been proposed with the aim of reproducing thermal damage from a macroscopical point of view. The principal objective is to be able to reproduce the effect of thermal cycles on the mechanical response of the material. The numerical model presented accomplishes this objective.

Some specific observations are presented below:

• Thermal damage is interpreted as due to interactions between the thermal expansion of the different minerals that compose the rocks. This effect can be adequately reproduced by the composite constitutive law as this model is based on the modeling interactions between two different components provided with their own stresses, strains and constitutive laws.

• The numerical model proposed is a modification from the constitutive law proposed by [START_REF] Vaunat | Bond degradation and irreversible strains in soft argillaceous rock[END_REF]. The principal change consisted in modifying the partition of external stresses and strains between the two material components; the original constitutive law considers a system in serial while the modification performed in this work considers a parallel system. It generates some important changes in the mathematical formulation that controls the response or the composite material.

• The theoretical results presented by Hashin & Shtrikman (1963) to define the bounds for the bulk modulus of composite materials has shown to be provide an adequate range to the measurements of La Roque Gageac limestone. The stiffness values measured in samples drilled from blocks and directly from the face of the cliff belong to the range defined with this method and considering typical stiffness values for calcite and quartz (two main components of this rock).

• The coefficients Xm and Xb are indicators of the microstructural arrangement of the rock, including mineral content and stiffness and grain arrangement, which controls in fine the stiffness of the composite material.

• The undamaged stiffness obtained experimentally on LaRG limestone samples presents a wide range of variations, between 7000 MPa and 26000 MPa. Taken into account that these samples are mainly composed by calcite and quartz, the constitutive model is able to reproduce this wide range of stiffness values, considering fixed values of the mineral stiffness (Km and Kb) and different microstructural features, mineralogical content and geometrical arrangement.

• The coefficients Xm and Xb have an important effect on the stiffness of the composite material, its initial value and evolution. If Xb is higher than Xm then the response of the composite material will be mainly controlled by the bond response, otherwise the stiffness of the composite material is governed by the matrix response.

• The bulk modulus bounds defined by Hashin & Shtrikman (1963) are used to restrict the microstructural coefficients Xm and Xb.

• The composite constitutive law has been implemented in finite element method code (CODE_BRIGHT). It has also been explicitly integrated within an algebraic expression with the aim to provide an verification for model implementation in the code. Results have shown a good agreement, validating the model implementation.

• Variations in temperature generate internal stresses (matrix and bond stress). When bonding stresses exceed the elastic limit defined by the loading function, the bond suffers damage and the loading function increases homothetically.

• If no bond rate dependency is considered and thermal cycles are imposed, almost all the damage will be developed for the first cycle of temperature. However, this characteristic does not correspond to what is observed in the experimental program, moreover, the thermal damage is a fracture propagation phenomenon. Thus, the response is time dependent. In order to reproduce this phenomenon, rate dependency has been introduced in the bond damage model

• The numerical model is used to reproduce the response obtained during the experimental program on 6 different samples: 2 obtained from blocks, 2 vertically drilled and 2 horizontally drilled from the face of the cliff. This modeling task is developed with two different approaches. In first place, each sample is modelled as a Gauss Point. Afterwards sample have been simulated as a boundary problem with an axisymmetric geometry. Both approaches give a good approximation of the response observed in the samples.

• For the numerical modeling of block samples, coefficient Xb and bond fraction has been taken higher than coefficient Xm and matrix fraction. Therefore, the evolution of the composite material bulk modulus is controlled by the bond bulk modulus. On the other hand, for the samples cored from the face of the cliff, values of Xm is taken higher than Xb, and matrix and bond fractions are the same. In that case, the composite material bulk modulus is governed by the matrix stiffness.

• Samples obtained from blocks present more damage at the boundary of the sample (higher damage multiplier and lower bulk modulus) while samples obtained from the cliff developed more damage at the center of the sample. It is because the thermal expansion coefficient considered for the block samples are very similar in axial and radial direction, while the thermal expansion coefficient obtained in axial direction is almost 3 times the one considered in radial direction for the cliff samples. This generates higher internal stresses in the center of the cliff samples and finally a higher damage multiplier.

• The damage multiplier obtained after more than 1000 thermal cycles is much higher for samples obtained from the cliff face (1.5) than for the block samples (0.13). This characteristic is related to the initial bulk modulus measured in the samples. It is higher in the block samples and experiences a lower reduction than in the cliff samples. In fact, higher values of final damage multiplier are linked to higher reductions in the bulk modulus.

• For all the tested samples, it has been observed a link between the evolution of damage multiplier thought thermal cycles and the bulk modulus.

• The results obtained with the numerical modeling of samples show a slight strains accumulation, contrary to what is observed in the experimental program results, where high values of strains accumulation are developed. Nevertheless, it is important to mention that with the numerical modeling the strains accumulate in negative direction (compression).

• The numerical model is able to reproduce the strains amplitude registered during a 24 hours interval (6 cycles for block samples and 6.5 for cliff samples). The stiffness variations computed by the numerical model are mainly related to stress accumulation.

• The modified composite material constitutive law has generally a good capability to reproduce the behavior of LaRG samples experimentally observed. It is an important tool to understand and explain the different experimental patterns observed in the samples. The phenomenon of thermal damage is evaluated from the case of La Roque Gageac cliff, a site that has suffered several rockfalls. The last event occurred in 2010, leads to the installation of an instrumentation system that follows displacements and temperature at the cliff, during more than 6 years. From these measurements, it is observed that there is a strong relationship between the deformations in the cliff and the temperature.

Conclusions

The deformations registered with the extensometers in the site are not reversible, in other words, there is an accumulation of deformations through the time (thermal cycles). In the same way, the jointmeters also show irreversible opening of the cracks with the thermal cycles. It is important to mention that according to the jointmeters measurements, the opening of the fissures is related to cold periods, explaining why most of the rockfalls observed on the site are registered during the winter.

The instrumentation observations show that temperature is a controlling factor on the deformation of the cliff. Based on this, an experimental program is developed on limestone samples from La Roque Gageac, showing that the temperature rates registered in the site are able to induce damage in this limestone and to generate the propagation of existent fissures.

Because of the heterogeneity of the limestone of the cliff, samples were extracted from different locations (within the cavern and at the cliff face). It should be remarked that differences on the mineralogical composition are observed regarding the original location of the samples, furthermore, this heterogeneity also generates differences in the response to the thermal cycles observed.

The behavior of the samples through thermal cycles is evaluated with the evolution of three principal measurements: deformations, elastic wave propagation velocity and uniaxial compressive strength. All these parameters show an evolution with the applications of thermal cycles.

In general, the samples present some strains accumulation. However, the behavior observed don't give any response that can be linked to point of extraction of the samples and/or the mineralogical composition or the sample direction regarding the cliff face. Except by the fact that samples obtained from the cliff face show highest thermal expansion in the perpendicular direction of the cliff face, which is the direction where the cliff suffered major decompression during the valley formation.

Furthermore, the elastic wave propagation velocity has shown to be a good parameter to evaluate the damage in rocks, because of its sensitiveness to the presence of internal microfissures, moreover, as it is a nondestructive measurement, its evolution can be followed for the whole experimental program. It is observed that thermal cycles generate a reduction in the elastic wave propagation velocities, especially evident in P waves (VP), that might be related to the presence of internal micro-fissures because any visible crack was observed in the face of the samples.

From the elastic wave propagation velocities, the dynamic strength parameters of the rock can be estimated. Samples from La Roque gageac shows a reduction in the dynamic bulk modulus is observed through the imposition of thermal cycles. Moreover, from VP values the uniaxial compression strength (UCS) can be estimated for samples obtained from inside the troglodyte cavern. The UCS also presents a reduction with the imposition of thermal cycles.

It was observed that atmospheric thermal cycles are not capable to generate any new macro-fissure in the limestone from La Roque Gageac, at least at photo scale. Nevertheless, the propagation of existent fissures is well identified.

The thermal damage on rocks is a topic of interest for several authors, in fact, many experimental programs have been focused on the effect of thermal damage in rocks from different lithologies. However, these studies are generally focused on high temperatures, such as those observed during fires or in underground waste disposal. The experimental program presented in this work shows that atmospheric thermal rates can also generate damage in the rock matrix.

From the results obtained in the experimental program, the phenomenological characteristics of the thermal damage in rocks are defined. They are used to establish the bases considered to proposed a constitutive law capable to reproduce the effect of thermal cycles on the mechanical response of the limestones from La Roque Gageac. With this purpose, a composite material constitutive law is prosed.

From a phenomenological point of view the thermal damage is related to the interaction between the thermal response of the minerals that composed the rock. This characteristic is taken into account by the numerical model as it considers that the material (rock) is composed of two materials and the interaction between them define the response of the composed material.

The constitutive model proposed considers parameters that characterized the microstructural arrangement as one of the principal characteristics that defined the response of the composite material. In fact, by considering typical strength values for quartz and calcite (the main minerals conforming La Roque gageac limestone) it is possible to reproduce all the varieties of responses observed in the experimental program, only by changing these microstructural parameters.

The composite material constitutive law has the ability to reproduce the behavior of the samples from LaRG evaluated in the experimental program. It is an important tool to understand and explain the different experimental patterns observed in the samples.

The experimental study and numerical model presented in this work shows that atmospheric thermal cycles appreciated in the region of La Rock Gageac may induce a damage in the rock matrix reducing its strength and helping to the propagation of the existent fissures, which confirm the idea that thermal cycles may be the triggering factor for the instabilities registered on the site.

Future development.

Based on the experimental results obtained during this thesis and the numerical model proposed some lines of future work are open, some of them are proposed bellow:

• From insitu observations it is observed that the water doesn't have a controlling role in the instabilities registered in the cliff of La Roque Gageac, in the same way, the experimental program showed that by isolating the thermal effect a considerable damage can be observed in the material. Nevertheless, it must be remembered that this cliff is subjected to the actions of the whole climatological conditions, for this reason it can be interesting to evaluate if the relative humidity plays a role in the thermal damage of this limestone.

• During the experimental program, it was observed that any new fracture was generated because of the thermal cycling, nevertheless, if there is a pre-existent crack these thermal cycles induce its propagation, that can be related to the reduction in strength of the rock mass. In order to have a quantitative evaluation of how thermal damage affect the propagation of fissures, it will be necessary to performed toughness tests on samples that were subjected to thermal cycles.

• The evaluation of thermal damage performed in this thesis is focused on the limestone from La Roque Gageac, it could be interesting to evaluate if thermal cycles produce the same effect with similar phenomenology in other types of lithologies, moreover, if the numerical model proposed in this work is also able to reproduce the response observed in those other lithologies.

• The numerical model proposed is capable to reproduce the effect of strength and stiffness reduction in the rock due to the thermal cycles, nevertheless the phenomena of propagation of fractures is not directly considered. This phenomenon could be included in the numerical modeling. In the finite elements scheme, this can be performed by two strategies. One option is to represent the existent crack as a weak zone and its propagation as the localization of damage through a non-local formulation [START_REF] Pijaudier | Non local Damage Theory[END_REF][START_REF] Mazars | From damage to fracture mechanics and conversely: A combined approach[END_REF][START_REF] Jirásek | Nonlocal models for damage and fracture: Comparison of approaches[END_REF]. Another options is to define the extent crack as a joint element and discretize the continuum medium with zero-thickness elements [START_REF] Carol | Normal/Shear cracking model: Application to discrete crack analysis[END_REF][START_REF] Segura | On zero-thickness interface elements for diffusion problems[END_REF]) thus the propagation of this fissure will be achieved by the opening and relative displacement of those elements

• Finally, it is important to be able to model the real case of the cliff of La Roque Gageac. This should be performed focusing in the joints response, as it is where the instrumentation system is installed, for this it will be necessary to have the crack propagation and response included in the model. Whit the numerical model of the case it is expected to estimate the real risk level. The increment of damage multiplier can then be obtained as:

Figure 1 -Figure 2 - 1 Figure 2 -Figure 2 - 3 Figure 2 - 4 Figure 2 - 5 Figure 2 - 7

 121223242527 Figure 1-1 Landslide-sequencing model (from Julian & Anthony (1996)) ……………………………... 3 Figure 2-1 Location of la Roque Gageac A) Position of La Dordogne department inside of France B) Situation of La Roque Gageac in the department ………………………………………………….. 12 Figure 2-2 situation of La Roque Gageac in relation to the Dordogne river and the cliff …….. 12 Figure 2-3 Geological situation of La Roque Gageac (obtained from the geological map from Sarlat La-Canéda -(BRGM, 1987)) …………………………………………………………………………….. 13 Figure 2-4 Example of a hydraulic deposit with different shapes and surfaces ……………………. 14 Figure 2-5 Location of principal events observed in the town (Virely and Guittard, 2010a) … 16 Figure 2-6 Source of rockfalls registered in 1920 (A) and in 1957 (B) ………………………………….. 16 Figure 2-7 Details of the collapse occurred in 2010 (from Ruiz, 2013) a) zone of failure b) view inside the cavern ………………………………………………………………………………………………………. 17

Figure 2 - 8

 28 Figure 2-8 First level of risk A) Remaining blocks inside the cavern B) Wire mesh installed ….18

Figure 2 - 9 Figure 2 -Figure 2 -Figure 2 -Figure 2 -Figure 2 -Figure 2 -Figure 2 -Figure 2 -Figure 2 -Figure 2 -

 292222222222 Figure 2-9 Flacking observed in the face of the cliff ……………………………………………………………..19 Figure 2-10 Stability condition of the beam A) General situation B) West support C) East support …………………………………………………………………………………………………………………….. 19 Figure 2-11 Principal discontinuities families observed (Ruiz, 2013) …………………………………….20 Figure 2-12 Supporting structures inside the cavern …………………………………………………………... 20 Figure 2-13 Meteorological station data …………………………………………………………………………….. 21 Figure 2-14 Cross-section scheme of the cliff with the localization of the extensometers (D1 and D2) …………………………………………………………………………………………………………………….. 22 Figure 2-15 Localization of the instrumentation devices (from Gasc-Babier, Virely and Guittard, 2015) ……………………………………………………………………………………………………………………….. 22 Figure 2-16 Temperature evolution at 2 and 6 meters from the Cliff face …………………………...23 Figure 2-17 Extensometers measures: a) D1 b) D2 ……………………………………………………………... 24 Figure 2-18 Possible discontinuity across the extensometers ……………………………………………. 25 Figure 2-19 Jointmeters displacement and temperature inside the cavern ………………………….25

Figure 2 -Figure 2 -Figure 2 -Figure 2 -

 2222 Figure 2-20 Joint behavior under low temperatures a) Horizontal crack b) Vertical crack (Ruiz Restrepo 2013) …………………………………………………………………………………………………………. 26

Figure 3 -

 3 Figure 3-49 Temperature -deformation evolution sample CH5 -Vertical deformation ……..

Figure 3 -

 3 Figure 3-50 Vertical direction accumulated deformation -samples group 1 ……………………….

Figure 3 -

 3 Figure 3-51 Radial direction accumulated deformation -samples group 1 ………………………….

Figure 3 -

 3 Figure 3-52 Vertical direction accumulated deformation -samples group 2 ……………………….

Figure 3 -

 3 Figure 3-53 Radial direction accumulated deformation -samples group 2 ………………………….

Figure 3 -

 3 Figure 3-54 Compressive elastic wave velocity evolution samples group 1 ………………………….

Figure 3 -

 3 Figure 3-55 Compression elastic wave propagation velocity evolution-samples group 2: a) samples drilled in vertical direction b) samples drilled in horizontal direction …………...

Figure 3 -

 3 Figure 3-56 Rock cores obtained in horizontal direction: a) Core 10 -samples CH1, CH2 and CH3 b) Core 11 -Samples CH4 and CH5 …………………………………………………………………….

Figure 3 -

 3 Figure 3-57 Shear elastic wave velocity evolution -samples group 1 …………………………………..

Figure 3 -

 3 Figure 3-58 Shear elastic wave propagation velocity evolution-samples group 2: a) samples drilled in vertical direction b) samples drilled in horizontal direction ………………………...

Figure 3 -

 3 Figure 3-59 Bulk modulus evolution -samples group 1 ……………………………………………………….

Figure 3 -

 3 Figure 3-60 Bulk modulus evolution -samples group 2: a) samples drilled in vertical direction b) samples drilled in horizontal direction …………………………………………………………………..

Figure 3 -

 3 Figure 3-61 Deformation -stress evolution ………………………………………………………………………….83 Figure 3-62 UCS and VP evolution ……………………………………………………………………………………….

Figure 3 -

 3 Figure 3-63 Young modulus evolution through thermal cycles …………………………………………….

Figure 3 -

 3 Figure 3-64 Relation between Est and Edyn ……………………………………………………………………... 86 Figure 3-65 Initial Vp and porosity relation ………………………………………………………………………….

Figure 3 -

 3 Figure 3-66 Bulk modulus variation (Calculated from Vp-Vs) with initial porosity (n) …………..

Figure 3 -

 3 Figure 3-67 deformations and VP evolution -sample B3 ……………………………………………………..

Figure 3 -

 3 Figure 3-68 deformations and VP evolution -sample CV4 ………………………………………………... 89 Figure 3-69 deformations and VP evolution -sample CH2 ………………………………………………….

Figure 3 -

 3 Figure 3-70 VP -UCS relation ……………………………………………………………………………………………….

Figure 3 -

 3 Figure 3-71 a) sample 1345, with speckle pattern b) sample 1346 without speckle pattern ..

Figure 3 -Figure 4 - 1 Figure 4 - 2 Figure 4 - 3 Figure 4 - 4 Figure 4 - 5 Figure 4 - 9 Figure 4 -Figure 4 -Figure 4 -Figure 4 -Figure 4 -Figure 4 -Figure 4 -Figure 4 -Figure 4 -Figure 4 -Figure 4 -Figure 4 -Figure 4 -Figure 4 -Figure 4 -Figure 4 -Figure 4 -Figure 4 -Figure 4 -Figure 4 -Figure 4 -Figure 4 -Figure 4 -Figure 4 -Figure 4 -Figure 4 -Figure 4 -Figure 4 -Figure 4 -Figure 4 -Figure 4 -Figure 4 -Figure 4 -Figure 4 -Figure 4 -Figure 4 -

 3414243444549444444444444444444444444444444444444 Figure 3-72 Digital image correlation results -sample 1345 ………………………………………………..

Figure 4 -

 4 Figure 4-53 Sample B3: Axial strains ………………………………………………………………………………...

Figure 4 -

 4 Figure 4-54 Sample B3: Radial strains ………………………………………………………………………………..

Figure 4 -

 4 Figure 4-55 Sample B3: Mean stress evolution ………………………………………………………………….

Figure 4 -

 4 Figure 4-56 Sample B3: Strains amplitude a) radial strains b) axial strains ……………………….

Figure 4 -

 4 Figure 4-57 Sample B7: Bulk modulus evolution ………………………………………………………………..

Figure 4 -

 4 Figure 4-58 Sample B7: Damage multiplier evolution ………………………………………………………..

Figure 4 -

 4 Figure 4-59 Sample B7: Damage multiplier evolution in the sample ………………………………....

Figure 4 -

 4 Figure 4-60 Sample B7: Axial and radial strains evolution ………………………………………………….

Figure 4 -

 4 Figure 4-61 Sample CV4: Bulk modulus evolution ……………………………………………………………..

Figure 4 -

 4 Figure 4-62 Sample CV4: Bulk modulus evolution in the sample ………………………………………..

Figure 4 -

 4 Figure 4-63 Sample CV4: Damage multiplier evolution ……………………………………………………..

Figure 4 -

 4 Figure 4-64 Sample CV4: Damage multiplier evolution in the sample ………………………………..

Figure 4 -

 4 Figure 4-65 Sample CV4: Axial strains evolution ………………………………………………………………..

Figure 4 -

 4 Figure 4-66 Sample Cv4: Radial strains evolution ……………………………………………………………...

Figure 4 -

 4 Figure 4-67 Sample CV4: Strains amplitude a) radial strains b) axial strains ……………………….

Figure 4 -

 4 Figure 4-68 Sample CV4: Mean stress evolution ……………………………………………………………….

Figure 4 -

 4 Figure 4-69 Sample CV5: Bulk modulus evolution ……….…………………………………………………….

Figure 4 -Figure 4 -

 44 Figure 4-70 Sample CV5: Damage multiplier evolution ………………………………………………………168 Figure 4-71 Sample CV5: Axial and radial strains evolution ……………………………………………….

Figure 4 -

 4 Figure 4-72 Sample CH2: Bulk modulus evolution …………………………………………………………….

Figure 4 -

 4 Figure 4-73 Sample CH2: Bulk modulus evolution in the sample ……………………………………….

Figure 4 -

 4 Figure 4-74 Sample CH2: Damage multiplier evolution …………………………………………………….

Figure 4 -Figure 4 -

 44 Figure 4-75 Sample CH2: Damage multiplier evolution in the sample ……………………………… Figure 4-76 Sample CH2: Radial strains evolution ……….…………………………………………………….

Figure 4 -

 4 Figure 4-77 Sample CH2: Axial strains evolution ……………………………………………………………….

Figure 4 -

 4 Figure 4-78 Sample CH2: Mean stress evolution ……………………………………………………………….

Figure 4 -

 4 Figure 4-79 Sample CH2: Strains amplitude a) radial strains b) axial strains ……………………….

Figure 4 -

 4 Figure 4-80 Sample CH3: Bulk modulus evolution …………………………………………………………….

Figure 4 -

 4 Figure 4-81 Sample CH3: Damage multiplier evolution …………………………………………………….

Figure 4 -

 4 Figure 4-82 Sample CH3: Strains evolution …………………………………………………………………….… Figure A-1. Pore size distribution for samples submitted to thermal cycles ………………………..

Figure A- 3

 3 Figure A-3 DIC analysis for sample C1 -face without speckle pattern ………………………………..

Figure A- 4

 4 Figure A-4 DIC analysis for sample C3 -face with speckle pattern ……………………………………..

Figure A- 5

 5 Figure A-5 DIC analysis for sample C3 -face without speckle pattern ………………………………..

  1 Density and elastic wave propagation velocity ........................................................ Table 3-2 Porosity for La Roque Gageac Limestone ................................................................ Table 3-3 Mechanical properties LaRG limestone ................................................................... Table 3-4 RMR evaluation ........................................................................................................ Table 3-5 Thermal expansion results ....................................................................................... Table 3-6 Geometry and characteristics cylindrical samples (CV: samples drilled in vertical sense from the face of the cliff, CH: Samples drilled in horizontal sense from the face of the cliff) ........................................................................................................................... Table 3-7 Geometry and characteristics prismatic samples .................................................... Table 3-8 Mineralogical analysis samples ................................................................................ Table 3-9 Number of cycles imposed and test performed for each sample ........................... table 3 10 Accumulated axial strains ....................................................................................... Table 3 11 Accumulated radial strains ..................................................................................... Table 3 12 Mean deformations amplitude for thermal-strains evolution ...............................

  ..................................................................................... Table 4 1 Parameters considered for the numerical model analysis ..................................... Table 4 2 Parameters used for the numerical modeling of samples ..................................... Table 4 3 Initial bulk modulus values (MPa) .......................................................................... 1. Introduction 1.1. General aspects. .......................................................................................................... 3 1.2. Aims of the research and thesis layout. ...................................................................... 7 1.1. General aspects.

Figure 1 -

 1 Figure 1-1 Landslide-sequencing model (from Julian & Anthony (1996))

Figure 2 - 1

 21 Figure 2-1 Location of la Roque Gageac a) Position of La Dordogne department inside of France b) Situation of La Roque Gageac in the departmentThis town is located at the right bank of the Dordogne river, in a regional area denominated The Périgord noir (Figure2-1b), at the toe of a limestone rock-cliff of more or less 100 meters high. In this cliff, there is an ancient troglodyte fort located at 50 m high (see Figure2-2). This cavern was used as a refuge during the hundred Years' War, which makes it one of the tourist attraction for the region.

Figure 2 -

 2 Figure 2-2 situation of La Roque Gageac in relation to the Dordogne river and the cliff

Figure 2 - 3

 23 Figure 2-3 Geological situation of La Roque Gageac (obtained from the geological map from Sarlat La-Canéda -(BRGM, 1987))

Figure 2 - 4

 24 Figure 2-4 Example of a hydraulic deposit with different shapes and surfaces

Figure 2 -Figure 2 - 5 Figure 2 - 6 Figure 2 - 7

 2252627 Figure 2-5 shows the location of three of the four principal events registered in the zone of interest: 1920: A block of several thousand of m 3 fell in the town from the east side of the cliff, without casualties. (Figure 2-6a) 1957: At the west of the village a large block fell and kill 3 people (Figure 2-6b) 1994: Several small blocks were detached from the cliff 2010: Collapse of the roof of the cavern (Figure 2-7)

Figure 2 Figure 2 - 8 Figure 2 -

 2282 Figure 2-8 First level of risk a) Remaining blocks inside the cavern b) Wire mesh installed

Figure 2 -

 2 Figure 2-12 Supporting structures inside the cavern

Figure 2 -

 2 Figure 2-13 Meteorological station data

Figure 2 -Figure 2 -

 22 Figure 2-14 Cross-section scheme of the cliff with the localization of the extensometers (D1 and D2)

Figure 2 -

 2 Figure 2-16 Temperature evolution at 2 and 6 meters from the Cliff face An extensometer measures a relative displacement between the points of interest (2, 4 and 6 m depth from the cliff face) and the end of the device, which is consider fixe. The negative displacements values observed represent a reduction in this distance, associated to the contraction of the rock mass, due to the temperature decrease. Hence, the positive values denote an expansion of the cliff and negatives the compression.

Figure 2 -

 2 Figure 2-17 Extensometers measures: a) D1 b) D2

Figure 2 -

 2 Figure 2-18 Possible discontinuity across the extensometers. The response recorded by the jointmeters in comparison with the external temperature is presented in Figure 2-19, Likewise with the extensometers, these devices results show a direct relation between deformations and temperature variations.

Figure 2 -

 2 Figure 2-19 Jointmeters displacement and temperature inside the cavern Nevertheless, the results obtain with the jointmeters have a different meaning of those obtained with the extensometers. While the extensometers, measure the deformation of the complete rock mass, the jointmeters measures the aperture of a discontinuity delimited by two different rock blocks.

Figure 2 -

 2 Figure 2-20 Joint behavior under low temperatures a) Horizontal crack b) Vertical crack[START_REF] Ruiz | Thermo-mechanical analysis of the stability of a rock-cliff under climatic actions[END_REF]

Figure 2 -

 2 Figure 2-21 Displacements vs temperature J-5

Figure 2 -

 2 Figure 2-22 displacement obtained with Jointmeters J-1 and J-2

Figure 2 -

 2 Figure 2-23 presents the displacements versus temperature evolution observed with the device J-2. It can be observed that in the period that includes February 2012 the cumulated displacements are much higher to what is observed in the other three periods evaluated when no extreme low temperatures where registered.

Figure 2 -

 2 Figure 2-23 displacements vs temperature J-2

Figure 3 - 1

 31 Figure 3-1 Discontinuity in sample 1018 Besides the mechanical tests performed within the framework of the DO-SMS project, some additional Brazilian tensile strength tests were performed on samples obtained from a third block (Block C) fell during the 2010's event and removed in February 2015. Including these last tests, the mean value of tensile strength is 5.39 MPa with a standard deviation of 0.50. The small standard deviation value obtained can be related to the probabilistic density associated with this type of test and the volume of the specimen associated, as illustrated in Figure 3-2.

Figure 3 - 2

 32 Figure 3-2 Tensile strength variation as function of type of test and specimen volume from (Hudson & Harrison 1997) Finally, triaxial test were performed over 3 samples under 2.5, 5 and 10 MPa confining pressure. The results obtained are presented on Figure 3-3 and Table 3-3. The response observed with these triaxial tests includes the transition between a fragile behavior to a ductile one, associated to the level of confining pressure (Hudson & Harrison 1997) . In fact, in Figure 3-3 it can be observed that samples under lower confining pressure show a fragile behavior, while sample under 10 MPa of confining pressure has a ductile response.

Figure 3 - 3

 33 Figure 3-3 Triaxial test results

Figure 3 - 4

 34 Figure 3-4 Samples obtained from blocks-February 2015The second set consist of samples directly drilled from the cliff face, nearby the troglodyte cavern, as illustrated in Figure3-5a. For this second group of samples, two types of specimens are extracted, half of the samples were drilled in vertical sense (Figure3-5b) and the other half in horizontal sense (Figure3-5c). It is important to mention that around the troglodyte cavern the cliff shows an almost horizontal stratification (see Figure3-5a), therefore horizontally drilled samples are parallel to the stratification planes while vertical samples are perpendicular to it.

Figure 3 - 5

 35 Figure 3-5 Sample extraction from the face of the Cliff-April 2016: a) location b) Samples obtained in vertical sense c) samples obtained in horizontal sense

Figure 3 - 7

 37 Figure 3-7 Infrared spectra for simple M1, M2 and M3, compared with typical calcite bands X-Ray diffraction.As from spectrometry test was not possible to determine the mineralogical composition of samples from the face of the cliff, X-Ray diffraction (XDR) analysis were performed on samples M4 and M5 belonging to group 2. It is important to mention that sample M4 was obtained from a sample drilled in the vertical direction and sample M5 in the horizontal direction both close to the face of the cliff. The XDR technique allows to obtain the mineralogical composition of the sample by collecting the X-ray reflected by the powder sample. These diffraction peaks are converted to d-spacing allowing to identify the minerals components because each mineral has a set of d-spacing characteristic. XDR results for samples M4 and M5 are shown in Figure 3-8 and Figure 3-9 respectively. It can be observed that both samples are composed almost in a 50-50 relation by quartz and calcite. More precisely sample M4 is composed in a 55% by quartz and 42% by calcite and sample M5 in a 47% quartz and 48% calcite.

Figure 3 -

 3 Figure 3-8 Rx-diffraction results for sample M4

FigureFigure 3 -

 3 Figure 3-10 Thin sections images under natural light

Figure 3 -

 3 Figure 3-11 Mineral identification sample M1-B

Figure 3 -

 3 Figure 3-12 Holes observed in thin section of sample M3

Figure 3 -

 3 Figure 3-13 Modal composition computed from image analyses of the thin sections.

Figure 3 -

 3 Figure 3-14 Shape ratio distribution for calcite type micrite

Figure 3 -

 3 Figure 3-15 Shape ratio distribution for calcite type sparite

Figure 3 -

 3 Figure 3-16 Shape ratio distribution for quartz

Figure 3 -

 3 Figure 3-17 crystals fissure density

Figure 3 -

 3 Figure 3-18 Climatic chamber employed The typical thermal cycle imposed lasts 3.5 hours and consists in: 45 minutes at 10ºC, then one hour to go from 10ºC to 50ºC, later, 45 minutes at 50ºC and finally, 1 hour to go from 50ºC to 10ºC, as illustrated in Figure 3-19. The external measurements (Vp-Vs and DIC) are performed each day during an interval of tree hours with a constant temperature of 10ºC.

Figure 3 -

 3 Figure 3-19 Thermal cycle imposed

  Figure 3-20 Thermal probe on sample B4: a) Sample face b) within the sampleThe measured temperatures are shown in Figure3-21. The temperature on the face and within the samples varies with a rate of 0.35ºC/min, much less than the rate of 1ºC/min, defined by[START_REF] Hall | Thermal fatigue and thermal shock in bedrock: An attempt to unravel the geomorphic processes and products[END_REF] as the threshold for destructive internal effects for thermal shock. In the same way, the maximal temperature difference between the face and the center of the sample is 2ºC.

Figure 3 -

 3 Figure 3-22 Thermal conductivity back analysis result

Figure 3 -

 3 Figure 3-23 Scheme strain gauges installation Measurements are performed with an amplifier type MJCplus, connected to a personal computer equipped with a data acquisition software. The measurements of voltage variation (mV) are registered each five minutes, during the whole experimental program. Based on the gauge characteristics (gage factor and nominal resistance) and voltage measurements deformation variations are calculated. As an example, Figure 3-24 presents the strain gauges glued on samples B1 and sample CH5.

Figure 3 -

 3 Figure 3-25 Elastic wave propagation velocities measurement equipment.

  Figure 3-26 Typical arrivals elastic wave time a) P-waves b) S-waves

Figure 3 -

 3 Figure 3-27 MTS rock test system 816

Figure 3 -

 3 Figure 3-28 Fissures observed in a former experimental study, Gasc-Barbier et al. (2014)

  -29). It allows to put the sample in the same position each time, having the camera immobile during the whole experimental program. To reduce the light noise, photos are taken in a dark room, employing just a fixed lamp as light source.

Figure 3 -

 3 Figure 3-29 Fixing device

Figure 3 -

 3 Figure 3-30 Vertical deformations for samples B1, B2 and B3After 948 cycles, sample B1 accumulates a vertical deformation of 6.8e-4 m/m, sample B3 behaves in the same way with an accumulated deformation of 8.5e-4 m/m. While, the accumulated vertical strains for sample B2 is only 3.6e-4 m/m but in the negative direction.

Figure 3 -

 3 Figure 3-31 Vertical deformation samples CV1, CV4, CH2 and CH5

Figure 3 -

 3 Figure 3-32 Radial deformation for samples B1, B2 and B3

  Figure 3-34 Radial deformation evolution for sample CH5: a) General tendency b) Daily response

  Figure 3-36 Daily deformations response for sample B2: a) vertical deformation b) Radial deformation

  -42. In order to compare the general behavior, only 6 typical cycles are considered: cycle 1, 50, 100, 200, 500 and 900. Results obtained for samples obtained from blocks are presented in Figure 3-43 to Figure 3-45, and for samples drilled from the face of the cliff in Figure 3-38 to Figure 3-41.

Figure

  Figure 3-42 Typical thermal cycle considered

3

  -50 and Figure3-51 present the accumulated strains through thermal cycles, in vertical and radial direction respectively for samples from group 1, and Figure3-52 Figure3-53, for samples obtained from the face of the cliff. It is important to mention that for these figures; a logarithmical scale is considered in the cycles axis only for visual purposes.

Figure 3 -

 3 Figure 3-50 Vertical direction accumulated deformation -samples group 1

Figure

  Figure 3-54 Compressive elastic wave velocity evolution samples group 1

  Figure 3-56 Rock cores obtained in horizontal direction: a) Core 10 -samples CH1, CH2 and CH3 b) Core 11 -Samples CH4 and CH5

Figure

  Figure 3-57 Shear elastic wave velocity evolution -samples group 1

Figure

  Figure 3-59 Bulk modulus evolution -samples group 1

Figure 3 -

 3 Figure 3-61, shows the deformation-stress evolution for samples that have known different number of thermal cycles. It is observed that samples submitted to a larger number of cycles, has a lower value of Uniaxial compressive strength, coupled with a lower bent in the strains-strength evolution.

Figure 3 -

 3 Figure 3-61 Deformation -stress evolutionThe variation of UCS through the thermal cycles, as the evolution of elastic wave propagation velocities (Vp final), is shown in Figure3-62. It can be observed a decrease in the elastic wave propagation velocity, coupled with a decrease in the rock strength.

Figure 3 -

 3 Figure 3-62 UCS and VP evolution

Figure 3 -

 3 Figure 3-63 Young modulus evolution through thermal cycles

Figure 3 -

 3 Figure 3-64 Relation between Est and Edyn

Figure 3 -

 3 Figure 3-65 Initial Vp and porosity relation Moreover, from Figure 3-66 it can be observed that samples obtained from blocks with higher porosity develop more damage (evaluated as the reduction of Vp) associated with thermal cycling. Vp values evolution measured on samples drilled from the face of the cliff decrease more for samples with lower porosity.

Figure 3 -

 3 Figure 3-66 Bulk modulus variation (Calculated from Vp-Vs) with initial porosity (n)

Figure 3 -

 3 Figure 3-67, presents the strains and Vp evolution of sample B3. It can be observed that vertical and radial strains increase for the whole experimental program, but Vp remains constant, nevertheless, despite this positive tendency vertical strains accumulate almost no deformation after 900 cycles.

Figure 3 -

 3 Figure 3-67 deformations and VP evolution -sample B3

Figure 3 -

 3 Figure 3-68 deformations and VP evolution -sample CV4

Figure 3 -

 3 Figure 3-69 deformations and VP evolution -sample CH2According to[START_REF] Najibi | Empirical relations between strength and static and dynamic elastic properties of Asmari and Sarvak limestones, two main oil reservoirs in Iran[END_REF] it is possible to obtain an empirical potential relationship between VP and Uniaxial compressive strength since the porosity and density has a direct effect in UCS. Having this in mind for samples mainly composed by calcite the uniaxial compressive strength can be predicted from the Vp value by the equation Eq. 3-12, where Vp

  Figure 3-71 a) sample 1345, with speckle pattern b) sample 1346 without speckle pattern All samples considered (1343 to 1348) presents similar results. To illustrates the results obtained with DIC, the analysis obtained for samples 1345 (with the paint treatment) and 1346 (without paint treatment) are presented in Figure 3-72 and Figure 3-73, respectively. It can be

Figure 3 -

 3 Figure 3-72 Digital image correlation results -sample 1345

Figure 3 -

 3 Figure 3-73 Digital image correlation results sample 1346

Figure 3 -

 3 Figure3-74 and Figure3-75 present the DIC results obtained for sample C2 for two different faces of the sample, one with a speckle pattern done with black and white paint (Figure3-74) and the other one under natural conditions (Figure3-75). Note that no identifiable crack can be observed from these analyses, similar results are observed by[START_REF] Leiss | Fabric anisotropy and its influence on physical weathering of different types of Carrara marbles[END_REF] on different types of marbles, subjected to heating until 50ºC.

  -74) than for the face under natural conditions (Figure 3-75). Similar results are obtained in sample C1 and C3, presented in Appendix B.

Figure 3 -

 3 Figure 3-74 DIC analysis for sample C2 -face with speckle pattern

Figure 3 -

 3 Figure 3-76 Fissure generated on samples B20.

Figure 3 -

 3 Figure 3-77 Fissure evolution for sample B20

Figure 3 -Figure 3 -

 33 Figure 3-79 Sample C4 Unlike the case of sample B20, is not possible to visually identify the propagation of fissures in samples C4, except by the fact that between test day 150 (cycles 972 to 978.5) and 155 (cycles 1004.5 to 1011) a small block was detached from the top of the sample, as is presented in Figure 3-80.

Figure 3 -

 3 Figure 3-81 Sample C4 condition at 730 and 820 cycles.

Figure 3 -

 3 Figure 3-82 DIC analysis for sample C4-face1 from 0 to 293 cycles

Figure 3 -

 3 Figure 3-84 DIC analysis for sample C4-face2 from 0 to 290 cycles

  These internal stresses may lead to local material weakening, fissuration or micro-cracking. As a matter of fact, Figure4-1 presents photomicrographs obtained by[START_REF] Malaga-Starzec | Microscopic and macroscopic characterization of the porosity of marble as a function of temperature and impregnation[END_REF] on a sample of marble before and after application of atmospheric temperature cycles. The marble is predominantly composed by calcite (>99%). After 50 thermal cycles, material experiences an important augmentation of porosity as the result of void development around grains boundary.

Figure 4 - 1

 41 Figure 4-1 Photomicrographs of calcite samples, CR0: Uncycled sample CO: Temperature cycled (Malaga-Starzecet al., 2006)

  1+e)b) The bond ratio eb gives the bond volume (Vb) over the volume of solid (Vs). It is related to the volumetric strain experienced by the bond by dε vb = The matrix void ratio em, gives the volume existing between matrix grains, that is the sum of void and bond volumes (Vv+Vb) over the volume of solid (Vs = Vm+Vb). The variation of em gives how close matrix grain go and is related to the volumetric strain of the composite material (dε vm ) by dε vm -

Figure 4 - 2

 42 Figure 4-2 Conceptual scheme of the composite material that represents an argillaceous rock, REB:Representative element boundary -From[START_REF] Vaunat | Bond degradation and irreversible strains in soft argillaceous rock[END_REF]

Figure 4 - 3

 43 Figure 4-3 New conceptual scheme for the composite material constitutive model

  a) Strains and stresses acting over the whole bond area labeled σ ij b and ε ij b . There are apparent quantities acting externally on the bonds. b) Strains and stresses acting over the unfissured bond area, labeled σ ̂ij b and ε ̂ij b . They are effective quantities, acting internally in the bonds.

Figure 4 -

 4 4 and Figure4-5 present the bounds obtained for the bulk modulus as a function of the porosity, for block and cliff drilled samples respectively. On the same figure, the bounds are compared with the values obtained from elastic wave propagation velocities (see section 3.5.2) measured on the tested samples. It can be observed that the experimental values lay between the Hashin & Shtrikman (1963) bounds. Therefore, this method appears to provide bounds for the bulk modulus that can be used for the rock under study.

Figure 4 - 4

 44 Figure 4-4 Hashin & Shtrikman bulk modulus bounds as a function of porosity for block samples

Figure 4 - 6

 46 Figure 4-6 Acceptable range of Xm-Xb coefficients for block samples with porosity equal to 0.2

Figure 4 -

 4 Figure 4-7 Thermo-Mechanical coupling in geological media The governing equations solved in the Finite Element scheme are the mass balance of the solid phase (Figure 4-43), the balance of energy (Figure 4-44) and of momentum (Figure 4-45) for the medium.

Figure 4

 4 Figure 4-8 CODE_BRIGHT model geometry considered Table 4-1 Parameters considered for the numerical model analysis Porosity () 0.2 Matrix linear thermal expansion -m (ºC -1 ) 6e-6 Cb 0.4 Bond linear thermal expansion -b (ºC -1 ) 6e-6 Km (GPa) 50 r1 0.01 Kb (GPa) 50 r0 (MPa) 1e-5

Figure 4 -

 4 Figure 4-9 and Figure 4-10 presents the results obtained with the CODE_BRIGHt model and the algebraic expression.It is important to mention that here in after the strains are draw with rock mechanics convention, in other words, positive values express expansion. Results obtained with both solutions are very similar, however, some small differences are observed at the end of the heating interval, possible related to the integration scheme of the finite element model.

Figure 4 - 9

 49 Figure 4-9 Strains evolution comparison between algebraic solution and CODE_BRIGHT model

Figure 4 -Figure 4 -

 44 Figure 4-11 Loading function evolution for an increment of temperature between 0ºC to 30ºCFigure4-12 presents the evolution of the loading function F during thermal cycles between 10ºC to 50ºC, without considering any rate dependency ( b =0). In this case, all the damage occurs during the first cycle, when the increase of bond mean stress brings the damage locus to the value of 6 MPa. Afterwards, temperature cycles cause changes of mean stress below 6 MPa, that is inside the elastic zone and no damage occurs.

Figure 4 -Figure 4 -

 44 Figure 4-12 Rate-independent loading function evolution for thermal cycling between 0ºC to 50ºC: a) initial condition; b) 1 st cycle, end of heating; c) 1 st cycle, end of cooling; d) 200 th cycle, end of heating; d) 200 th cycle, end of cooling

Figure 4 -

 4 14 shows theHashin & Shtrikman (1963) bounds envelopes for the considered ratios of Km/Kb. Values considered for Xm/Xb ratio is in all cases located within the admissible zone.

Figure 4 -

 4 Figure 4-14 Coefficients Xm and Xb considered for the evaluation of different Km -Kb ratio

Figure 4 -

 4 Figure 4-15 Evaluation of bulk moduli with Xm/Xb=0.5: total strains evolution

Figure 4 -

 4 Figure 4-16 Evaluation of bulk moduli for Xm/Xb=0.5: internal stresses evolution

Figure 4 -

 4 Figure 4-17 Evaluation of bulk moduli for Xm/Xb=0.5: matrix and bonding mechanical strains evolution

Figure 4 -

 4 Figure 4-18 Evaluation of bulk moduli for Xm/Xb=0.5: bulk modulus and damage multiplier evolution Xm/Xb=1 (Xm=0.4 Xb=0.4).This case addresses a composite material made of two components having a symmetric effect on the rock behavior. This is obtained by setting equal volumetric content (cb = cm), equal thermal expansion coefficient and equal values for coefficients Xm and Xb. In this case, because the definition of Xm and Xb (equations 4-25 and 4-26), the total strains for the bond and matrix are equal and must thus correspond only to the thermal strains (Figure4-19), In other words, each component can expand freely under the thermal load without violating nor the stress and strain partition, neither the bond and matrix constitutive laws.

Figure 4 -

 4 Figure 4-19 Evaluation of bulk moduli for Xm/Xb=1: total strains evolutionTherefore, no internal stresses are generated, as can be observed in Figure4-20, In other words, because of the specific internal microstructure considered for the composite material, matrix and bond can develop their complete thermal expansion without interacting each other. Consequently, no damage will be registered in the bond (damage multiplier equal to zero) and

Figure 4 -Figure 4 -

 44 Figure 4-20 Evaluation of bulk moduli for Xm/Xb=1: internal stresses evolution

Figure 4 -

 4 Figure 4-22 Evaluation of bulk moduli for Xm/Xb=2: total strains evolution

Figure 4 -

 4 Figure 4-23 Evaluation of bulk moduli for Xm/Xb=2: mechanical strains evolution

Figure 4 -

 4 Figure 4-24 Evaluation of bulk moduli with Xm/Xb=2 -Internal stress evolution

Figure 4 -

 4 Figure 4-25 Evaluation of bulk moduli for Xm/Xb=2: bulk modulus and damage multiplier evolution

Figure 4 -

 4 Figure 4-26 Evaluation of thermal expansion coefficient for Xm/Xb=0.5: total strains evolution

Figure 4 -

 4 Figure 4-27 Evaluation of thermal expansion coefficient for Xm/Xb=0.5: internal stress evolution

Figure 4 -

 4 Figure 4-28 Evaluation of thermal expansion coefficient for Xm/Xb=0.5: bulk modulus and damage multiplier evolution Xm/Xb=1 (Xm=0.4 Xb=0.4).

  5 and mb=2 in Figure 4-29).

Figure 4 -

 4 Figure 4-29 Evaluation of thermal expansion coefficient for Xm/Xb=1: total strains evolutionConversely to what was observed in the analysis of the bulk modulus, if Xm is equal to Xb some mechanical strains and internal stresses are going to be developed in the components, except for the case where m/b be equal to Xm/Xb, as presented in Figure4-30.

Figure 4 -

 4 Figure 4-30 Evaluation of thermal expansion coefficient for Xm/Xb=1: internal stress evolution

Figure 4 -

 4 Figure 4-31 Evaluation of thermal expansion coefficient for Xm/Xb=1: bulk modulus and damage multiplier evolution Xm/Xb=2 (Xm=0.8 Xb=0.4).

Figure 4 -

 4 Figure 4-32 Evaluation of thermal expansion coefficient for Xm/Xb=2: total strains evolution

Figure 4 -

 4 Figure 4-33 Evaluation of thermal expansion coefficient for Xm/Xb=2: internal stress evolution

Figure 4 -

 4 Figure 4-34 Evaluation of thermal expansion coefficient with Xm/Xb=2 -bulk modulus and damage multiplier evolution

Figure 4 -

 4 Figure 4-35 Coefficients Xm and Xb considered for the evaluation of different bond content Xm/Xb=0.5 (Xm=0.125 Xb=0.25).

Figure 4 -

 4 Figure 4-36 Evaluation of bond content for Xm/Xb=0.5: total strains evolution

Figure 4 -

 4 Figure 4-38 Evaluation of bond content for Xm/Xb=0.5: bulk modulus and damage multiplier evolution

Figure 4 -

 4 Figure 4-39 Evaluation of bond content for Xm/Xb=1: total strains evolution

Figure 4 -

 4 Figure 4-40 Evaluation of bond content for Xm/Xb=1: internal stress evolution

Figure 4 -

 4 Figure 4-41 Evaluation of bond content for Xm/Xb=1: bulk modulus and damage multiplier evolution Xm/Xb=2 (Xm=0. 5 Xb=0.25).

Figure 4 -

 4 Figure 4-42 Evaluation of bond content with Xm/Xb=2 -Total strains evolution

Figure 4 -

 4 Figure 4-43 Evaluation of bond content for Xm/Xb=2: internal stress evolution By comparing Figure 4-43 and Figure4-44 it can be observed that the maximal bond stress is obtained when the bond content is 0.6, as usual, this is related to a higher increase of damage multiplier, therefore to a higher reduction of bond bulk modulus. Nonetheless, the composite material bulk modulus is not much affected by this reduction. This can be associated to two different phenomena: first, the value of Xm considered is larger than Xb, which generates that K be principally controlled by Km, additionally K is less affected by Kb as higher is Cb. This last condition explains why even if the damage multiplier for Cb=0.2 increases the less, the reduction of the bulk modulus of the composite is the highest.

Figure 4 -

 4 Figure 4-44 Evaluation of bond content for Xm/Xb=2: bulk modulus and damage multiplier evolution

Figure 4 -

 4 Figure 4-45 Coefficients Xb and Xm considered for the numerical modeling of sample CH3

Figure 4 -

 4 46 and Figure4-47 present the results obtained with this first approximation for the block and cliff samples respectively. It can be observed that the results obtained with the GP modeling follow closely the evolution observed with the experimental results.

Figure 4 -

 4 Figure 4-46 Gauss point modeling for block samples

Figure 4 -

 4 Figure 4-48 Geometry and mesh for the axisymmetric model

Figure 4 -

 4 Figure 4-49 presents the comparison of the bulk modulus evolution obtained with the axisymmetric model at points A, B, with the Gauss Point model and the experimental data.

Figure 4 -

 4 Figure 4-49 Sample B3: Bulk modulus evolution

Figure 4 -

 4 Figure 4-50 Sample B3: Bulk modulus evolution in the sample

Figure 4 -

 4 53 and Figure 4-54 presents the comparison between measurements and results obtained with the asymmetrical model for axial and radial strains respectively. The numerical model results considered are obtained at point B, because experimental measurements were performed at the boundary of the sample.

Figure 4 -Figure 4 -

 44 Figure 4-53 Sample B3: Axial strains

Figure 4 -Figure 4 -

 44 Figure 4-55 Sample B3: Mean stress evolutionNevertheless, from Figure4-56 it can be observed that strains response obtained for a 24 hours period is very similar to the measurements, validating the thermal expansion coefficients considered.

Figure 4 -

 4 Figure 4-57 presents the bulk modulus evolution obtained with the axisymmetric model, Gauss Point model and measured data. Note that as was obtained for sample B3, the bulk modulus reduces less within the sample (point B) than at its face (point A), and the results obtained with the Gauss Point model are closer to the results obtained at point B, because the temperature evolution imposed for the gauss point model is the one registered at the boundary of the sample. Nevertheless, results obtained with the axisymmetric model follows closely the data obtained in the experimental program.

Figure 4 -Figure 4 -

 44 Figure 4-57 Sample B7: Bulk modulus evolution Figure 4-58 presents the evolution of the damage multiplier obtained at points A and B, note that the evolution is very similar to the one observed for sample B3 (Figure 4-52), the damage multiplier is much higher at the boundary of the sample (point B) than within it (point A). In fact, from Figure 4-59 it can be observed that the evolution of damage multiplier thought the thermal cycles for sample B7 is almost equal to the one obtained for sample B3 (Figure 4-52). Similar results are obtained in the evolution of the bulk modulus.

Figure 4 -

 4 Figure 4-58 Sample B7: Damage multiplier evolution

Figure 4 -

 4 Figure 4-59 Sample B7: Damage multiplier evolution in the sample

Figure 4 -

 4 Figure 4-61 Sample CV4: Bulk modulus evolutionMoreover, the bulk modulus tends to a stable value after about cycle 500, Figure4-64 shows the stiffness evolution in the sample. Note that the stiffness values from cycle 500 and cycle 1000 are very similar and almost constant for the whole sample, while for cycle 100 it is lower closer to the symmetry axis.

Figure 4 -

 4 Figure 4-62 Sample CV4: Bulk modulus evolution in the sample

Figure 4 -

 4 Figure 4-63 Sample CV4: Damage multiplier evolution

Figure 4 -

 4 65 and Figure4-66 presents the axial and radial strains obtained in the experimental program and the modeling of sample CV4. As observed in the former modeling a slight strains accumulation is obtained with the axisymmetric model, -1.10e-4 m/m in axial and -7.25e-5 m/m in radial direction after 1000 thermal cycles. Moreover, the strains amplitudes (Figure4-67) are adequately reproduce by the model.

Figure 4 -Figure 4 -

 44 Figure 4-64 Sample CV4: Damage multiplier evolution in the sample

Figure 4

 4 Figure 4-68 Sample CV4: Mean stress evolution

Figure 4 -

 4 Figure 4-69 Sample CV5: Bulk modulus evolutionFrom Figure4-69, it can be remarked that similar to sample CV4, the bulk modulus reduces faster within the sample (point A) than at the boundary of the sample (point B). Furthermore, the bulk modulus evolution obtained with the Gauss Point model is placed above the evolution obtained with the axisymmetric model at point A and B. It is also important to point out that this parameter show almost any evolution after 500 cycles.

Figure 4 -Figure 4

 44 Figure 4-70 Sample CV5: Damage multiplier evolution

Figure 4 -Figure 4 -

 44 Figure 4-72 Sample CH2: Bulk modulus evolution

Figure 4 -

 4 Figure 4-74 Sample CH2: Damage multiplier evolution

Figure 4 -

 4 76 and Figure4-77 presents the results obtained with the numerical modeling for radial and axial strains respectively, at the boundary of the sample, compares with the measurements results.

Figure 4 -Figure 4 -

 44 Figure 4-76 Sample CH2: Radial strains evolution

Figure 4 -

 4 Figure 4-80 Sample CH3: Bulk modulus evolution

Figure 4 -

 4 Figure 4-81 Sample CH3: Damage multiplier evolution

Figure 4 -

 4 Figure 4-82 Sample CH3: Strains evolution

  a parameter of viscosity. Note that the dependency of F db on λ db through r is already considered in ∂res F db ∂r db (r and λ db are considered independent variables in the derivation). tangent matrix of the damage model must be provided to construct the tangent matrix of the composite model that will enter in the global tangent stiffness matrix of the Finite Element scheme. It is obtained by stating the consistency condition on the damage loading function

Table 3 13

 3 Thermal expansion coefficient, from strains measurements -Samples group1 ...

Table 3 14

 3 Thermal expansion coefficient, from strains measurements -Samples group2 ...

Table 3 15

 3 Vp measurements without and with strain gauges ...............................................

Table 3 16

 3 Bulk modulus evolution ..........................................................................................

Table 3 17

 3 Uniaxial compressive strength results ....................................................................

Table 3

 3 18 calculated Young modulus .

  3.1. Introduction ............................................................................................................... 3.2. La Roque Gageac's limestone general characterization. .......................................... 3.3. Samples description................................................................................................... 3.4. Experimental layout ................................................................................................... 3.5. Mechanical response .................................................................................................

	3.6. Generation and propagation of fissures. ...................................................................
	3.7. Experimental work conclusions ...............................................................................

Table 3 -

 3 1 Density and elastic wave propagation velocity

	Sample	Block	Density (kg/m 3 )	Vp (m/s)	Vs1 (m/s)	Vs2 (m/s)
	998	A	2434	5420	3014	3016
	999	A	2399	5046	2814	2870
	1000	A	2444	5029	2875	2856
	1013	A	2317	4757	2736	2840
	1014	A	2437	3269	2675	2602
	1015	A	2340	4684	2731	2811
	1016	A	2400	4786	2537	2612
	1017	A	2290	4681	2561	2477
	1018	A	2423	5241	2675	2720
	1019	A	2367	4884	2572	2663
	1020	B	2373	4431	2427	2462
	1021	B	2366	4078	2376	2411
	1022	B	2377	4153	2571	2563
	Mean		2382	4651	2659	2683
	St deviation		45	545	171	178

Table 3 -

 3 2 Porosity for La Roque Gageac Limestone

	Sample	d (kg/m 3 )	s (kg/m 3 )	nc Connected porosity	nt Total porosity	Wsat (%)
	Block1	2356	2868	0.111	0.179 4.7
	Block2	2329	2835	0.113	0.179 4.8

Table 3

 3

			-3 Mechanical properties LaRG limestone
	Sample	Block	c (MPa)	tb (MPa)	c (MPa) 1 3	E (GPa)	 
	998	A			80.2	5.0	42.58	0.15
	999	A			81.68	10.0	46.97	0.12
		A			65.00	2.5	38.35	0.15
		A		6.18
		A		5.49
		A		5.10
		A		4.16
		A		5.44
		A	47.44
		A	47.28
		A	53.73
		A	59.71
		A	35.44
		A	31.30
		A	63.29
		B	51.68
		B	43.20
		B	40.65
		C		5.63
		C		4.91
		C		5.86
		C		5.50
		C		5.58
		C		5.51
	Mean	C	47.37	5.39			42.63	0.14
	St deviation		9.63	0.50			3.52	0.01

Table 3 -

 3

	4 RMR evaluation
	Property	Value	Note
	Strength	c=47.4 MPa	4
	RQD	75 -90%	17
	Discontinuities separation	200 -600 mm	10
	Discontinuities roughness	Roughness surface	25
	water	dry	15
	RMR		71

Table 3

 3

			-5 Thermal expansion results
	Prismatic samples P1 hor 1 P1 vert	P1 hor 2	P2 hor 1 P2 vert P2 hor 2
	10 -6 ºC -1 )	3.1	14.2	7.0	6.2	4.8	1.6
	
	Cylindrical samples	C1 vert		C1 hor	C2 vert		C2 hor
	10 -6 ºC -1 )	6.1		9.7	7.6		10.6

Table 3

 3

-6 presents the characteristics of the cylindrical samples and Table 3-7 of the prismatic ones. It is important to mention, that prismatic sample from group 2 (C4) is a sample with natural fissures, removed during a technical visited realized in April 2016.

Table 3 -

 3 6 Geometry and characteristics cylindrical samples (CV: samples drilled in vertical sense from the face of the cliff, CH: Samples drilled in horizontal sense from the face of the cliff)

	Sample	Group	Diameter (mm)	Height (mm)	Weight (gr)	VP (m/s)	Porosity (%)
	B1	Group 1-Block 1	51.03	101.07 473.00		15.86
	B2	Group 1-Block 1	50.79	100.96 500.00		10.13
	B3	Group 1-Block 1	50.85	101.13 493.00		11.75
	B4	Group 1-Block 1	50.93	100.99
	B5	Group 1-Block 1	51.02	100.96 474.52		15.47
	B6	Group 1-Block 1	50.40	100.86 458.29		16.29
	B7	Group 1-Block 1	50.61	101.31 457.44		17.49
	B8	Group 1-Block 1	51.03	100.92 475.08		15.38
	B9	Group 1-Block 1	50.05	101.02 457.96		15.29
	B10	Group 1-Block 2	49.71	112.04 511.48		13.52
	B11	Group 1-Block 2	50.03	112.07 511.02		14.73
	B12	Group 1-Block 2	50.21	112.16 514.57		14.80
	B13	Group 1-Block 2	49.84	112.03 512.04		13.88
	B14	Group 1-Block 2	49.43	112.07 501.35		14.29
	B15	Group 1-Block 2	49.22	112.07 498.09		14.13
	B16	Group 1-Block 2	49.31	112.06 504.65		13.29
	B17	Group 1-Block 2	49.83	112.02 511.27		13.95
	B18	Group 1-Block 2	49.64	112.07 510.63		13.45
	B19	Group 1-Block 2	49.26	111.98 495.98		14.56
	B20	Group 1 -Block1	50.66	96.60	472.72
	CV1	Group 2	50.77	111.31 505.71		18.40
	CV2	Group 2	50.73	111.17 490.70		20.87
	CV3	Group 2	50.73	111.23 503.15		18.91
	CV4	Group 2	50.68	111.01 498.78		20.76
	CV5	Group 2	50.75	111.33 497.87		19.91
	CH1	Group 2	50.89	111.14 509.87		18.29
	CH2	Group 2	50.85	111.02 5106.6		17.93
	CH3	Group 2	50.90	111.08 514.87		17.46
	CH4	Group 2	50.88	111.09 498.57		20.03
	CH5	Group 2	50.82	111.11 493.33		20.70

Table 3

 3

		-7 Geometry and characteristics prismatic samples
	Sample Group	Width (mm)	Length (mm)	Height (mm)	Weight (gr)
	C1	Group 1	81.50	61.02	83.56 983.04
	C2	Group 1	81.44	61.72	81.11 920.08
	C3	Group 1	82.09	61.62	80.87 973.69
	C4	Group 2 82.80* 72.90* 74.80* 953.46

*Approximated measurements because the irregular shape of the sample (Natural sample)

Table 3

 3

			-8 Mineralogical analysis samples
	Sample	Group	Thin sections Spectrometry Rx

Table 3 -

 3 9 Number of cycles imposed and test performed for each sample

	Sample	Cycles applied	Def	Vp-Vs	UCS	DIC
	B1	948	X
	B2	948	X	X
	B3	948	X	X
	B4	Measure of temperature inside the sample
	B5	828

table 3

 3

		-10 Accumulated axial strains
	Samples obtained from blocks (group 1) Samples obtained from the cliff face (group 2)
	Sample	Axial accumulated strain (mm/mm)	Sample	Axial accumulated strain (mm/mm)
	B1	2.7e-4	CV1	-2.6e-4
	B2	-8.1e-4	CV4	2.5e-4
	B3	9.0e-4	CH2	-26.0e-4
			CH5	-8.6e-4

Table 3

 3

	Sample	Axial accumulated strain (mm/mm)	Sample	Axial accumulated strain (mm/mm)
	B1	-0.9e-4	CV1	0.26e-4
	B2	9.0e-4	CV4	2.5e-4
	B3	2.7e-4	CH2	-3.3e-4

-11 Accumulated radial strains Samples obtained from blocks (group 1) Samples obtained from the cliff face (group 2)

Table 3 -

 3 12 Mean deformations amplitude for thermal-strains evolution

	Sample	Vertical (m/m)	Radial (m/m)
	B1	2.87e-5 17.5e-5
	B2	14.1e-5 17.3e-5
	B3	6.84e-5 10.3e-5
	CV1	46.9e-5 7.14e.5
	CV4	12.3e-5 1.77e-5
	CH2	28.6e-5 2.86e-5
	CH5	38.4e-5

Table 3 -

 3 13 Thermal expansion coefficient, from strains measurements -Samples group1

		Sample B1	Sample B2	Sample B3
	Cycle	Vertical (ºC -1 )	Radial (ºC -1 )	Vertical (ºC -1 )	Radial (ºC -1 )	Vertical (ºC -1 )

Table 3 -

 3 14 Thermal expansion coefficient, from strains measurements -Samples group2

		Sample CV1	Sample CV4	Sample CH2	Sample CH5
	Cycle	Vertical (ºC -1 )	Radial (ºC -1 )	Vertical (ºC -1 )	Radial (ºC -1 )	Vertical (ºC -1 )	Radial (ºC -1 )

Table 3 -

 3 15 Vp measurements without and with strain gauges

	Sample	Vp-without gages (m/s)	Vp-with gages (m/s)
	B2	4651	4637
	B3	4470	4480
	CV1	3247	3227
	CV4	2842	2810
	CH2	3037	3042
	CH5	2830	2827

Table 3

 3

		-16 Bulk modulus evolution
	Sample	Initial K (MPa)	Final K (MPa)	Variation (%)
	B2	26500	26024	1.7
	B3	24300	23500	3.0
	B5	20600	19400	6.0
	B6	20300	19460	4.5
	B7	19330	18025	6.7
	CV1	13638	1068	23.2
	CV2	9630	6700	30.4
	CV3	7160	5694	20.5
	CV4	10470	9020	13.8
	CV5	8590	5630	34.4
	CH1	10610	7814	26.4
	CH2	10360	8145	21.3
	CH3	14030	9950	29.1
	CH4	9389	7840	16.5
	CH5	8539	7870	7.80
	3.5.3. Uniaxial compressive strength.

Uniaxial compressive strength tests are performed on 15 samples that have been submitted to a different number of thermal cycles (see Table

3-9

). The initial and final elastic compression waves propagation velocities are considered as comparison parameters. Table

3

-17 presents the results obtained for each sample considered, related to its initial and final Vp value.

It can be observed, that samples with a higher Vp value (final Vp) present a superior strength (UCS) value, even if the sample has been submitted to thermal cycles. For instance, sample B13 has an initial Vp value of 4485 m/s, after 420 cycles of temperature the elastic wave propagation velocity decrease until 4450 m/s. Nevertheless, its UCS value (47.77 MPa) is similar to the value obtained in sample B18 (47.54 MPa), that has a similar initial Vp value, but has suffered no thermal cycles.

Table 3

 3

		-17 Uniaxial compressive strength results
	Sample	Cycles imposed	Initial VP (m/s)	final VP (m/s)	UCS (MPa)
	B5	828	4342	4175	41.23
	B6	828	4282	4151	38.03
	B7	828	4177	3981	37.56
	B8	624	4407	4181	39.08
	B9	624	4404	4181	39.93
	B10	624	4321	4130	40.01
	B11	420	4373	4198	41.96
	B12	420	4402	4154	41.72

Table 3

 3

			-18 calculated Young modulus
	Sample	Etan (GPa)	Esec (GPa)	E50 (GPa)	E80 (GPa)
	B5	13.86	12.03	17.58	17.50
	B6	13.37	11.74	15.30	15.08
	B7	13.24	11.68	15.29	14.84
	B8	13.79	11.74	15.70	15.78
	B9	13.40	12.02	15.81	15.76
	B10	14.50	12.76	16.94	16.73
	B11	15.10	12.09	15.03	15.59
	B12	15.13	13.11	16.90	17.33
	B13	16.21	13.86	17.95	18.47
	B14	15.69	13.60	17.39	15.59
	B15	15.59	13.67	17.38	17.57
	B16	16.00	14.09	17.73	18.20
	B17	14.30	12.84	16.66	16.30
	B18	16.22	14.08	17.78	18.30
	B19	15.41	13.83	17.59	17.60

Table 4 -

 4 3 presents the comparison between experimental values of rock modulus and the values obtained by the model after back-analyzing Xm and Xb. Figure 4-45 presents an example of the check of the Hashin-Shtrikman condition for Xb and Xm, realized with the values obtained for sample CH3.

Table 4

 4

		-3 Initial bulk modulus values (MPa)
		Experimental results Numerical modeling
	B3	24300	24020
	B7	19330	19360
	CV4	10470	10750
	CV5	8590	8300
	CH2	10610	10200
	CH3	10360	13570

  5.1. Concluding remarks. ................................................................................................ 181 5.2. Future development. ............................................................................................... 183

	5.1. Concluding remarks.

Numerical modeling

Acknowledgements

Acknowledgements

Additionally, [START_REF] Carol | On the formulation of anisotropic elastic degradation. I. Theory based on a pseudo-logarithmic damage tensor rate[END_REF] introduce a loading function (F) to define an elastic domain in the stress space F( b , b ) in which the stiffness remains undamaged. Once the loading surface is reached (F=0), the onset of degradation takes place. In the case of the isotropic damage, the loading function F is defined as:

Eq. [4][5][6][7][8][9][10][11][12][13][14] where u is the bond current (secant) elastic energy and r(L) is the damage locus. r evolves generally with the intensity of damage L. Three different laws have been implemented for this evolution:

Linear: r = r 0 + r 1 𝐿

Eq. [4][5][6][7][8][9][10][11][12][13][14][15] Exponential: r = r 0 e 𝑟 1 𝐿 Eq. [4][5][6][7][8][9][10][11][12][13][14][15][16] Logarithmical: r = r 0 + r 1 ln ( 1 + 𝐿 1 + L )

Eq. [4][5][6][7][8][9][10][11][12][13][14][15][16][17] Where, r0 and r1 are model parameters and indicate the initial damage locus and the damage rate, respectively.

A rate dependency is introduced following the theory of Perzyna (1966) proposed for elastoplastic models. In this theory, the plastic multiplier is expressed as a function of the distance between the current stress point and the inviscid plastic locus. This concept can be extended to damage theory by considering that damage multiplier L is a time-dependent function of the distance between the damage locus and the current stress state. The expression for this function is given by:

where dt is the time increment,  b the bonding viscosity and  are the Macauley brackets.

Matrix and bond coupling.

The constitutive law of the composite material is obtained by coupling the ones of the bonds and the matrix. This coupling should insure that the external strain applied on the composite material must be compatible with the strains of the bonds and the matrix, while external stresses must be in equilibrium with the local ones.

Compatibility between external and local strains is obtained from the volume conservation; same results can be obtained from equations of mass conservation of bonds, matrix and medium. The volume conservation states that the total volume V is the sum of volume of matrix (Vm), voids (Vv) and bonds (Vb):

Eq. 4-19

Experimental results modeling.

The modified composite material constitutive law presented before and implemented in CODE_BRIGHT is used to reproduce the response of cyclic thermal tests performed in the laboratory and presented in section 3. Six tests have been modelled: two of them were performed on samples extracted from fall blocks (B3 and B7) and four on samples obtained from boreholes drilled from the cliff face (CV4, CV5, CH2 and CH3). CV4 and CV5 have been cored in the vertical direction, CH2 and CH3 in the horizontal one.

The parameters considered for the numerical modeling of these samples are listed on Table 4-2, it can be observed that rate dependency considered for all the samples is the same (1000).

Volumetric fraction of voids (porosity), matrix (cm) and bonds (cb) vary according to the sample. Matrix is associated to the hardest mineral phase (quartz) while bonds are associated to the calcite phase, more susceptible to experience degradation. In block samples, quartz has a volumetric fraction close to 10%. In cored samples, quartz and calcite appeared to have a similar

Appendix A. Mercury intrusion porosimetry.

In order to evaluate if atmospheric thermal cycles may have any identifiable effect on the porosimetry of La Roque Gageac limestone, the porosity of four small samples submitted to 0, 276, 420 and 630 cycles between 10ºC to 50ºC were analyzed with the mercury intrusion porosimetry technique.

Figure A-1 presents the pore size distribution obtained for the 4 samples, it can be observed that all samples present similar response, therefore there is any remarkable response on the porosity when samples are subjected to thermal cycles. For these samples, a predominant pore size of 1.56 m is observed. It can be observed that any new fissure is observed at least at photo scale, similar to what was observed on sample C2 and presented in section 3.6.2. It is important to point out that the results obtained with the three samples are comparable, in terms of strains magnitude and distribution, which corroborates that the small strains observed can be attributed to variations on the illumination pattern. Moreover, with the three samples evaluated it is observed that for this strain magnitude the face with the speckle pattern presents less noise.

C.1 Rate equations for general bond damage model.

Implementation is based on the implicit integration scheme proposed by Simo and Hughes, (1998) for elastoplastic models. The scheme is generally used to integrate all model equations: damage rate equations (bond model), elastoplastic rate equations (matrix model) and bondmatrix coupling equations (stress and strain partition and Xb dependency). In order to use this scheme, damage model is described by the general rate equations proposed by Carol et al. (1996) who showed that they are similar to those of elastoplasticity. They include the elastic law, the flow rule, the hardening law and the Kuhn-Tucker condition. where dλ db is equal to the parameter dL defined by [START_REF] Carol | On the formulation of anisotropic elastic degradation. I. Theory based on a pseudo-logarithmic damage tensor rate[END_REF].

The strain entering in the mechanical model is defined as the difference between the total (ε ij totb ) and the thermal strains (ε ij tempb ).

Bond strain ε ij b andnbond stress σ ij b are related by the secant stiffness and compliance matrix:

The secant matrix D ijkl eb is related to the initial (undamaged) secant matrix D ijkl eb0 by the damage law:

Finally, the damage loading function is defined in terms of energy:

and r(λ db ) is the bond damage locus Within this formalism, the flow rule is directly the strain tensor:

The hardening

The evolution of the damage locus can be defined by three different expressions:

Linear: r = r 0 + r 1 λ db in which case, h db = r 1

Exponential: r = r 0 e r 1 λ db in which case, h db = r 1 r

Logarithmical r = r 0 + r 1 ln (

C.2 Discrete equations for the damage model.

Elastic law.

The scheme used to solve the discrete equation is not directly based on the linearization of the rate equations. The presence of an exponential in the damage law makes indeed inaccurate its approximation by a linear function, leading to small steps and a high number of sub-increments in the integration procedure. Alternatively, the discrete stress increment has been directly expressed as the difference between its value at the beginning and the end of the step:

))

The previous equations implicitly state the following discrete flow rule

It can be noted that the classical discrete linearization of the damage law is recovered by developing in series the term within parenthesis up to order 1. However, as previously indicated, this linearization has not been considered for accuracy issues.

Noting:

the bond trial stress increment and the stress increment tensor, the system of discrete equations becomes:

The discrete expression of the hardening law simply reads:

Kuhn-Tucker condition.

The Kuhn-Tucker condition is given by

and

Stress point algorithm for damage model.

The stress point algorithm is based on an implicit scheme. The Newton-Raphson procedure is used to solve the non-linearity of the system.

Newton-Raphson procedure

The NR procedure aims at cancelling the residues of the elastic law (res ij σ b ), hardening rule (res r db ) and Kuhn-Tucker condition (res F db ). There residues read: