Simulation of an experimental fire in an underground limestone quarry for the study of Paleolithic fires

Abstract Numerous fire marks occur on the walls of the Chauvet-Pont d’Arc cave. Dating indicated that some of the fires were contemporary to the Aurignacian. Violent thermal shocks were observed in surprisingly narrow areas of the cave. This raises numerous archaeological questions about the function of the fires; the answers depend on the location of the hearths, and the intensity of the fires. Numerical simulation was used here to provide information about the behaviour of fires in such confined spaces. An underground non-archaeological site, in a limestone quarry, was equipped to monitor fires in an environment similar to that of the Megaceros gallery of the Chauvet-Pont d’Arc cave. The fire and the movement of heat and smoke in the quarry were simulated by the open source code “Fire Dynamics Simulator (FDS)”. Results were validated on wall temperatures recorded behind and above the fire. The thermo-mechanical impact of the fire on the rock was simulated with CAST3M software, providing the most probable zones for limestone spalling due to thermal gradients. The validated approach will, in a forthcoming study, be applied to the Chauvet-Pont d’Arc cave, in which coupled simulations in the air and in the rock should indicate the location of the hearths and the intensity of the fires that generated the marks.


Background
In the Chauvet-Pont d'Arc cave, from the earliest investigations, marks of rubification and shedding related to violent thermal shocks associated to combustion areas were noted on the walls and ceiling of the Megaceros gallery and the gallery of the Crosshatches [1] [2]. Research begun since 2008 on the taphonomy of the decorated walls allowed the identification of thermal alterations in the entrance sector, in the Chamber of the Bear Hollows and in several places of the Megaceros Gallery [3] [4]. They are revealed by a change of colour of the limestone and spalling of the rock surface [5]. The thermoluminescence study and the comparison with a thermal referential indicate a range of temperatures between 300 and 375 C [6].
The dating of charcoal concentrations located close to walls or ceilings heated in the Entrance sectors and in the Gallery of Crosshatches [7] as well as the presence in the Chamber of the Bears and in the Megaceros Gallery of drawings related to the ancient phase of human frequentation [8] [9] and posterior to the thermal marks, indicates that some fires might have been contemporary to the Aurignacian. Recent thermoluminescence dating results obtained on heated calcite fragments sampled at the wall of the Megaceros gallery and at the ceiling of the Paleolithic entrance confirmed that fires were lit during the first paleolithic occupations of the cave [10]. Nevertheless, it cannot be totally excluded that others were the works of the Gravettians.
The function of the fires has not yet been clearly identified but several hypotheses can be considered: lighting, beacons, torch rekindling, pigment production, protection from animals, heat or light production without practical need other than the one linked to the symbolism [11]. In the Chauvet-Pont d'Arc cave, few hearths have been found. They have been moved by bears or men, deliberately or not. Their absence leads to a lack of information, especially concerning their function. Locations of the fire marks in very narrow areas of the cave raises questions about lighting fires without being injured. Fires release large amounts of toxic gases into the atmosphere, and the evacuation of such gases constitutes a problem in confined areas. Simulations of a fire in a geometry similar to the Chauvet-Pont d'Arc cave may give additional information about the production of toxic gases, the temperatures close to the hearth and the temperature of the rock.
Various works on simulations of fires in underground or confined areas have provided valuable information on air flows and smoke evacuation. Computational Fluid Dynamics has been used to simulate spontaneous combustion of coal and evaluate oxygen levels in mines [12]. Catastrophic tunnel fires have placed a focus on fire safety issues concerning road and rail links, generating numerous simulation works [13] [14]. Fires occurring in subways have also been thoroughly studied from a numerical point of view [15]. The thermomechanical impact of fires on the structure of a car park has equally been studied [16].
Works on toxic levels of gases emitted during combustion have given valuable information about the toxicity of the fumes and the thermal hazards. Purser [17] proposed the definition of a dose (Fractional Effective Dose) to predict when an incapacitating or lethal quantity of fumes has been received. Following this research, the toxicology of the combustion gases in confined areas has been applied to aircraft [18] and more recently to tunnels [19].
Simulations of fires in an archaeological context have so far focused on the ground supporting the fire. 3D simulations of the impact of outdoor fires on the ground have given information about fire duration or for instance reuse of archaeological hearths (Pincevent, Seine-et-Marne, France) [20]. More recently, heat transfer simulations on the ground of an experimental fire gave information on the duration of the fire starting from the temperature of the sediment [21] [22].
As far as the authors know, no works have simulated fires in the archaeological context of a rock art cave.

General objectives
The simulation of the impact of a fire in the air and on the rock of an underground area is a part of a project dedicated to studying fire in a monitored underground area. We aim to reach a better understanding of the phenomena involved in the combustion in such confined areas.
The project originates from archaeological observations of the walls of the Chauvet-Pont d'Arc cave. From the very first scientific studies of these caves, areas on the walls affected by a rubification due to heat from fires were identified in various zones. Arguments based on relative chronology and indirect dating suggested Aurignacian frequentation was responsible for the state of the walls. Thermoluminescence analysis of rubified rock (limestone) samples indicated that temperatures of at least 300 C were reached [5].
Starting from these archaeological observations, several questions arose. What types of hearths can produce marks (rubification, flakes, soot deposit) similar to those observed on the Chauvet-Pont d'Arc walls? How hot were their flames? What were the functions of the fires? Did the hearths occupy any particular positions? To answer these questions, several research objectives had to be defined. First of all, the marks observed in Chauvet-Pont d'Arc had to be experimentally reproduced by calibrating the fire correctly. The marks result from the acquisition of a pink to grey colour, flaking of the limestone, a temperature at the surface of the wall of at least 250 C, or soot deposit on the walls. Secondly, our aim was to characterize the effect of heat on the walls by a transformation gradient approach at macroscopic and microscopic scales. To achieve this, thermal and hygrometric measurements on the surface and in the rock were made. Also the different textures were studied as a function of the thermal history and the depth from the surface. Thirdly, the residues resulting from combustion for each type of experimental hearth were studied. Residues include charcoal, ash, or solid particles ejected into the atmosphere of the cave. Finally, the consequences of fires on the environment of the cavity were estimated by simulation taking a narrow gallery with dimensions similar to those of the Megaceros Gallery of the Chauvet-Pont d'Arc cave. Heat balance (distribution and variation of the temperatures), smoke flow, and air renewal were calculated and compared to experimental data. Once the model had been validated, it was used to simulate the impact of fires in archaeological sites such as the Chauvet-Pont d'Arc cave.
We aimed to be able to provide an accurate description of a fire in an archaeological context underground. For obvious reasons, the site chosen for the experiments was free of archaeological interest. We used an underground limestone quarry presenting morphological similarities with the burned galleries in Chauvet-Pont d'Arc cave.
The archaeological objective of the project was to obtain indications on the position of the paleolithic population in the caves near the fires. From the position of the fire, simulation should indicate the zones in which it would be impossible to remain without serious intoxication from the combustion gases.
The experimental process is detailed first. A description is given of the experimental site, the instrumentation used, and the way the fire burned. Then, fire simulations are presented in the site. An open source code (Fire Dynamics Simulation e FDS) was used to calculate the temperatures in the air and on the wall surface, as well as the level of toxic gases. Finally, we discuss the thermo-mechanical impact of the fire on the wall.

The experimental site
The experimental fires were set up in an abandoned underground limestone quarry situated in Lugasson in the Gironde district, France. The quarry presents several advantages: it is located far from any site of archaeological interest so fires would not cause any damage, its shape and volume are relatively easy to simulate (to mesh) and it is located close to the project members' laboratories.
The fires were made in an L shaped chamber with only one of the extremities open to the outside. They were placed against a wall located at the back of the chamber. The opening of the room was tall enough to enable us to maintain the fire. A photogrammetry technique was used to visualize the inner volume of the experimental room with a resolution of 1 cm. A view of the opening as seen from outside and a 3D representation of the volume are presented in Fig. 1.

Metrology
Numerous instruments were set up in the experimental cave during the fire tests. In this article we only present the temperatures of the chamber atmosphere, of the surface of the wall directly exposed to the fire and of the bulk of the limestone at various depths. Measurements of air velocity in the upper part of the entrance of the cave are also presented. These data provided information about the thermo-aeraulic behaviour of the cave during a fire. The simulation results in x3.5 are compared to the measured data for validation of the numerical model.

Instrumentation in the volume of the room
The objective of the instrumentation in the volume of the chamber was to measure the way the air temperature varied during the fire. In particular, measurements were made close to the cave's entrance. For this purpose, type K thermocouples were placed along two vertical poles placed at the entrance to the cave. On each pole, 3 thermocouples were distributed vertically in order to assess the vertical temperature gradient during the fire.
In addition, a hot-wire anemometer was placed on one of the poles, 1.67 m from the ground, in order to measure the air velocity during the fire. Being placed in the upper part of the cave's entrance, velocities measured mainly correspond to the departure of hot smoke. Fig. 2 shows the distribution of these different sensors in the cave and their position in relation to the fire.

2.2.2.
Instrumentation of the wall directly exposed to fire The objective here was to measure the thermal impact at the surface of the wall directly exposed to the fire and in the depth in the first centimetres of the limestone.  Six type-K thermocouples were fixed on the surface of the wall exposed to the fire. Four of them were vertically distributed 50, 90, 120 and 170 cm from the ground (Fig. 3).
Note that these 4 thermocouples were not perfectly aligned vertically. Significant temperature variations can then be observed according to the horizontal position of the sensors (a horizontal offset of few centimetres can result in temperature differences of the order of 100 C during a fire due to flame turbulence, see x3.5.4). A thermocouple (T ce in Fig. 3) was fixed to the ceiling 30 cm from the vertical wall. The temperature at the ceiling is an important parameter to study because it indicates whether the intensity of the experimental fires is sufficient to cause rubification of the ceiling limestone (250 C), as observed in some parts of ceilings in the Chauvet-Pont d'Arc caves. Finally, a thermocouple was placed on the ground 11 cm from the vertical wall, under the clay sole that supports the hearth (T fl in Fig. 3).
Other type-K thermocouples were introduced into the limestone of the wall. They were placed in holes drilled parallel to the surface exposed to the fire (Fig. 3, photograph). The sensitive end of the thermocouple was placed as near as possible to the vertical axis of the fire. Thermocouples were placed at 5 heights off the ground: 20, 60, 90, 120 and 160 cm. For each height 5 thermocouples were placed between 3 and 15 cm from the surface exposed to the fire.

Fire characteristics
The fire was made in the experimental cave on May 15th, 2012. The hearth was placed on a 6-cm thick clay sole (Fig. 4).
The wood used was Scots pine (Pinus sylvestris) as it was probably the same as that used as fuel in the Chauvet-Pont d'Arc Cave [23]. Four-kg wood bundles of the deadwood picked up from the ground in the Gironde area were prepared in advance. Each bundle consisted of pine branches of different diameters (between 15 and 50 mm) and length 40 cm.
The bundles of wood were placed on the hearth according to two successive protocols: -The first protocol, called "slow", lasted 1 h 30min, with a feeding rate of 3 bundles per hour. Each bundle was placed like a "tepee" around 10 cm from the wall. -The second protocol, called "fast", lasting 3 h and 45 min, with a feeding rate of 6 bundles per hour. Each bundle was laid on the residual embers, directly against the wall.
Finally, the fire was fed for more than 5 h, with a total mass of burned wood of around 100 kg. After introduction of the last bundle, measurements were continued until extinction of the hearth and natural cooling of the cave. Fig. 3. Left, cross-sectional view of the wall exposed to fire (see axis A-A in Fig. 2). Vertical distribution of the surface thermocouples and deep thermocouples (into the limestone). Right, picture of the wall before the fire.

Temperature measurements
Fig . 5 shows the temperatures measured at the surface of the wall directly exposed to the fire. Temperature oscillations correspond to the introduction of the wood bundles. The two protocols are clearly distinguished in the graphs: "slow" regime from 0 to 2 h, "fast" regime from 2 to 6 h.
A vertical thermal gradient is clearly observed, with surface temperatures decreasing from the ground up to the ceiling. At 50 cm from the ground (in almost direct contact with flames), the surface temperature reached around 600 C; at 170 cm from the ground, the temperature did not exceed 175 C. Fig. 6 shows the temperatures measured at the ceiling and under the clay sole.
The temperature at the ceiling remained lower than 200 C, i.e. lower than the temperature necessary to cause rubification of the limestone. Post-fire observations in the experimental cave confirmed the absence of limestone rubification of the ceiling. The fact that rubification of 3 m high ceilings occurred in certain areas of the Chauvet-Pont d'Arc cave shows that the experimental fire in the present study was probably smaller than the Paleolithic fires. The temperature measured under the sole remained lower than 100 C. This is explained by the insulating role played by the clay sole. Indeed the sole was initially quite wet, then the heat of the fire was attenuated by evaporation of the water in the clay. Fig. 7 shows the air temperatures measured at the experimental cave's entrance.
A vertical thermal gradient is clearly observed and up to 1.3 m, temperatures did not exceed 35 C. Temperatures in the upper section (1.47e1.67 m from the ground) reached 60 C. This can be explained by the circulation of the hot smoke in the upper section of the entry (see x 3.5).

Simulation of the fire
In order to better analyse the thermo-aeraulic behaviour of the gases into the quarry during the fire, a simulation has been carried out with the software Fire Dynamics Simulator (FDS). Moreover, the temperature simulated on the walls will be used as thermal  . Surface temperature on the wall directly exposed to fire versus time. In the graph caption, the number after the letter "T" is the thermocouple height from the cave floor.
boundaries for the thermo-mechanical simulation (x 4).

A brief description of FDS (Fire Dynamics Simulator)
The Fire Dynamics Simulator (FDS) used was developed by the National Institute of Standards and Technology (NIST), USA. The code was precisely described by McGrattan et al. [24]. It resolves the simulation of fire-generated heat and smoke transport in a confined area. Rehm and Baum [25] developed the equations to describe fire-driven fluid flows and proposed an approximate form of the Navier-Stokes equations appropriate for low-Mach number applications. The flow is considered as non-compressible (hydrodynamic) and non-acoustic (a low-Mach number flow).
The combustion model is based on the mixing-limited, infinitely, fast reaction of "lumped species". The chemistry of the fire is simplified as the reaction of fuel with air, giving products. Whereas the fuel is a single species, the air and products are referred to as "lumped species". "Lumped species" represents the mixture of gas species transported together and that react together. From the point of view of the numerical model it can be dealt with as a single species.
Toxicological aspects of a fire are also included in the code. The Fractional Effective Dose index [17] is commonly used to measure human incapacitation due to the combustion gases. It corresponds to the sum of the integration over time of different gas concentrations, here O 2 , CO and CO 2 concentrations. The values of the index show the zones where it is still possible to breathe normally, the zones where it becomes dangerous, and where the risk of death is high.

Description of the geometry
The shape of the quarry was reproduced by 24 obstacles defining the morphology of the site (Fig. 8). The elementary volume  considered was a cube with 10 cm edges. The whole simulation domain had 162000 nodes.
Concerning the boundary conditions, the heat diffusion into the walls are taken into account. For that, the thermal and physical properties of the limestone (see 3.3) are considered in the simulation. An open condition was attributed to the entrance.

Thermal characteristics of the limestone
The apparent density of the limestone was measured by hydrostatic weighing of samples of a few cubic centimetres taken from a wall of the experimental cave. This technique also gave an assessment of the water porosity of the material. Thermal conductivity and specific heat were measured on cylinders (diameter 45 mm, height 60 mm) also taken from a wall in the experimental cave. The method known as "Transient Plane Source" [26] was used. The average characteristics obtained are presented in Table 1.
The limestone of the experimental cave had a high porosity and low density. This results in a relatively low thermal conductivity.
For the fire-related simulations, the thermal properties of the limestone were considered constant regardless of the temperature. For the thermo-mechanical study (see x4), the properties of the limestone were temperature-dependent.

Fire characteristics
The location of the fire is described in Fig. 2. The fire was considered here as a constant gas burner. In that case, the fire is basically modelled as the ejection of gaseous fuel from a solid surface, with a specified Heat of Release Rate (HRR) in units of kW. The HRR is the product of the fire surface and the Heat Release Rate Per Unit Area, HRRPUA, in units of kW/m 2 . In this study, we assume that the fire surface keeps constant in time. This is not strictly correct as the morphology of the heart continuously changes during the fire but in this way, the HRRPUA becomes the only input parameter of the simulations. The fire is modelled as a parallelepiped with dimensions 0.68 Â 0.3 Â 0.1 m, and we consider that only the top and the sides of the prism are burning (i.e. 0.40 m 2 of heat source).
In order to reproduce the experimental conditions in which we fuelled the fire every 10 min with a new bundle of wood, the HRRPUA was not considered constant over the fire time. During the first minute of the simulation, the HRRPUA reached 100% of its maximum value, then during the next 8 min the value remained constant, and in the last minute it decreased to 80% of the initial value. 30 successive ramps occurred for 5 h, and then there was a return to HRRPUA ¼ 0 for an hour (see Fig. 9).  Since the HRR was not measured during the experiment (it could have been estimated by continuously weighing the burning wood [27]), it has been estimated by successive simulations with different maximum values for HRRPUA. For instance, Fig. 10 shows the comparison between measured and calculated temperatures on the wall close to the fire, depending on three HRRPUA maximum values (187, 250 and 375 kW/m 2 ). It is well known that the input value of the HRRPUA has a strong influence on the simulations results [28] and it is also demonstrated here. Indeed, the temperatures on the wall close to the fire are almost proportional to the HRRPUA maximum value.
The best estimation of the measured temperatures was obtained with a maximum value of HRRPUA of 250 kW/m 2 (see x 3.5.4). With a fire surface of 0.40 m 2 , it corresponds to a maximum HRR of 100 kW. By considering an effective heat of combustion for pine (with 30% of water content) equal to 12 MJ/kg [29], and a burning rate of 8 g/s, a power of 96 kW is obtained, which is close to the value used in the simulations.

Results and discussions
Simulation results are presented and compared to values recorded during the experimental fires. Data from the simulations are provided for the environment near the fire. In fact, simulation data close to the fire were validated by experiment. This is due to the fact that the simulations do not take variations of the outside weather into account although they are known to have a major impact on the exchanges.

Smoke distribution
Smoke was concentrated at the ceiling of the quarry (Fig. 11). The products of combustion (gases, soot, smoke) are released above the fire. As they are hotter than the ambient air, thus lighter, they remain in the upper parts of the quarry chamber. This is observed both numerically and experimentally.

Temperature distribution
Temperature had a cone-shaped distribution (Fig. 12).
The temperature gradient was very high, both horizontally and vertically. Temperatures in the lower part of the cone reached 300 C, i.e. where rubification occurred in our experiments (Fig. 12). Simulation provided similar results.

Velocity distribution
The velocity vectors were coloured depending on the air temperature (Fig. 13).
On being heated by the fire, the air becomes lighter and moves along the ceiling of the quarry. When it is becomes cooler on making contact with the ambient air, it falls again as it becomes heavier. The thermal gradient is very high both horizontally and vertically. Above the fire, on a horizontal slice, temperature shows a huge spatial variability. On both sides of the fire, on vertical slices the thermal gradient is high between the middle and the top of the room.

Comparison between numerical and experimental values of temperature
During the experiments, temperatures were recorded at several heights, 0.36 m, 0.88 m, 1.14 m, and 1.56 m above the hearth on the wall behind the fire (Fig. 3). Numerical and experimental values of temperature were compared at the locations of the thermocouples in Fig. 14.
On the whole, numerical simulations reproduce the global thermal behaviour with similar mean values. However, simulation curves are smoother than experimental ones. During the experiments, we observed the oscillations of the fire, visible on the values of the temperatures reported in Fig. 14. They were much greater in the experiments than in the simulations. This is due to the boundary conditions. Experimentally, they changed during the 5 h, a weather station placed outside the quarry showed the variation of meteorological conditions. Numerically, we chose an open boundary condition, the flows were balanced with the hot air leaving the quarry acting as the driving force for the cooler air coming in. Fire oscillation was thus weaker numerically than experimentally. At 36 and 88 cm above the hearth, the numerical temperature was higher as the heat was less diffuse. Then, at 1.14 m, simulation   underestimated the temperatures. The gradient is high along a horizontal slice through the fire, where temperatures can differ strongly within a couple of centimetres. The positions of the thermocouples in the fire were properly referenced, but it would be safe to consider an incertitude of 3 cm between their location and that in the simulator.
This incertitude is sufficient to explain the difference between the numerical and experimental values at 1.14 m.

Fractional effective dose index around the fire
For this value we have no experimental data, the simulation here is helpful to provide information on the toxicity of the air in the environment of the fire and on the possibility to breathe in the quarry during combustion (Fig. 15).
In FDS, the FED is based on the concentrations of O 2 , CO and CO 2 . The highest values of the FED correspond to the location of the smoke, which is coherent; the highest concentration of toxic gases is generally within the smoke. Concerning the values of the FED, when equal to 1, it means 50% of people will be able to either escape or survive post-exposure [18]. The safest places are on the ground while it is extremely dangerous to stand up near the fire. From an archaeological point of view, this index is useful to isolate zones from where it was possible to feed the fire in a more complex environment, as shown by Lacanette et al. [30] for the Chauvet-Pont d'Arc cave.

Objectives
The objective of the thermomechanical model was to estimate the likelihood of a Paleolithic fire (or several) on the basis of thermal (rubification of limestone) and mechanical marks (superficial spalling of limestone) in the Chauvet-Pont d'Arc cave. For this  purpose a thermo-mechanical model was coupled with the Computational Fluid Dynamics code previously described in x3.1.
Here, we present a validation of the approach. From the numerical thermal data obtained in x3.5.4 (the temperature field at the wall's surface), the diffusion of heat into the limestone is simulated and the rubification of the limestone analyzed. In a second step the mechanical behaviour of the wall under the effect of thermal diffusion is studied and the risk of limestone spalling analyzed.

Thermo-mechanical model
A 2D transient thermo-elastic model with weak coupling between thermal and mechanical considerations is used. During the iterations, the heat diffusion equation (Fourier) equation is solved, and then thermal deformations of the solid are calculated. Numerical resolution is performed by the finite elements software "CAST3M" [31]. The thermal and mechanical behaviour of limestone are considered to be isotropic. The simulations were carried out over a period of 6 h (initial temperature in limestone of 20 C), including a short cooling phase of the rock.

Geometry
The cavity was modelled in two dimensions, i.e. only a representative cross-section of the wall was meshed (see cross section "A-A" on Fig. 16).
The meshed area is large enough to be free of thermal side effects. The first centimetres of the wall were more finely meshed due to the high thermal gradients occurring. Part of the ceiling and the floor of the cavity were also meshed. The area was meshed with rectangular surface finite elements, with a total of 12 700 finite elements. A view of the mesh is presented in Fig. 16.

Thermal properties
The thermal properties of the limestone at the initial temperature (20 C) are those presented in x3.3. For thermo-mechanical simulations, the variation of the thermal properties with temperature was taken into account. Based on the results of Bazant et al. [32], Vosteen et al. [33] and Wang et al. [34], the values reported in Table 2 were used.
The thermal conductivity of the limestone decreased due to drying during heating and to the formation of micro-cracks (thermal degradation). The density of the limestone decreased very weakly because the temperatures studied remained below that of limestone decarbonation i.e. above 650 C. Finally, based on the method proposed by Bonacina et al. [35], a peak of specific heat (7232 J kg À1 . K À1 ) was introduced between 100 and 120 C to take into account the energy consumed by the vaporization of the initial moisture of limestone.
Simulations were conducted to examine the sensitivity of the calculated temperatures depending on the thermal properties of limestone, by varying the initial thermal diffusivity given in Table 2. The relative change of the thermal diffusivity with temperature remains the same as that given in Table 2. Fig. 17 presents the evolution of the "temperature ratio" defined as the temperature calculated with a modified thermal diffusivity divided by the temperature calculated with the initial thermal diffusivity ( Table 2). This ratio is calculated for a thermal diffusivity 20% higher and 20% lower than the initial value. The calculations are performed at 1 and 3 cm deep in the wall close to the fire, at 88 cm from the ground.
Under the thermal exposure of this study, uncertainty of thermal diffusivity involves weak modification of the temperature in areas of interest. Indeed, the limestone rubification is observed after a sufficiently long time (around one hour of fire, see Fig. 18), and essentially in the first 3 cm of rock (Fig. 19). At 88 cm from the ground, and after 1 h of fire, the temperature changes only more or less 5e10% for an uncertainty of thermal diffusivity of more or less 20% (Fig. 17).

Mechanical properties
Compression tests were carried out on cylindrical cores (diameter 45 mm, length 60 mm) drilled from a wall of the experimental cave. Four compression tests were performed at room temperature. An additional compression test was performed on a specimen heated to 300 C (at 5 C/min) then naturally cooled to room temperature. 300 C corresponds to the maximum heating capacity of the lab's furnace but also matches the order of magnitude of the temperature reached in the first centimetres of the limestone during the experimental fire (see Fig. 18).
Experimentally measured values of the modulus of elasticity were taken into account in the calculation model. The value at 500 C is considered equal to the value measured on the specimen heated to 300 C. Indeed, according to Mao et al. [36], the mechanical properties of limestone are little changed below the    . 17. Evolution of the temperature ratio at 1 (left) and 3 (right) cm deep in the wall close to the fire (at 88 cm from the ground), depending on the change of the initial thermal diffusivity. Temperature ratio is defined as the temperature calculated with more or less 20% of thermal diffusivity divided by the temperature calculated with the initial thermal diffusivity (Table 2). decarbonation temperature (650 C). For these simulations, we make the assumption that Poisson's ratio did not vary with temperature (ʋ ¼ 0.3). Finally, during the heating phase of the specimen at 300 C, the thermal expansion of the limestone was measured, allowing us to take into account the variation of the coefficient of thermal expansion with temperature. The value of the coefficient of expansion at 500 C is considered equal to that measured at 300 C. The average mechanical properties used for the simulations are reported in Table 3.

Thermal boundary conditions
Dirichlet type thermal boundary conditions are used, i.e. temperatures are imposed at the borders of the mesh. The variation of imposed temperatures was taken directly from FDS simulations (see x3). Depending on the border, imposed temperatures were as follows: -On the sole (Fig. 16, segment AB), a constant temperature of 100 C was applied. -The vertical wall (Fig. 16, segment BC) was partitioned into 18 segments: from the ground a first 30-cm segment (which corresponds to the position of the fire in the CFD simulation) then 10-cm segments until the ceiling. On each segment, the temperature simulated by the CFD model was imposed at the centre of the segment. -The ceiling (Fig. 16, segment CD) is partitioned into ten 10-cm segments. Temperatures simulated by the CFD model were imposed at the centre of each segment. -For all the other borders of the 2D mesh, no thermal boundary conditions were applied. Fig. 18. Simulated (grey curves) and experimental (black curves) temperature profiles in the first 10 cm of the limestone wall exposed to the fire (a e 20 cm from the ground, b e 90 cm from the ground, c e 160 cm from the ground) and at the ceiling (d e 30 cm from the wall, e À 70 cm from the wall).

Mechanical boundary conditions
In order to take into account the massiveness of the wall, segments EF, FG, GH and HA (Fig. 16) are totally restrained. The selfweight of the meshed rock and the dead weight of the rock above the experimental cave were not taken into account in the simulations. In this way, only the mechanical effects caused by fire are simulated.

Thermal results
Temperature profiles simulated in the first 10 cm of the limestone after 1, 2.5 and 5 h of fire are reported in Fig. 18. Profiles are given for the wall at 20, 90 and 160 cm from the ground and for the ceiling at 30 and 70 cm from the wall. The temperature profiles measured in the wall during the experimental fire are also presented.
Overall, the simulated temperatures were higher than the temperatures measured during the fire. This difference can come from the thermal boundary conditions. Indeed, temperatures simulated by the CFD model are higher than temperatures measured on the wall surface during the experimental fire (see x3. 5.4). The difference can also be explained by the uncertainty about the actual position of thermocouples in the wall and by the significant temperature difference observed between two points on the wall located at the same height from the ground but a few centimetres apart on the horizontal axis. Thermo-mechanical simulations being performed on a 2D model, this gap cannot be taken into account. However, the order of magnitude of thermal gradients over the first 3 cm of the limestone is well transcribed.
In the graphs of Fig. 18 we plot the rubification temperature of limestone that we consider equal to 250 C [6]. This allows us to estimate estimating the thickness of rubified limestone at different moments of the fire, and more particularly the residual rubification at the end of the fire. Table 4 presents the simulated thicknesses of rubified limestone after 1, 2.5 and 5 h of fire.
A front view of the wall exposed after the fire and evacuation of residual coals and ash in the hearth is shown in Fig. 19. Fig. 19 also presents the location of the small cores drilled for measurement of the rubification thickness. Simulated rubification thicknesses of limestone were greater than the measured thicknesses. This is mainly explained by simulated temperatures being higher than experimental temperatures (see x3.5.4). Nevertheless, the coupling between the CFD model of fire and the thermomechanical model can be adapted to the study of Paleolithic fires of the Chauvet-Pont d'Arc Cave. The comparison between the simulated thicknesses of limestone rubification and thicknesses Fig. 19. Left: picture of the wall directly exposed to fire after removal of the ashes, coals and half of the clay sole (Credit: C. Ferrier). The location of the drilled cores of limestone samples is indicated. Right, pictures of the drilled cores of limestone (Credit: A. Brodard). The depth of rubification is given depending on the height above the floor of the chamber.

Mechanical results
Fig . 20 shows the evolution of the mechanical stresses induced in the limestone due to thermal diffusion.
We present the normal stress along the vertical axis of the wall (y-axis in Fig. 15). Compressive stresses are considered negative. Under the effect of heat diffusion in the wall, and due to the fact that the coefficient of expansion of limestone is positive (dilatant material), it can be noted that: -Thermal expansion of the lower part of the wall, which receives more heat, is greater than in upper part of the wall.
-At the same height, the limestone close to the surface exposed to fire (hottest) expands more than limestone in deeper areas (colder).
The consequence of the restrained thermal expansions (along the height and thickness of the wall) is the generation of compressive stresses along the vertical axis of the wall (y-axis in Fig. 15). Note that in a real wall, compressive stresses along the zaxis are probably also generated (not simulated in our 2D working model). Thus the first few centimetres of the real wall are actually in a biaxial compression state.
As reported in Fig. 20, compressive stresses in the first 5 cm of limestone increased steadily until the fire cooled down. In particular, in the first centimetres, compressive stresses increased very quickly during the first hour of the fire. Compressive stresses decrease with the depth of rock due to the decrease of the thermal gradients. Similarly, compressive stresses decrease up the height of Fig. 20. Vertical compressive stress versus time at 1, 3 and 5 cm into limestone at the wall exposed to the fire (a e 20 cm from the ground, b e 90 cm from the ground, c e 160 cm from the ground) and at the ceiling (d e 30 cm from the wall, e À 70 cm from the wall). the wall and are very weak at the ceiling.
In the lower part of the wall, the compressive stress reaches 12 MPa at 1 cm in depth and 8 MPa at 3 cm in depth. Rupelien limestone compressive strength was measured at room temperature on cylindrical specimens (diameter 45 mm, length 60 mm) taken from a wall of the experimental cave and after heating at 300 C. The values obtained are respectively 1.96 MPa and 2.33 MPa. We therefore find that the compressive stresses induced by the fire are higher than the compressive strength of limestone. Consequently, fires such as the one modelled can cause localized damage to the limestone, especially in the first 3 cm of exposed wall, and especially at the base of the wall (where the heat was greatest during the experiment). These numerical results are in good accordance with observations of the wall after the experimental fire: flakes around 2 cm thick were observed coming off the lower part of the wall (Fig. 21).
They were probably induced by mechanical buckling of the first centimetres of limestone, i.e. where compressive stresses are high. These assumptions are generally accepted for explaining granitic rock spalling during "blast hole drilling" [37] or concrete spalling due following exposure to fire [38]. Spalling probably induced by Paleolithic fires is also observed in the Chauvet-Pont d'Arc cave [5].

Conclusions and perspectives
The simulation approach was validated on a fire in an underground quarry in order to apply it to the Chauvet-Pont d'Arc cave. The aim was to study a geometry similar to that of the Megaceros Gallery, where numerous thermal marks are found on the walls. The main objective was to couple a Computational Fluid Dynamics code and a thermo-mechanical code. With an archaeological objective, the simulation of a fire includes the effects of the combustion in the air and its impact on the rock. Wood combustion in a volume like the one of the quarry chamber was simulated with the software Fire Dynamics Simulator (FDS). The Heat Release Rate of the fire was set with respect to experimental temperatures, as was its behaviour over the time of the fire. In the air, variations of temperature, velocities, and carbon monoxide levels were calculated during the combustion. Temperatures at the wall were validated with the data acquired during the experiments and the differences were discussed. Simulation was useful to go further with analysis of the distribution of the combustion gases. The Fractional Effective Dose delimited the zones where it was less dangerous to stay during the combustion. Up to 1 m high, it was possible to circulate around the fire, making it feasible to add fuel to the fire during the combustion provided the people moved in this lower zone. Besides, the experimental quarry was chosen for its dimensions similar to the Chauvet-Pont d'Arc Megaceros gallery, and the height of the smoke ceiling corresponds to the soot marks around the gallery [11].
Starting from the variation of the simulated wall temperatures, the thermo-mechanical impact on the rock was then simulated with the CAST3M software. Heat diffusion analysis in the rock enabled the thickness of rubified limestone to be simulated. Secondly, thermo-mechanical analysis was used to localise the most probable zones for limestone spalling due to thermal gradients.
In the experimental quarry, post-fire observations and numerical simulations showed that areas of rubification and spalling of limestone were located at the bottom part of the wall where thermal gradients were the steepest. It is thermal gradients induced by fire that are likely to have been the cause of limestone spalling, with flakes around 2 cm thick. These observations were correctly simulated.
In parallel, observations in the Chauvet-Pont d'Arc cave showed that rubification and spalling of limestone mainly occurred in the upper parts of the walls which clearly indicates that the Paleolithic fires were larger than the experimental fires made in this study. In addition, the flakes attributed to fire in the Chauvet Cave are a lot thinner (0.5 cm thick on average) than those found after the experimental fire in the quarry (2 cm thick). This is likely due to higher thermal gradients in the rock during the Paleolithic fires, inducing higher compressive stresses that are more confined to the surface of the rock exposed to fire. Coupled with the combustion model, the thermomechanical tool is able to quantify the thermal gradients in the rock and therefore allows estimation of flake thickness depending on the fire intensity. Therefore, the coupled model developed is useful for judging the intensity of the fires used in the Paleolithic caves.
Our ultimate goal is archaeological, now that our simulations have been validated by the experiments in the quarry, the approach will be applied in the Chauvet-Pont d'Arc cave. The aim is to provide hypotheses on the location of the hearths, starting from the positions of the fire marks on the walls. Several simulations will be performed to find the original locations of the fireplaces and the intensity of the fires having led to the thermal marks observed on the walls. Indeed, the fire modifies the air dynamics and the simulation will allow refinement of the position of the fires. The fire marks are not necessarily located exactly above the position of the hearths since, depending on the presence of draughts, they can be shifted.
The "CarMoThap" project (Characterization and Modelling of the Thermo-alterations and the residues of combustion on the walls, funded of the Aquitaine Region) will begin during the second half of 2015. This project focuses on the study of thermal marks found in Paleolithic caves, especially the Chauvet-Pont d'Arc cave, and for which the function of the fires has not yet been identified. The "CarMoThaP" program is based on an interdisciplinary approach, based on data acquired by experimental fires and results of aero-thermo-mechanical simulations. One of the projects is to make new fires in the underground quarry in Lugasson, sufficiently large to generate thermal marks in the upper part of the walls, such as those observed in the Chauvet-Pont d'Arc cave. The coupled numerical tool presented in this paper will help in the design of the experimental protocol, in terms of the position and intensity of the fires to be built. The coupled approach with experiment and simulation mutually nourishing each other is complementary to the acquisition of observation data in the cave in order to better understand the function of the fires.